Light Dose and Fluorescence Imaging Depth in Dual-Wavelength PDT: a Numerical Study for Various Photosensitizer Distributions in a Layered Biotissue
Daria Kurakina*, Ekaterina Sergeeva, Aleksandr Khilov, and Mikhail Kirillin
A.V. Gaponov-Grekhov Institute of Applied Physics of the Russian Academy of Sciences, 46 Ulyanov str., Nizhny Novgorod 603950, Russian Federation *e-mail: [email protected]
Abstract. Dual-wavelength photodynamic therapy (PDT) in combination with fluoescence monitoring has great potential in medical applications, in particular, in the treatment of skin diseases. In this study we conducted numerical Monte Carlo experiments in a two-layer skin model to simulate the absorbed light dose distribution for PDT with chlorin-based photosensitizers as well as to estimate the fluorescence imaging depth depending on the probing radiation wavelength and the photosensitizer (PS) in-depth distribution. The obtained dependencies provide a detailed overview of the light absorption profile and fluorescence signal formation for different PS distribution in layered biotissues, and can be used in various fluorescence based light dosimetry techniques. © 2024 Journal of Biomedical Photonics & Engineering.
Keywords: chlorin-based photosensitizers; photodynamic therapy; fluorescence imaging; light transport in biotissue; Monte Carlo simulations.
Paper #9186 received 25 Oct 2024; revised manuscript received 27 Nov 2024; accepted for publication 29 Nov 2024; published online 27 Dec 2024. doi: 10.18287/JBPE24.10.040318.
1 Introduction
Photodynamic therapy (PDT) is a modern non-invasive therapeutic method used for treatment of a wide range of tumor and non-oncologic pathologies via externally stimulated cytotoxic effects, tumor vascular damage, as well as immune responses [1-3]. PDT procedure is based on the administration of a photosensitizer (PS) and its accumulation in the target tissue followed by illumination by light of an appropriate wavelength causing production of cytotoxic agents [1-3]. The success of a PDT procedure depends on local concentration of cytotoxic agents which is affected by the PS concentration, light transport in the tissue and tissue oxygenation [4]. According to modern PDT protocols, only the PS dose administered to a patient and the light energy delivered to the surface of the biological tissue during the procedure are specified at the stage of the treatment planning [4]. This limited set of parameters does not consider such individual features as heterogeneous PS accumulation and non-uniform distribution of optical radiation in biological tissues of various localizations. Thus, in regard with the current trends towards personalized medicine, individualization
of the PDT procedure scheme at the stage of its planning is an urgent problem.
For PSs with high absorption in different bands of visible spectrum, PDT can be performed at several wavelengths. Chlorin-based photosensitizers are a wide class of PDT agents with two pronounced absorption peaks located in the blue and red bands of the visible spectrum (402 and 662 nm, respectively). Optical properties of biotissues at these therapeutic wavelengths differ significantly [5]: for most biotissues subjected to PDT (superficial and mucous biotissues) both the absorption and scattering coefficients at X = 402 nm are several times higher than those at X = 662 nm [5-7], which leads to different penetration depth of therapeutic light for these two wavelengths and thus, to the difference in the PDT impact area. Therefore, development of multi-wavelength PDT protocols may benefit from combining the effects of exposure at several wavelengths to ensure the effect over wide depth.
Treatment of skin diseases is an essential part of PDT applications with specific requirements for light dosimetry. Skin is a layered structure with epidermal, dermal and hypodermal layers that are characterized by
specific values of thickness and optical properties. These characteristics can significantly deviate for diseased versus normal skin [8]. In the case of non-uniform tissue, determination of light fluence distribution in biotissue is a challenge. Light dosimetry problem for PDT has been studied previously using analytical or numerical approaches, however, they are mainly limited by the case of uniform biotissue [9, 10]. Direct interstitial measurements of light fluence rate [11, 12] are invasive and low informative in terms of providing volumetric information. Monte Carlo simulation technique [13] is the most promising approach for modeling the propagation of light in biotissue with complex structure. The simulation consists in tracing the large number of random photon trajectories for a given probing configuration in the medium with predefined optical properties followed by statistical processing of the obtained data. This approach makes it possible to estimate PDT impact volume for the specified optical properties of the medium, tissue geometry, and the procedure mode, by calculating the absorbed light dose distribution. In Ref. [14-15] numerical Monte Carlo simulations of fluence in multilayer skin are carried out for one [14] or several [15] wavelengths followed by its implementation into PDT pharmacokinetics mathematical model to evaluate the production of singlet oxygen. In Ref. [16], using Monte Carlo simulation of a PDT procedure, the effect of tumor density heterogeneity on the distribution of light and photodynamic dose is investigated. In paper [17], based on Monte Carlo simulation of the irradiation distribution, dose-volume histograms are obtained which demonstrate the fraction of the organ volume receiving particular level of light dose.
In order to improve the efficiency of a PDT procedure and to obtain additional information about the tissue in the course and in the follow-up of the procedure, optical imaging techniques are used [3, 18]. Fluorescence properties of most PSs make it possible to employ fluorescence imaging (FI) techniques to monitor PS distribution in tissues [19-21]. Broad-beam excited FI allows for monitoring the accumulation of PS in biological tissue prior to the PDT procedure, as well as evaluating the PS photobleaching during PDT [22, 23] that can serve as a dosimetric indicator of the procedure effectiveness [24-28]. In monitoring the PDT procedure with the chlorin-based PSs characterized by two absorption peaks, dual-wavelength fluorescence imaging allows to analyze information from different depths simultaneously due to different localization of fluorescence excited at 402 nm and at 662 nm [19-21], while the thickness of the PS distribution can be estimated by the ratio of fluorescence signals at two excitation wavelengths [29, 30].
At the same time, the question arises whether the depth from which the fluorescence is registered is the depth where it has been emitted. Due to isotropic fluorescence emission and further stochastic behavior of light in biotissue, the depth from which the fluorescence is registered in explicit dosimetry PDT systems is challenging to predict and it cannot be measured directly.
The maximum depth from which the fluorescence photons are detected (further called "FI depth") depends on the optical properties of the tissue and on the way of PS administration, being an important characteristic for further optimization of the fluorescence monitoring technique. Monte Carlo approach is widely utilized for modeling fluorescence distribution in biological media [5, 31-35]. It is a convenient tool for tracing individual trajectories of the detected fluorescence photons and it can be used to analyze their depth distribution. In Ref. [35] numerical investigation of fluorescence probing depth dependence on the source-detector separation is conducted for a configuration of a surface contact probe equipped with multiple detector fibers. However, studies of signal formation in fluorescence imaging with several excitation wavelengths are quite limited [19-21, 36, 37].
The aim of this work is to perform numerical experiments in a two-layer biotissue mimicking skin with epidermal and dermal layers to simulate the distribution of the absorbed light dose for PDT with chlorin-based photosensitizers, as well as to estimate the FI depth in dependence on the probing radiation wavelength and the PS distribution pattern.
2 Materials and Methods
2.1 Monte Carlo Simulations 2.1.1 Distribution of the Excitation Light Absorption
In this study we used our previously reported code [20] which was adopted and employed for simulation of signal formation in a system for widefield surface FI [29]. The feature of this code is that all photons are processed in parallel in contrast to a conventional MC algorithm with consequent calculations. Biotissue is considered as a slab of scattering and absorbing medium which is illuminated by a wide beam from its upper boundary z = 0 along z direction. The slab consists of two layers with differing optical properties: thin topical layer mimicking epidermis and thick layer mimicking dermis. To calculate a three-dimensional distribution of the absorbed dose A(x,y,z,Xex) of light with the wavelength Xex, at each n-th step of modeling and for each photon the absorbed fraction of the photon weight was added to the corresponding voxel of the absorption map array:
An(x,y,z,Aex) = An-1(x,y,z,Aex) +
+ YW _^a(x,y,Z,ÀeX)_
1 Va O, 7, Z, Aex) + Vs O, 7, Z, Kx) '
(1)
where Wi is the current weight of the i-th photon located in a voxel with the center coordinates (x, y, z) at step n, Ha (x, y, z, Xex) and ns (x, y, z, Aex) are local absorption and scattering coefficients of the medium in this voxel. Eq. (1) implies that all the photons are processed in parallel and at each step of the simulations the absorbed weight fractions of all the propagating photons are
summed up within a given voxel. Local optical properties account for contribution of a background biological tissue (epidermis or dermis) and of a photosensitizer. We assumed that the background biotissue was transversely homogeneous with parameters HaSSUe(z,Xex) and ^tissue (z,xex), while non-uniformly distributed PS was manifested by the additive term to the local absorption coefficient as in Refs. [29, 38, 39]:
ßa(x,y,z,Äex) = ^PS(x,y,z,Äex) + ^Jissue(z,Xex); (x, y, z, Aex) = ßls'ssue (z , iex).
(2)
Here ^pas(x,y,z,Xex) = oaPS(Xex)CPS(x,y,z) is the absorption coefficient of PS with the absorption cross-section oaPS{Xex) and volume concentration distribution Cps(x,y,z).
To simulate the formation of fluorescence signal, we employed previously developed dual-stage algorithm [20, 29, 40]. At the first stage of Monte Carlo simulation, the map of absorbed light dose A(x, y, z, Xex) was calculated for the considered optical properties, particular PS distribution in the tissue and the irradiation wavelength. Further this map was recalculated to generate the light dose density Aps(x,y, z,Xex) absorbed directly by the PS using the ratio of local PS absorption coefficient to the total one:
Aps(x,y,z,Aex) = = A(x,y,z,Aex)
ß%S(x,y,z,Aex)
(3)
ßFas(xy,z,Äex)+ß^e(z,Äexy
The values of A (x, y, z, Xex) and Aps (x, y, z, Xex) are measured in photons. To transform them into volumetric distributions of fractional absorbed dose one should normalize the corresponding characteristics by the value of A0 = Nph AV [photonsmm3] where Nph is the number of excitation photons and A V = AxAyAz is the volume of an elementary voxel with the sizes Ax, Ay, and Az.
For PDT planning, it is more practical to analyze one-dimension in-depth distribution of energy absorbed by a PS only (Eq. (3)) that allows to evaluate the impact depth for PDT at the particular irradiation wavelength. It can be calculated by transversal integration of the normalized PS absorption dose within an elementary layer at the depth z:
Aps(z,Xex) = J
Aps(x,y,z,Aex)
■ex) / ,
A0
x,y
AxAy =
NphAzßps(z, Äex) + (z, Äex)
~Y^^A(x,y,z,Xex)
(4)
and is measured in inverse millimeters. Another key parameter that cannot be measured directly is the dimensionless fraction of energy Aps(Z,Xex) absorbed by the PS in a layer of a thickness Z measured from the surface of the tissue which is calculated by the integral:
Aps(Z,Aex) =YjAps(z,Aex)Az.
(5)
2.1.2 Distribution of Fluorescence
At the second stage the PS absorption map Aps (x, y, z, Xex) served as a distribution of fluorescence emission sources. The emission of fluorescence photons was assumed isotropic with quantum yield equal to 1 for simplicity. As we are interested not in absolute values of fluorescence but in weights distribution, the magnitude of quantum yield is not critical. In the evaluation of absolute fluorescence, real value of quantum yield amounting 0.1-0.2 [41] for chlorin-based PSs can be easily introduced by the corresponding multiplicative factor. Surface distribution of fluorescence intensity Ifi(x,y,z = 0,Xex,Xem) [1/mm2] can be calculated as:
If i(x,y,Z = 0,Xex,Xem) =
= W^S X 1em(x',У',Z',X,У,Z = 0,Xex,Xem), ^
ph x',y',z'
where Iem(x',y',z',x,y,z = 0,Xex,Xem) is the weight of fluorescence photons with the wavelength Xem that are launched isotopically from the center of an elementary voxel with coordinates (x',y',z') having the initial weight equal to Aps (x', y', z', Xex) and leave the biotissue through an elementary pixel AS = AxAy with coordinates (x, y) at the tissue surface z = 0.
Along with the fluorescence signal from a distributed source, the following parameters for all detected fluorescence photons were recorded: the maximum depth zjnax reached by j-th fluorescence photon in the medium ("FI depth") and the photon weight Wj when it exits the biotissue. The weight distribution over the FI depth Wd(zmax) for photons detected at the boundary z = 0 from all directions and collected at the whole sample surface was binned according to the following Eq.:
Wd(Zmax) =^Wj {zJ™*€(zmax,zmax+Az]}, ] = 1
zmax=kAz, k = 0,...Nk-1,
(7)
where each k-th bin contains the summed weight of the photons with FI depth within the bin edges
+ Az]
N is the total number of
fluorescence photons that have reached the tissue surface, Az is the bin size, Nk is the total number of bins. The effective FI depth was calculated as the median value zmax of the corresponding distribution Wd(zmax).
2.2 Simulation Parameters
We aim to model light distribution during a PDT procedure of various skin diseases. There are two ways in PS administration in this task. Topical application is used in lesions that affect mainly epidermis so therapeutic impact is localized above basal membrane which is almost impenetrable for a PS. In the therapy of malignant lesions deeper impact is required so intravenous injection of a PS is implemented according to the protocol. We further considered three models of in-depth PS concentration
1
x,y
z=0
distribution CPS(z) in biological tissue mimicking skin (Fig. 1) similar to those previously discussed in paper [29]: (i) boxcar PS distribution in the upper tissue layer of a given thickness dE (a simplified topical application model assuming the PS is uniformly distributed over the whole epidermis thickness dE); (ii) uniform PS distribution in a subsurface tissue covered by a PS-free layer of a given thickness dE (an intravenous injection model assuming PS accumulation below epidermis due to blood flow); (iii) exponential in-depth decay of PS concentration in a two-layer tissue (a refined topical application model assuming realistic PS penetration through upper layers, i.e. stratum corneum and epidermis).
In distributions (i) and (ii), the PS absorption coefficient corresponding to the constant level of the drug concentration was denoted as nPs(Aex). In distribution (iii), the PS absorption coefficient profile was set in accordance with the following Eq.:
№p,exp Aex)
Mo(A-ex)
d
l/e
exp
V d1/e)
(8)
where M0(Aex) = oPs(Aex) • Z is a dimensionless parameter equal to the product of the PS absorption cross-section (^ps(Aex)) and the PS surface density (Z) at the moment of its application on the tissue, d1/e is a spatial scale of PS concentration decrease with depth. This model is more accurate approximation for the topical application case [42, 43] which allows for tracking the dynamics of the PS penetration if the dependence of d1/e on time is given.
In order to study spatial distribution of the absorbed light and the FI depth distribution in the task of dual-wavelength PDT for all PS concentration models (further denoted as Case (i), Case (ii) and Case (iii), correspondingly), we conducted simulations of light propagation at two excitation wavelengths
Aex = 405 nm and 660 nm typically used in PDT with chlorin-based PSs. Fluorescence from chlorin-based PS in dual-wavelength systems is typically detected at the wavelengths red shifted from the emission peak at 665 nm to avoid a cross-talk with the reflected red-excitation light [20]. In our simulations the PS fluorescence was modeled for the emission wavelength of 760 nm equal to that in the custom-designed dual-wavelength FI setup [20].
Total amount of Nph = 107 excitation photons was used in the simulations. The sample was illuminated from the top normal to the surface by a circular beam with radius of 9 mm and uniform intensity distribution. Transverse dimensions of the sample were set to 20 mm x 20 mm, its thickness was 20 mm, and the size of an elementary voxel equaled 0.25 mm x 0.25 mm x 0.01 mm. The value of dE in the distributions (i) and (ii) varied in the range of 0.05-1 mm in accordance with the data on epidermis thickness in norm and pathology in different skin types [44-47], the value of d1/e in the distribution (iii) varied in the range of 0.1-1 mm. Optical properties of epidermis and dermis were set in accordance with the data in Ref. [8] and are presented in Table 1 for the considered excitation and emission wavelengths (405 nm, 660 nm, and 760 nm). The values of the PS absorption coefficient nPs (A) in Cases (i) and (ii) corresponded to the constant concentration level of Revixan® Derma of 0.1% vol that was achieved by in vivo topical administration of 1 |g/g [48, 49]. The value of M0 (Aex) was chosen so that the surface absorption coefficient npsexp(0,Aex) = M0(Aex)/d1/e was of the same order as the corresponding nps(A) for the first two cases. The value of scattering anisotropy factor was set as g = 0.8 for all wavelengths as we assumed that the effects of small-angle scattering related to the shape of a scattering phase function were insignificant.
Fig. 1 Considered spatial PS distribution models with different in-depth profiles: topical application model (Case (i)), intravenous injection model (Case (ii)), refined topical application model (Case (iii)). Here dE is the thickness of epidermis, d1/e is the 1/e scale of exponential PS concentration distribution.
Table 1 Optical properties of epidermal (E) and dermal (D) skin layers and chlorin-based PS used in simulations.
Wavelength, nm Ma5, mm-1 M0 MEa, mm-1 «E, mm 1 mm-1 mm 1 g
405 0.10 0.020 1.25 50 0.9 38 0.8
660 0.02 0.004 0.26 22 0.14 14 0.8
760 0 0 0.2 14 0.13 12 0.8
3 Results and Discussion
3.1 Absorbed Light Dose Distribution
Fig. 2 shows log10 of normalized absorbed light dose distribution A(x,y = 0, z, Aex)/ A0. The maps are plotted for a boxcar PS distribution (i) with the thickness of dE = 0.2 mm and two irradiation wavelengths Aex of 405 nm and 660 nm. Due to significant dependence of skin optical properties on the wavelength (see Table 1), for Aex = 405 nm the light is primarily absorbed in superficial tissue layer, while for Aex = 660 nm deeper layer is involved. Central part of the maps represents distribution of absorbed dose over depth which can be used for planning PDT procedures with chlorin-based PS in skin.
Fig. 3a shows in-depth profiles of the PS-absorbed dose Aps(z,Aex) calculated with Eq. (4) for different thicknesses dE of epidermis containing uniformly distributed PS and for two excitation wavelengths Aex of 405 nm and 660 nm (Case (i)). Smaller values of dE correspond to normal thickness of epidermis while larger ones, to the cases of its pathological thickening, e.g., for actinic keratosis [44, 45]. Fig. 3b depicts the
dependences of absorbed energy fraction Aps(Z,Ae
calculated with Eq. (5) for the same PS distribution. Faster accumulation of total light dose delivered with blue light compared to that for red light is shown.
Fig. 2 Logio of normalized volume density of the absorbed energy fraction A(x,y = 0,z,Aex)/A0 in the biotissue containing PS uniformly distributed in epidermis (dE =0.2 mm) for irradiation wavelengths of (a) 405 nm and (b) 660 nm. White dashes show the boundary of the PS containing layer.
(a)
@405 nm, d = 0.1mm @405 nm, d = 0.25 mm @405 nm, d = 0.5 mm @405 nm, d = 1 mm nm, d = 0.1 mm nm, d = 0.25 mm nm, d = 0.5 mm nm, d = 1mm
0.4 0.6 0.8
z, mm
(c)
5-, ■
4-
(D II
3-
K
2-
K
1-
(b) —
0.06-
@405 nm, d = 0.1mm @405 nm, d = 0.25 mm @405 nm, d = 0.5 mm @405 nm, d = 1 mm
nm, d = 0.1 mm nm, d = 0.25 mm nm, d = 0.5 mm nm, d = 1mm
co"
Q- W
K
0.04-
0.02-
0.00
0
0.0 0.2 0.4 0.6 dE, mm
0.8
1.0
Fig. 3 (a) In-depth distribution of energy Aps(z,Aex) absorbed by a PS uniformly distributed in a topical (epidermal) layer of a thickness dE for the irradiation wavelengths of 405 nm and 660 nm; (b) fraction of energy Aps(Z,Aex) absorbed by a PS in a superficial layer of tissue versus the layer thickness Z. (c) Ratio of total light doses absorbed at Aex = 405 nm to that absorbed at Aex = 660 nm as the function of the PS layer thickness.
uaPS (1=405nm)/ii:S(1=660nm)
For the given PS concentration total absorbed energy fraction does not exceed units of percent. Similar to the plots in Fig. 3a, Fig. 3b shows that the dose accumulation stops at the depth Z equal to the thickness of the PS-containing layer dE. Note that the ratio of total doses absorbed at different wavelengths for the same PS distribution (Fig. 3c) is not equal to the ratio of the PS absorption coefficients at these two wavelengths
.^%s(405 nm)
(_
PaS(660 nm)
= 5) and decreases with the increase of the
PS layer thickness.
In Fig. 4 density profile of the PS-absorbed light dose and absorbed energy fraction are shown for Case (ii) mimicking intravenous injection with a uniformly distributed PS in a subsurface slab of biotissue covered by a PS-free layer (epidermis) of thickness dE. The absorption profile (Fig. 4a) shows a jump at the depth equal to dE while the curve shape is determined by the optical properties of dermal tissue. Fig. 4b demonstrates the dependences of absorbed dose Äps(Z,Äex) for
(a)
N
0.20
0.15
0.10
0.05
0.00
@405 nm, dE = 0.1mm @405 nm, dE = 0.25 mm @405 nm, dE = 0.5 mm @бб0 nm, dE = 0.1 mm @бб0 nm, dE = 0.25 mm @бб0 nm, dE = 0.5 mm @бб0 nm, dE = 1mm
Case (ii) showing the behavior of the curves for 405 nm and 660 nm similar to Case (i). Absorbed dose of light at the wavelength of 660 nm is comparable or even exceeds that at the wavelength of 405 nm except the case of the smallest dE value.
The same set of curves is plotted for the case of exponential PS distribution in the two-layer biotissue (Case (iii)) for different values of d1/e (Fig. 5). The increase of d1/e mimics the dynamics of in-depth redistribution of a PS after topical application. During this process the fraction of the PS in superficial layers decreases and its fraction in deeper layers grows, however total amount of PS is constant. Despite that, total absorbed light doses that correspond to the asymptotic values in Fig. 5b decrease with d1/e, which is caused by light attenuation. Also, the accumulated doses reach the asymptotic values faster (Fig. 5b) as compared to Case (i) with uniform PS distribution within epidermis (Fig. 3b).
(b)
N
S
CL W <
0.0б
0.04
0.02
0.00
@405 nm, dE = 0.1mm @405 nm, dE = 0.25 mm 0.5 mm 0.1 mm
1.0 z, mm
1.0 Z, mm
Fig. 4 (a) In-depth distribution of light dose Aps(z,Xex) absorbed by a PS uniformly distributed below a topical (epidermal) layer of a thickness dE for the irradiation wavelengths of 405 nm and 660 nm; (b) fraction of light dose (Z, Xex) absorbed by a PS in a superficial layer of tissue versus the layer thickness Z.
(a)
(b)
0.5 0.4 ^ 0.3
N
< 0.2
nm, d1/e = 0.1mm nm, d1/e = 0.25 mm nm, d = 0.5 mm nm, d = 0.1 mm nm, d = 0.25 mm nm, d,,„ = 0.5 mm
0.4 0.б z, mm
N
CL W <
0.0б
0.04
0.02
0.00
@405 nm, d1/e = 0.1mm @405 nm, d1/e = 0.25 mm
■ @405 nm, d1/e = 0.5 mm @бб0 nm, d1/e = 0.1 mm @бб0 nm, d1/e = 0.25 mm
■ @бб0 nm, d1/e = 0.5 mm
0.4 0.б Z, mm
Fig. 5 (a) In-depth distribution of light dose Aps(z,Xex) absorbed by a PS distributed exponentially with different values of d1/e for the irradiation wavelengths of 405 nm and 660 nm; (b) fraction of light dose Aps(Z,Xex) absorbed by a PS in superficial layer of tissue versus the layer thickness Z.
3.2. Fluorescence Imaging Depth
Based on the PS fluorescence excitation at two wavelengths from different spectral ranges, dual-wavelength FI provides the capability to measure the PS localization depth using the ratio of the detected fluorescence signals [19, 20, 29, 38]. In order to analyze the fluorescence signal formation in dual-wavelength FI, one should figure out from which volume the signal is detected for the considered excitation wavelength. This characteristic is statistically related to the maximum depth (FI depth) reached in the medium by the photons which contribute to the detected fluorescence.
Fig. 6 shows distribution of fluorescence photon weights over maximum depths Wd(zmax) calculated with Eq. (7) for all three models of PS distribution (Cases (i)—(iii)) at two excitation wavelengths. The widths of the distributions for Aex = 405 nm and Aex = 660 nm noticeably differ which is explained by a much smaller number of absorbed photons that generate fluorescence at a large depth for blue light compared to red light. As a result, the distribution maxima for Aex = 405 nm are located closer to the tissue surface than those for Aex = 660 nm.
(a)
8X10"4-!
. 6X10-4-
=3
ra
ro 4X10-4
I
2X10-4-0
(c)
4x10-4-=3 3x10-4-
£ 2X10-M
1x10-4 -
(b)
(e)
0.05 mm 0.1 mm 0.25 mm 0.5 mm 1 mm
1.0
zmav, mm
0.1 mm 0.25 mm 0.5 mm 0.75 mm
1.0 zm„„, mm
0.1 mm 0.25 mm 0.5 mm 1 mm
3x10-
2x10-
■ 2x10-
1x10-
5x10-
(d)
2.5x10n
(f)
0.05 mm 0.1 mm 0.25 mm 0.5 mm 1 mm
1.0 z„a„, mm
0.1 mm 0.25 mm 0.5 mm 0.75 mm
0.1 mm 0.25 mm 0.5 mm 1 mm
Fig. 6 Distribution of fluorescence photon weights over FI depth Wd(zmax) normalized by the number of launched photons for excitation wavelengths of 405 nm (left column) and 660 nm (right column) and for different PS distributions: (a), (b) Case (i); (c), (d) Case (ii); (e), (f) Case (iii). Colored lines correspond to different values of the epidermis thickness dE (a)—(d), and the 1/e scale of exponential PS concentration distribution d1/e (e)—(f) indicated in the plots.
^ mm
zmax, mm
^ mm
In Figs. 6a,b for the Case (i), a sharp drop is observed in the distribution profiles at values of zmax = dE when the PS-containing layer thickness dE is small compared to the transport length of fluorescence light
(1 - g760)) amounting about 0.3—0.4 mm in epidermis and dermis. This fact is explained by the presence of two diffuse fluxes of fluorescence photon: a forward flux directed inside the tissue and a backward flux directed to the front boundary of the tissue where the detector is located. Photons from the forward flux can redirect into the backward flux as the result of backscattering and vice versa. For the backward flux, the maximum depth reached by most fluorescence photons corresponds to the depth at which fluorescence emission occurred, however, there is a fraction of those travelled forward and then backscattered. For the forward flux, the maximum depth of a photon is defined by the depth at which it has changed its direction towards the detector as a result of backscattering. Thus, photons from the backward flux mostly contribute to the bins with lower index corresponding to small values of m a x . The bins with larger index (i.e, larger zmax values) collect also from the forward flux due to backscattering. The overall distribution depends on the absorption profile within PS
(a)
1.0-,
E 0.8E
"S. 0.6-1
<u
"O
.1 0.4-
"O
<1J
0.2-
0.0
l = 405 nm IQR, l = 405 nm l = 660 nm IQR, l = 660 nm
0.0
0.2
0.4
0.6
0.8
dE mm
(c)
layer and the attenuation of fluorescence photons energy along their paths. The presence of a sharp drop at m a x = dE is caused by disappearance of the straight backward flux contribution since only the photons with at least one backscattering can form the fluorescence signal with zmax > dE. The number of such photons decreases with the increase of m a x which causes a decaying tail in the dependence. When the PS layer thickness exceeds lt760, the peak at the boundary zmax = d is smoothened due to pronounced multiple scattering when the forward and the backward fluxes are mixed.
For Case (ii) of uniform PS distribution below epidermis, the dependencies in Figs. 6 c,d show the opposite trend: minimum value of zmax is determined by the epidermis thickness dE . Multiple scattering is responsible for long tails in these dependencies whereas slower decay is observed for 660 nm, similar to Case (i) (Fig. 6b). In Case (iii) of exponentially decaying PS concentration, the distributions Wd(zmax) for various values of parameter d1/e in Figs. 6 e,f appear smooth with the maxima that are located close to the boundary for Aex = 405 nm (Fig. 6e) and are shifted to larger depths for Aex = 660 nm (Fig. 6f).
(b)
2.0-
E E
1.5-
<u
"O c
xs <u
1.0-
0.5-
0.0
l = 405 nm IQR, l = 405 nm l = 660 nm IQR, l = 660 nm
1.0
0.0
0.2
0.4
0.6
0.8
1.0
dE, mm
1.0
E 0.i
c
<u 0.6
"O a ro
0.4
0.2
0.0
l = 405 nm IQR, l = 405 nm l = 660 nm IQR, l = 660 nm
0.0
0.2
0.4 d.
1/e1
0.6 mm
0.8
1.0
Fig. 7 The dependence of the median FI depth zmax (solid lines) on the epidermis thickness dE ((a) Case (i) and (b) Case (ii)), and on the effective thickness d1/e ((c) Case (iii)) for excitation wavelengths of 405 and 660 nm. Filled areas between dashed lines indicate the interquartile range (IQR) between the 25th and the 75th percentiles of zmax distribution shown in Fig. 6.
Fig. 7 demonstrates the dependence of the median depth zmax on the characteristic scale of the PS distribution for Cases (i)-(iii). In Fig. 7a describing Case (i), the dependence of zmax on the thickness dE is not linear and it tends to an asymptotic value with the increase of dE due to decrease in the contribution of fluorescence photons from deep layers (see Figs. 6 a,b).
For X„x = 405 nm the median zl
reaches its
asymptotic value for the values of dE~ 0.75 mm, while for Xex =660 nm the median z££0x keeps growing. The ratio z^ax/zmax is close to 1 at small values of dE, however it increases with the epidermis thickness and for dE =1 mm it reaches 1.6 which characterizes the difference in the effective imaging depths in FI at these wavelengths. Fig. 7b shows the dependence of the median zmax on the epidermis thickness for Case (ii). It is almost linear for both considered wavelengths in the specified range of dE values. This is presumably due to the major role of the straight backward flux in the formation of fluorescence signal. However, in the whole range of epidermis thickness values the median value Zmax exceeds dE because it is influenced by the photons emitted from the back of the dermis, not only by those from the interface. Similar to Fig. 7a, the plot in Fig. 7c describing Case (iii) shows saturating growth of the median zmax with the effective PS layer thickness d1/e, whereas z^Oxl^max > 1 for a given value of d1/e. Still, the ratio z^ax/Zmlx is slightly less in Case (iii) compared to Case (i) owing to higher contribution of fluorescence coming from smaller depths.
3.3 Fluorescence Imaging Depth and Fluorescence Signal Ratio
Dual-wavelength FI has been reported as a tool for estimating the PS localization depth [29]. This
estimation is based on the monotonous dependence of
,660 1 fl
the fluorescence signal ratio Rx = Hjos on the
'fi
characteristic thickness of a PS-containing region. If60 and Ifl°5 are the fluorescence signals either registered from the tissue upon excitation at the wavelengths of 660 nm and 405 nm, correspondingly, or evaluated numerically with Eq. (6). However, for the development of reconstruction techniques it is reasonable to evaluate the relation between Rx and the FI imaging depth at both probing wavelengths, since this dependence should be more direct. The performed simulations allowed to plot the dependence of the ratio Rx on the median depth ratio z^O/Zm£x which characterizes the difference in the imaging volumes at two wavelengths employed in dual-wavelength FI. Fig. 8 demonstrates that the dependence of Rx on Zmax/Zmax is close to linear one for both cases mimicking topical administration case (Cases (i) and (iii)). Such a behavior is caused by a growing difference in FI depth at 660 nm versus 405 nm with the increase of the PS localization depth owing to the difference in the absorption dose in-depth distribution (see Figs. 3a,
5a). In Case (ii) the dependence monotonously decreases since both z^0 and z^x linearly grow with the thickness dE of the upper layer, however their difference becomes smaller (see Fig. 7b), while the ratio Rx increases due to larger attenuation of blue light compared to that of red light in the PS-free epidermal layer. The dependence is shown in semilogarithmic scale and can be attributed by an inverse power law.
—■— PS in episermis, Case (i) —a- PS in dermis, Case (ii) —•— exponential distribution, Case (iii)
U 0.4 -
~660 /~405 ^max ^max
Fig. 8
Dependence of fluorescence signal ratio Rx = Ifi60/lfi05 on the ratio z66a0x / zZ5x of median FI depths for different cases of PS distribution in two-layer biotissue.
4 Conclusion
In this study we have conducted numerical Monte Carlo modeling of light characteristics essential for planning and monitoring of PDT with chlorin-based PS. We compared the effect of different PS distribution models in skin on the profiles of absorbed light dose and on the FI depths in a two-layer biotissue mimicking epidermis-dermis geometry when excitation is performed at the wavelengths close to the maxima in the absorption spectrum of chlorin-based PS (405 nm and 660 nm). Profiles of PS-absorbed light dose evaluated with account of layered biotissue structure allow for more accurate estimation of the PDT procedure efficiency compared to the light dosimetry models developed for uniform biotissue. We studied fluorescence photon distribution versus the depth from which the signal is detected for different PS distributions in biotissue and for optical characteristics at particular probing and emission wavelengths. Dependence of the median FI depth on the characteristic depth of the PS distribution differs sufficiently for topical and intravenous PS administration cases: the former demonstrates asymptotic trend for both boxcar and exponential distribution of PS concentration, while the latter in nearly linear. Ratio of fluorescence signals at wavelengths of 660 nm and 405 in dual-wavelength FI as the function of the ratio of median FI depths at these wavelengths also shows opposite trends for topical and intravenous administration models: the former demonstrates slow monotonous increase while the latter
can be described by a decaying inverse power law. The obtained dependencies provide a detailed overview of the light absorption profile and fluorescence signal formation for different PS distribution in layered biotissues, which has practical value in various fluorescence-based light dosimetry techniques assisting dual-wavelength PDT.
Acknowledgements
The study has been supported by the Russian Science Foundation (Project 24-15-00175, https://rscf.ru/project/ 24-15-00175/ ).
Disclosures
The authors declare no conflicts of interest.
References
1. J. F. Algorri, M. Ochoa, P. Roldán-Varona, L. Rodríguez-Cobo, and J. M. López-Higuera, "Light technology for efficient and effective photodynamic therapy: a critical review," Cancers 13(14), 3484 (2021).
2. R. R. Allison, K. Moghissi, "Photodynamic therapy (PDT): PDT mechanisms," Clinical Endoscopy 46(1), 24 (2013).
3. J. P. Celli, B. Q. Spring, I. Rizvi, C. L. Evans, K. S. Samkoe, S. Verma, B. W. Pogue, and T. Hasan, "Imaging and photodynamic therapy: mechanisms, monitoring, and optimization," Chemical Reviews 110(5), 2795-2838 (2010).
4. B. C. Wilson, M. S. Patterson, "The physics, biophysics and technology of photodynamic therapy," Physics in Medicine & Biology 53(9), R61-R109 (2008).
5. V. V. Tuchin, Tissue Optics: Light Scattering Methods and Instruments for Medical Diagnostics, 3rd Ed., Society of Photo-Optical Instrumentation Engineers (SPIE), Bellingham, WA, USA (2015). PDF ISBN: 9781628415179.
6. A. N. Bashkatov, E. A. Genina, V. I. Kochubey, and V. V. Tuchin, "Optical properties of human skin, subcutaneous and mucous tissues in the wavelength range from 400 to 2000 nm," Journal of Physics D: Applied Physics 38(15), 2543-2555 (2005).
7. S. L. Jacques, "Optical properties of biological tissues: a review," Physics in Medicine & Biology 58(11), R37-R61 (2013).
8. E. Salomatina, B. Jiang, J. Novak, and A. N. Yaroslavsky, "Optical properties of normal and cancerous human skin in the visible and near-infrared spectral range," Journal of Biomedical Optics 11(6), 064026 (2006).
9. S. L. Jacques, "How tissue optics affect dosimetry of photodynamic therapy," Journal of Biomedical Optics 15(5), 051608 (2010).
10. J. S. Dysart, M. S. Patterson, "Photobleaching kinetics, photoproduct formation, and dose estimation during ALA induced PpIX PDT of MLL cells under well oxygenated and hypoxic conditions," Photochemical & Photobiological Sciences 5(1), 73-81 (2006).
11. A. Dimofte, J. C. Finlay, and T. C. Zhu, "A method for determination of the absorption and scattering properties interstitially in turbid media," Physics in Medicine & Biology 50(10), 2291-2311 (2005).
12. T. Johansson, M. Soto Thompson, M. Stenberg, C. A. Klinteberg, S. Andersson-Engels, S. Svanberg, and K. Svanberg, "Feasibility study of a system for combined light dosimetry and interstitial photodynamic treatment of massive tumors," Applied Optics 41(7), 1462 (2002).
13. L. Wang, S. L. Jacques, and L. Zheng, "MCML—Monte Carlo modeling of light transport in multi-layered tissues," Computer Methods and Programs in Biomedicine 47(2), 131-146 (1995).
14. N. Lopez, R. Mulet, and R. Rodríguez, "Tumor reactive ringlet oxygen approach for Monte Carlo modeling of photodynamic therapy dosimetry," Journal of Photochemistry and Photobiology B: Biology 160, 383-391 (2016).
15. B. Liu, T. J. Farrell, and M. S. Patterson, "Comparison of photodynamic therapy with different excitation wavelengths using a dynamic model of aminolevulinic acid-photodynamic therapy of human skin," Journal of Biomedical Optics 17(8), 0880011 (2012).
16. C. L. Campbell, K. Wood, C. T. A. Brown, and H. Moseley, "Monte Carlo modelling of photodynamic therapy treatments comparing clustered three dimensional tumour structures with homogeneous tissue structures," Physics in Medicine & Biology 61(13), 4840-4854 (2016).
17. J. Cassidy, V. Betz, and L. Lilge, "Treatment plan evaluation for interstitial photodynamic therapy in a mouse model by Monte Carlo simulation with FullMonte," Frontiers in Physics 3, (2015).
18. B. C. Wilson, M. S. Patterson, and L. Lilge, "Implicit and explicit dosimetry in photodynamic therapy: a New paradigm," Lasers in Medical Science 12(3), 182-199 (1997).
19. A. V. Khilov, D. A. Kurakina, I. V. Turchin, and M. Y. Kirillin, "Monitoring of chlorin-based photosensitiser localisation with dual-wavelength fluorescence imaging: numerical simulations," Quantum Electronics 49(1), 63-69 (2019).
20. A. V. Khilov, M. Y. Kirillin, D. A. Loginova, and I. V. Turchin, "Estimation of chlorin-based photosensitizer penetration depth prior to photodynamic therapy procedure with dual-wavelength fluorescence imaging," Laser Physics Letters 15(12), 126202 (2018).
21. A. V. Khilov, D. A. Loginova, E. A. Sergeeva, M. A. Shakhova, A. E. Meller, I. V. Turchin, and M. Yu. Kirillin, "Two-wavelength fluorescence monitoring and planning of photodynamic therapy," Sovremennye Tehnologii v Medicine 9(4), 96 (2017).
22. D. J. Robinson, H. S. De Bruijn, N. Van Der Veen, M. R. Stringer, S. B. Brown, and W. M. Star, "Fluorescence photobleaching of ALA-induced protoporphyrin IX during photodynamic therapy of normal hairless mouse skin: the effect of light dose and irradiance and the resulting biological effect," Photochemical & Photobiological Sciences 67(1), 140-149 (1998).
23. J. Moan, "Effect of bleaching of porphyrin sensitizers during photodynamic therapy," Cancer Letters 33(1), 45-53 (1986).
24. S. Anbil, I. Rizvi, J. P. Celli, N. Alagic, and T. Hasan, "A photobleaching-based PDT dose metric predicts PDT efficacy over certain BPD concentration ranges in a three-dimensional model of ovarian cancer," SPIE Proceedings 8568, 85680S (2013).
25. J. Tyrrell, S. Campbell, and A. Curnow, "Protoporphyrin IX photobleaching during the light irradiation phase of standard dermatological methyl-aminolevulinate photodynamic therapy," Photodiagnosis and Photodynamic Therapy 7(4), 232-238 (2010).
26. M. A. Scott, C. Hopper, A. Sahota, R. Springett, B. W. McIlroy, S. G. Bown, and A. J. MacRobert, " Fluorescence photodiagnostics and photobleaching studies of cancerous lesions using ratio imaging and spectroscopic techniques," Lasers in Medical Science 15(1), 63-72 (2000).
27. A. Johansson, F. Faber, G. Kniebühler, H. Stepp, R. Sroka, R. Egensperger, W. Beyer, and F. Kreth, "Protoporphyrin IX fluorescence and photobleaching during interstitial photodynamic therapy of malignant gliomas for early treatment prognosis," Lasers in Surgery and Medicine 45(4), 225-234 (2013).
28. M. T. Jarvi, M. S. Patterson, and B. C. Wilson, "Insights into photodynamic therapy dosimetry: simultaneous singlet oxygen luminescence and photosensitizer photobleaching measurements," Biophysical Journal 102(3), 661-671 (2012).
29. M. Kirillin, A. Khilov, D. Kurakina, A. Orlova, V. Perekatova, V. Shishkova, A. Malygina, A. Mironycheva, I. Shlivko, S. Gamayunov, I. Turchin, and E. Sergeeva, "Dual-wavelength fluorescence monitoring of photodynamic therapy: from analytical models to clinical studies," Cancers 13(22), 5807 (2021).
30. M. Kirillin, D. Kurakina, A. Khilov, A. Orlova, M. Shakhova, N. Orlinskaya, and E. Sergeeva, "Red and blue light in antitumor photodynamic therapy with chlorin-based photosensitizers: a comparative animal study assisted by optical imaging modalities," Biomedical Optics Express 12(2), 872 (2021).
31. M.-A. Mycek, B. W. Pogue (Eds.), Handbook of biomedical fluorescence, 1st ed., CRC Press, Boca Raton, FL (2003). ISBN: ISBN: 9780203912096.
32. J. Swartling, A. Pifferi, A. M. K. Enejder, and S. Andersson-Engels, "Accelerated Monte Carlo models to simulate fluorescence spectra from layered tissues," Journal of the Optical Society of America A 20(4), 714 (2003).
33. Y. H. Ong, J. C. Finlay, and T. C. Zhu, "Monte Carlo modeling of fluorescence in semi-infinite turbid media," SPIE Proceedings 10492, 104920T (2018).
34. Y. H. Ong, M. M. Kim, J. C. Finlay, A. Dimofte, S. Singhal, E. Glatstein, K. A. Cengel, and T. C. Zhu, "PDT dose dosimetry for Photofrin-mediated pleural photodynamic therapy (pPDT)," Physics in Medicine & Biology 63(1), 015031 (2017).
35. Y. H. Ong, T. C. Zhu, "Estimation of fluorescence probing depth dependence on the distance between source and detector using Monte Carlo modeling," SPIE Proceedings 11628, 116280X (2021).
36. D. Y. Churmakov, I. V. Meglinski, S. A. Piletsky, and D. A. Greenhalgh, "Analysis of skin tissues spatial fluorescence distribution by the Monte Carlo simulation," Journal of Physics D: Applied Physics 36(14), 1722 (2003).
37. A. Ustinova, I. Bratchenko, and D. Artemyev, "Monte Carlo simulation of skin multispectral autofluorescence," Journal of Biomedical Photonics & Engineering 5(2), 020306 (2019).
38. A. V. Khilov, E. A. Sergeeva, D. A. Kurakina, I. V. Turchin, and M. Yu. Kirillin, "Analytical model of fluorescence intensity for the estimation of fluorophore localisation in biotissue with dual-wavelength fluorescence imaging," Quantum Electronics 51(2), 95-103 (2021).
39. R. M. Valentine, C. T. A. Brown, H. Moseley, S. Ibbotson, and K. Wood, "Monte Carlo modeling of in vivo protoporphyrin IX fluorescence and singlet oxygen production during photodynamic therapy for patients presenting with superficial basal cell carcinomas," Journal of Biomedical Optics 16(4), 048002 (2011).
40. M. Y. Kirillin, D. A. Kurakina, V. V. Perekatova, A. G. Orlova, E. A. Sergeeva, A. V. Khilov, P. V. Subochev, I. V. Turchin, S. Mallidi, and T. Hasan, "Complementary bimodal approach to monitoring of photodynamic therapy with targeted nanoconstructs: numerical simulations," Quantum Electronics 49(1), 43-51 (2019).
41. E. Zenkevich, E. Sagun, V. Knyukshto, A. Shulga, A. Mironov, O. Efremova, R. Bonnett, S. P. Songca, and M. Kassem, "Photophysical and photochemical properties of potential porphyrin and chlorin photosensitizers for PDT," Journal of Photochemistry and Photobiology B: Biology 33(2), 171-180 (1996).
42. L. O. Svaasand, P. Wyss, M.-T. Wyss, Y. Tadir, B. J. Tromberg, and M. W. Berns, "Dosimetry model for photodynamic therapy with topically administered photosensitizers," Lasers in Surgery and Medicine 18(2), 139149 (1996).
43. L. O. Svaasand, B. J. Tromberg, P. Wyss, M.-T. Wyss-Desserich, Y. Tadir, and M. W. Berns, "Light and drug distribution with topically administered photosensitizers," Lasers in Medical Science, 11, 261-265 (1996).
44. D. A. Lintzeri, N. Karimian, U. Blume-Peytavi, and J. Kottner, "Epidermal thickness in healthy humans: a systematic review and meta-analysis," Journal of the European Academy of Dermatology and Venereology 36(8), 1191-1200 (2022).
45. I. M. Heerfordt, C. V. Nissen, T. Poulsen, P. A. Philipsen, and H. C. Wulf, "Thickness of actinic keratosis does not predict dysplasia severity or P53 expression," Scientific Reports 6(1), 33952 (2016)
46. D. Yudovsky, L. Pilon, "Rapid and accurate estimation of blood saturation, melanin content, and epidermis thickness from spectral diffuse reflectance," Applied Optics 49(10), 1707 (2010).
47. V. V. Tuchin, "Tissue optics and photonics: biological tissue structures," Journal of Biomedical Photonics & Engineering 1(1), 3-21 (2015).
48. A. R. Kamuhabwa, R. Roelandts, and P. A. De Witte, "Skin photosensitization with topical hypericin in hairless mice," Journal of Photochemistry and Photobiology B: Biology 53(1-3), 110-114 (1999).
49. M. Gallardo-Villagran, D. Y. Leger, B. Liagre, and B. Therrien, "Photosensitizers used in the photodynamic therapy of rheumatoid arthritis," International Journal of Molecular Sciences 20(13), 3339 (2019).