Introduction and Context
Modeling of absorption, diffusion, and fluorescence pheno-mena in a biological tissue allows the study of the distribution of light that propagates in it, while taking into account its intrinsic optical properties.1, 2 These physico-optical characteristics depend chiefly on the concentrations, sizes, morphologies, spatial orientations, biochemical composition, and metabolic activity of the tissue compounds as well as their structural organization (cell and conjonctive tissue layers varying according to the type of organ, for instance).3 These elements are fundamental characteristics of the gold standard histological diagnosis of healthy and pathological states of biological tissues. Among the main light-absorbing tissues (chromophores), melanin, hemoglobin, and water are to be found.4 The elastic scattering of light in tissues is due to variations of the indices of refraction between intra- and extracellular elements (nuclear cells, cell organelles, etc).5 Finally, multiple fluorescing macromolecules also exist naturally in living tissues such as proteins (collagen, elastin, keratin), enzymes and coenzymes (various forms of flavins, NADH), or porphyrins.1, 6 These intrinsic fluorophores have relatively different but overlapping excitation and emission spectra cove-ring spectral bandwidths of and , respectively. These fundamental elements are targeted by optical methods for in-vivo diagnosis and biological tissue characterization.7, 8, 9
For a given tissue (especially epithelial ones), the presence and concentration of many fluorophores vary with the depth, and according to pathological changes, including cancer. Several studies show that a classification of the states of healthy, hyperplastic, and dysplastic tissues is made possible by using different wavelength bands, the shape of spectra, and the relative amplitudes of emission intensities (including the contrast between porphyrins and NADH/FAD, i.e., versus ),10, 11 or by in-depth access to the intrinsic fluorescence spectra (the fluorescence spectra acquired at the surface are modulated by the absorption from chromophores or fluorophores and the diffusion of the medium).12, 13, 14 For example, Biswal 15 show that the profile of autofluorescence spectra obtained with excitation is clearly altered by the characteristics of porphyrin absorption dips at 540 and . These studies also highlight that the measured bulk tissue autofluorescence results from the superposition of emissions of these molecules, whose absorption (excitation) and emission spectra overlap. For example, the absorption spectrum of elastin overlaps with the second peak of the NADH absorption spectrum and a portion of the one of collagen, as well as those of porphyrins and the first peak FAD , so that the emission spectra of fluorescent proteins (elastin, collagen, keratin, ) and NADH respectively overlap with the absorption spectrum of porphyrins and the second absorption peak of FAD.16 The consideration and analysis of these interactions come together to properly model the interactions of absorption/multiple excitation and spectral emission of several fluorophores.
Consequently, knowledge and/or identification of tissue optical properties is required in many clinical applications such as in-vivo diagnosis and noninvasive follow up, photodynamic therapy (PDT), or laser surgery.17 A large part of photodiagnostic methods target cancerology applications and especially studies of epithelial tissues (hollow organs, skin) where 85% of the cancers occur. The implementation of spectroscopic methods in situ can easily be performed by means of multiple-fiber probes (i.e., called “point spectroscopy” or “optical biopsy”).6, 18 The diagnostic efficiency of such approaches has already been validated separately for autofluorescence, diffuse reflectance, and elastic scattering spectroscopies on numerous tissues (skin,18, 19 cervix,20 oral mucosa,21, 22 bladder,23 and breast24). However, it can still be improved by combining those methods, as shown for the cervix,8, 25 breasts,26 head and neck,14 bronchial tree,27 and oral cavity.28
Models of light propagation in tissues are of importance and needed in a number of works where the extraction of local intratissular information based on global surface measurements is a challenge (for example, fluence spatial distribution, oxygen or chromophore concentration, intrinsic fluorescence quantification, and optical coefficient identification).2, 3, 9, 29, 30, 31, 32
The Monte Carlo method has been widely used in the past in many studies for simulating the propagation of light in turbid and complex media, such as biological tissues.33 It consists in propagating a packet of elementary energy (photon or group of photons) step by step, with random sampling of the length and direction of the displacement steps, as well as the proportion of energy absorbed in each elementary voxel of the medium, with reference to probability distribution functions depending on the optical coefficients of the medium.
The advantage of this statistical simulation method lies in first providing accurate results, whatever the interfiber or source-detector distances, and second in providing the detailed spatial distribution of the fluence while taking into account a large number of parameters: the intensity spectrum of the excitation source, the locations and geometries of both illumination and light collection, the local coefficients of absorption and diffusion, the refractive indices and dioptres, etc.34
Wilson and Adam35 and Flock 36 were among the first to introduce the use of the Monte Carlo method for modeling light-tissue interactions and propagation in multilayer structures. They validated their algorithm on highly diffusing phantoms at a unique wavelength of . In 1993, Wu, Feld, and Rava37 proposed a monolayer theoretical model coupling absorption, diffusion, and fluorescence, and taking the quantum yield as well as an isotropic emission of fluorescence into account. An experimental validation performed on arterial tissues highlighted that intrinsic fluorescence information may be extracted from global colocalized measurements of fluorescence and diffuse reflectance.
Between 1992 and 1995, Wang, Jacques, and Zheng34, 38 developed Monte-Carlo for multilayered media (MCML), a reference program in C-code for simulation in multilayer models, freely available on the web, and widely used since then for simulating light absorption and light diffusion between 350 and , in various experimental configurations39, 40 and on different types of tissues.3, 5, 33, 41
Monte Carlo methods are widely used to model the fluorescence emitted by the tissues to study the effects of excitation/emission geometries on fluorescence, as well as those of the absorption and diffusion properties of the medium, and to study the links between bulk measured fluo-rescence and intrinsic fluorescence at different depths of tissue.42 They have been applied by Keijzer 43 and Welch 44 to study the human aorta autofluorescence, depending in particular on excitation and emission geometries and on the fluorescence quantum yield; by Zonios 45 to connect the spectra measured in clinics and the human colon tissue microstructure; by Qu 12 to measure the optical properties to study the distortion spectra of autofluorescence depending on the optical properties of lung tissues; and by Drezek 13 to determine the relative contribution of two fluorophores present in the epithelial (NADH) and stromal (collagen) layers of a bilayer tissue model. In 1997, Zeng 2 simulated fluorescence excitation and emission in two steps: Monte Carlo calculation of the fluence distribution at one single excitation wavelength of , then fluorophore excitation with reference to that distribution and final calculation of the propagation for fluorescence emission (seven-layer model, bandwidth). Following the same approach, Welch in 199744 applied Monte Carlo to a multilayer model for simulating the propagation of an excitation light at and the emission of fluorescence light at to study the factors that influence the shape of the autofluo-rescence spectra measured. They showed that 90% of the emitted fluorescence measured comes from thicknesses of tissue up to (at ). This algorithm was then modified by Drezek in 2001, and applied to cervical tissues (monoexcitation at ).13 In 2001, Müller 29 validated a correction method of the absorption and diffusion contributions to recover an intrinsic fluorescence spectrum, taking into account the spectral dimension of the optical coefficients. However, the algorithm proposed was not validated for fluorophores with high quantum efficiency and high absorption, and secondary emission of fluorescence was not implemented. In 2003, Liu, Zhu, and Ramanujam33 were the first to describe precisely an algorithm of statistical simulation of fluorescence with one excitation wavelength and one emission wavelength. They also gave the expressions of the optical coefficients of the model with reference to absorption and diffusion by fluorescence. A rigorous experimental validation was carried out with fluorescence and diffuse reflectance on phantoms for an excitation wavelength of and an emission one at . Ma 32 compared the results obtained with a simplified analytical model of fluorescence and a simplified method of Monte Carlo for simulating time-resolved fluorescence. They observed that the Monte Carlo method gives better results when fluorophores lie in areas close to the source and to the detection fibers, as well as for short interfiber distances. A Monte Carlo code was developed by Vishwanath and Mycek46 in 2005 to simulate excitation and time-resolved fluorescent light propagation at multiple wavelengths in a bilayered turbid model. Experimental validation was conducted on two-phase (liquid and solid) layered tissue-simulating phantoms. In 2006, Chandra 47 obtained good agreement in comparing reflectance and fluorescence spectroscopic measurements made on cartilage tissue samples with simulation results using their time-resolved photon propagation code for multifluorophore emission.
The Monte Carlo technique developed by Churmakov 48, 49 simulates the spatial distribution of fluorescence excitation by combining a scattering-center matrix (calculating scattering paths) and an absorption calculation based on the path length applied to a multilayer model of human skin. In 2008, Palmer and Ramanujam50 experimentally validated a new extraction method of intrinsic fluorescence with a Monte Carlo model based on the work of Swartling 51 and Graff 52 The model is validated using liquid synthetic phantoms with a single fluorophore in light-induced fluorescence spectroscopy ( , ) and diffuse reflectance spectroscopy . The proposed model optimizes the fast-track approach based on the reciprocity principle to determine the absorption energy density grid, applying scaling relations to the two grids according to the set of optical parameters and using quasidiscrete Hankel transform for the convolution of the two grids (faster). Liebert 53 also developed a new Monte Carlo algorithm that simulates photon fluorescence generation and diffusion in a multilayer time-resolved model. The principle of absorption and scatte-ring of a packet of photons computes the probability of absorption fluorescence as a conversion (at an excitation wavelength) of a specific percentage of this packet; the emitted photons are then propagated at the emission wavelength. This methodology implies an isotrope distribution, values of reduced scattering coefficients identical at both emission and excitation wavelengths, and does not take the reabsorption into account.
Thus, one of the current difficulties resides in the coupled simulation of absorption, diffusion, and multiple fluorescence in a multilayer medium, over a broad spectral band. Furthermore, these works also confirm the necessity for a rigorous experimental validation of the simulation algorithms prior to their application to biological tissues.
Several of the algorithms that have been developed mostly consider multilayer models with a maximum of one fluorophore per layer and without taking into account any effect of reabsorption or re-emission. The work presented in this work is an incremental contribution to the continuation of the works of Liu, Zhu, and Ramanujam33 and Churmakov 48, 49 in particular by providing a solution taking into account both the spectral dimension (full spectra of absorption and emission of the fluorophores) and reabsorption/re-emission by several fluorophores whose absorption and emission spectra overlap.
Consequently, this work presents the principle of an algorithm developed for simulating the spectral absorption and emission of several fluorophores in a multilayer model. The experimental validation of the program implemented in C is performed with reflectance data measured at different distances between excitation and collection fibers on various mono- and multilayer phantoms containing absorbing, diffusing, and fluorescing compounds.
Materials and Methods
Monte Carlo Simulations
Basic principle and fundamental algorithm
The simulation of spatially resolved light transport in multilayered models is already well documented in the literature.34, 54, 55 We used these bases to develop our algorithm for multiple fluorescence.
Our modeling configuration of a semi-infinite medium with layers parallel to each other is considered, with one excitation optical fiber and several emission or collecting fibers (Fig. 1 ). Each layer is described by the following para-meters: the thickness (mm), the absorption coefficient , the scattering coefficient , the refractive index , and the anisotropy factor .
In short, the principle of random migration of photons consists in displacing a large number of photons or photon packets (usually a few hundreds of thousands at least) in directions and over distances drawn at random with reference to density probability functions characterizing the medium .
First, the source photons with an initial “weight” equal to unity are injected at origin coordinates from a random position inside the diameter of the excitation fiber. Each new 3-D direction of propagation is determined by statistical sampling of:
• the deflection angle in , i.e., the angle between the previous and new direction lines in the same plane, is defined by its cosine based on the Henyey-Greenstein scattering phase function56 as:
with being a random variable uniformly distributed between 0 and 1.
• the azimuthal angle (angle with the previous plane) in is defined as:
Algorithms previously developed for fluorescence simulation
Liu, Zhu, and Ramanujam33 experimentally verified Monte Carlo modeling of fluorescence and diffuse reflectance measurements in turbid, tissue-like phantom models. By consi-dering the optical parameters of a fluorophore (absorption coefficient of the fluorophore) and (fluorescence quantum yield), the authors defined the absorption probability for fluorescence as follows:and for measurements at a unique emission wavelength of .
Churmakov 48 simulated the spatial fluorescence distribution within the human skin. In their simulation algorithm, each photon packet produces only one fluorescent event. Their simulation results suggest that the distribution of autofluo-rescence is significantly suppressed in the near-infrared spectral region. In 2004, Chang 57 compared the results of Monte Carlo simulations with an analytical model of fluorescence in a two-layer model applied to normal and preneoplastic epithelial tissue. Their results show that the analytical model provides a good description of fluorescence for the optical properties of normal and precancerous cervical epithelium and stroma. Finally, Matuszak, Sawow, and Wasilewska-Radwanska58 implemented the simulation of fluorescence spectra in two steps: the simulation of the excitation spectrum to obtain the intensity of light inside cells containing a fluorophore, and then the simulation of emitted photons (emission spectra).
Simulation algorithm developed for multiple fluorescence (proposed solution)
Each layer of the medium (Fig. 1) is considered homogeneous with a thickness , a refractive index , an anisotropy factor , and absorption and scattering coefficients, respectively and . In our model, each layer can also contain one or more fluorophores expressed by the following complementary parameters: the absorption coefficient , the fluorescence quantum yield , the emission spectra , and the molecular concentration of the fluorophores.
The implementation of multiple fluorescence involves three principal steps: first the absorption of a photon by a fluorophore as a function of its absorption spectra and concentration in the medium, second the choice of an emission wavelength with reference to its emission spectra, third the generation of the fluorescence photon in a new direction. These various stages are detailed in the flowchart of Fig. 2, which is part of the main classical Monte Carlo algorithm.
Considering an absorbing and diffusing medium containing fluorophores with absorption coefficients and associated quantum yields , the stepsize is calculated in the following way:is a random number generated from a uniform probability distribution.
At each absorption event, a first test is carried out to decide whether it is statistically due to a fluorescent compound (Fig. 2). The probabilities that a photon is absorbed by a “pure” absorber and that a photon is absorbed by a fluorophore can respectively be defined as:: , or by the fluorophore if:
In the first case (photon absorbed by a “pure” absorber), the photon’s weight becomes:being the photon weights before and after absor-ption.
In the second case (photon absorbed by a fluorophore ), a second test is performed to decide whether conversion into fluorescence happens, depending on the fluorophore quantum yield (Fig. 2). Again, the principle consists in generating a random number uniformly distributed between 0 and 1.
• If this random number is higher than the quantum yield value , then the incident photon is not converted into fluorescence. It continues its way until the next event, after having updated its energy level with reference to as follows:
• If , then the incident photon is fully absorbed and a new photon of fluorescence is generated at the same coordinates in the medium. The first propagation direction of the new emission photon is randomly sampled between 0 and .
The emission wavelength of the new photon is also randomly sampled according to a probability distribution corres-ponding to the emission spectrum of the fluorophore, whatever the absorption wavelength. For example, if a fluorophore has an emission peak at a given wavelength, then the highest the probability emission of a fluorescence photon would be at this wavelength.
The values of the emission spectrum for each fluorophore being known, the following probability of emission of a fluorescence photon at a wavelength is associated:and accounts for the total number of points (wavelengths) used to discretize the emission spectrum.
A photon of fluorescence can cause another event of absorption fluorescence (the “cascade” effect). In addition, an incident photon absorbed at a wavelength cannot be converted into an emission photon with a wavelength lower than .
In the program based on the flowchart in Fig. 2, a specific file was created for each fluorophore in each layer, in which the values of the absorption and emission spectra, quantum yield, and absorption coefficient of fluorophores are stored. The simulation is made in the ascending way, i.e., from the lowest wavelength to the highest one, with sampling steps.
Finally, the function of cumulative distribution such asis used to choose the wavelength of emission in the following way: as a uniformly distributed random number.
Instrumentation Setup for Experimental Validation
We used a spatially resolved reflectance spectroscopy configuration. It consists in a laser diode (Laser 2000, France) at for the fluorescence excitation of all phantoms except phantom (see Sec. 2.3) for which a deuterium-tungsten halogen light source (DH2000, Ocean Optics, Dunedin, Florida) bandpass filtered at was used. Backscattered fluorescence emission spectra were measured using an imaging spectrograph (iHR320, HoribaJobin Yvon, France) with a back-illuminated charge-coupled device (CCD) detector, and a diffraction grating ( , blaze wavelength ). A multiple fiber optic probe (37 fibers with core diameter and 0.22 numerical aperture) coupled to a specific in-line fiber bundle at the entrance of the spectrograph allows multitrack operations, i.e., the simultaneous acquisition of up to 13 spectra correspon-ding to 13 different distances between illumination and collection fibers . Results shown in the present work are given for only one distance between illumination and acquisition fibers , but the results obtained were verified and validated for all other distances.
To improve signal-to-noise ratio, three spectra were automatically acquired for averaging. All experimental spectra were preprocessed to be free of any spectral distortions induced by instrumentation (spectral correction in wavelength and intensity). To improve reproducibility, a tripod was used to maintain the tip of the fiber optic probe in gentle contact with the surface of the phantom and perpendicular to it.
Phantoms were prepared with specific fluorophores, absor-ption, and scattering properties and dimensions, all of which are detailed hereafter. Liquid phantoms were made of absorber, scatterer, and one or two fluorophores, ethanol being used as thinner for the liquid phantoms since emission and absorption spectra of the fluorophores used here are given with ethanol as a thinner.38 The relevance of the monofluorophore phantoms was to check the good correlation between measurements and simulation results. The relevance of multifluorophore phantoms was to check if the cascade effect was correctly taken into account. Then, multilayer solid phantoms were prepared in two and three layers with one or two fluorophores.
India ink (Super Black India ink, Lefranc/Bourgeois, France) was selected as the absorber. The absorption coefficient of pure India ink is available in Ref. 33, and the absor-ption coefficient of India ink diluted with water (Kabi-Fresenius) at 0.011% was measured using a spectrophotometer (Lambda EZ210, in DCPR Department of Physics Chemistry of the Reactions, Nancy, France).
Intralipids (Sigma Aldrich, Saint Louis, Missouri) were chosen as scatterer. They are often used as scatterers because they do not have a strong absorption band in the visible. They are made of intravenously administered nutrients consisting of an emulsion of phospholipid micelles and water, easily avai-lable and relatively inexpensive.36, 59, 60 The optical coefficients of intralipids are available in Ref. 36. Figure 3 represents the absorbance and the scattering coefficient as a function of wavelength for pure India ink and intralipids, respectively.
Three fluorophores were selected for our experiments with reference to their distinct or partially overlapping emission spectra: eosin Y ( , Sigma Aldrich), fluore-scein ( , Fluka) and cryptocyanine ( , Sigma Aldrich). Indeed, eosin Y dissolved in ethanol has a quantum yield of 0.67 at , fluorescein a quantum yield of 0.79 at , and cryptocyanine a quantum yield of 0.007 at .38 Moreover, fluorescein is referred to several articles for the creation and validation of equipment.61, 62
For all the phantoms, the fluorophores were systematically weighed to obtain a quantity of of product. All simulations were performed with 1,000,000 photons launched at 40 different wavelengths between 400 and with steps (i.e., 40 million photons per spectrum), to ensure convergence of the results.
Four liquid phantoms were prepared: ( of fluorescein, of ethanol), ( of eosin, of ethanol), ( of cryptocyanine, of ethanol), and ( of fluorescein, of cryptocyanine in of ethanol as thinner). All solution phantoms were poured into a black container comparable to one expanded polystyrene container, until reaching a height of . The phantom mixing two fluorophores was prepared to test the cascading effect itself (or secondary fluorescence emission due to another fluorescence photon).
Solid mono- and multilayer phantoms
We used a solution of agarose (Agarose L, VWR Prolabo, France) at 5% ( of agarose powder and of distilled water) that has the property to become liquid at and solid around . Contribution of light scattering and fluorescence from agarose was checked to be negligible compared to the contributions of fluorophores at the wavelengths of excitation chosen for our experiments. Mono- and multilayer thicknesses were measured on sample sections with a Vernier micrometer (error: ).
A first set of bilayered phantoms with one fluorophore ( , , and ) was prepared, having at the bottom (second layer, thickness ) of fluorescein and of agarose powder in distilled water, and at the top (first layer, thickness ) either of India ink diluted at 0.011% , or of intralipids at 10% , or a mixture of both . For this set of phantoms, the first layer was kept in solution to prevent the optical properties from being modified by heating. A complementary single-layer phantom was also made with fluorescein, of agarose powder, and of distilled water (thickness ) to compare the spectral res-ponses obtained for the previous bilayered phantoms with the ones of the second layer alone.
A second set of two- or three-layered phantoms was made with two different fluorophores. One bilayer phantom was composed of of eosin, agarose, and of water on the top layer (thickness ), and of of fluorescein, of agarose powder, and of distilled water at the bottom layer (thickness ). The trilayered phantom was composed of a top layer of (thickness: ) India ink diluted at 0.011%; a middle layer of (thickness ) of eosin, of aga-rose powder, and of distilled water; and a bottom layer of (thickness ) of fluorescein, of agarose powder, and of distilled water.
Table 1 recapitulates the characteristics and composition of the six types of multilayer phantoms tested.
Characteristics of multilayer phantoms.
|First layer||India ink diluted at 0.011%||intralipids at 10%||solution 1||eosin|
|solution 2||agarose powder||solution 1|
|(Solution 1)||(Solution 2)||(Solution 4)|
|distilled water||—||Solution 3||Solution 4|
|Third layer||—||—||—||Solution 3||—||Solution 3|
To quantify the difference between experimental and corresponding simulated spectra, we calculated, for all phantoms, the mean of normalized differences as follows:representing the number of curve points, the experimental reflectance spectrum, and the simulated reflectance spectrum.
Results and Discussion
Primary Validation by Comparison to Simulation Results from Literature
The proper operation of our program, including our fluo-rescence algorithm, was first verified for absorption and diffuse reflectance without fluorophores. Therefore, our results were compared to other Monte Carlo simulation results referenced in the literature with identical model parameter-ization.34, 54, 63
The total diffuse reflectance of a monolayer slab of a turbid medium was computed with the following optical properties: absorption coefficient , scattering coefficient , anisotropy factor , relative refractive index (i.e., refractive-index-matched boundary), and thickness . Ten Monte Carlo simulations with 50,000 launched photons per wavelength were completed. The average and standard errors of the total diffuse reflectance are summarized in Table 2 together with the results from Van de Hulst,63 Prahl,54 and Wang, Jacques, and Zheng.34 All results agree with each other.
Comparison of published simulation results of average and standard errors of total diffuse reflectance Rd obtained for a slab of thickness e=0.02cm with a matched boundary ( μa=10cm−1 , μs=90cm−1 , g=0.75 , and n=1 ), and for a semi-infinite medium with a mismatched boundary ( μa=10cm−1 , μs=90cm−1 , g=0 , and n=1.5 ).
|Model||Source||Rd Average||Rd Error|
|Slab with a matched boundary||Ref. 63||0.09739||—|
|Semi-infinite medium with a mismatched boundary||Ref. 64||0.26000||—|
In the same way, we computed the average and standard errors of the total diffuse reflectance for a semi-infinite turbid medium that has a mismatched refractive index with the ambient medium and optical properties: , , (isotropic scattering), and . Ten Monte Carlo simulations with 5000 photons per wavelength were completed. Here as well, our results are in good correlation with the ones from Giovanelly,64 Prahl,54 and Wang, Jacques, and Zheng,34 as shown in Table 2.
As a primary validation, those results confirm that our program with fluorescence algorithm added operates properly in the absence of fluorophore in the simulation model.
Validation with Phantoms in Solution
All experiments were performed under monochromatic excitation with the laser diode source at , except the phantom performed with a a deuterium-tungsten halogen light source bandpass filtered at .
The correlation between simulated and measured emission spectra for each fluorophore was first verified using mono-fluorescence phantoms in solution (phantoms , , and ). Experimental and simulation curves are shown in Fig. 5 . The spectra were normalized with reference to the maximum peak value of the spectra. Experimental and simulated spectra were normalized in the same way. A relatively good correlation was obtained, with error values of 2.2, 1.9, and 2.5% respectively for phantoms , , and .
To test the proper operation of secondary fluorescence absorption (i.e., cascading fluorescence emission due to the absorption of another “already emitted” fluorescence photon), results in single excitation (phantoms , ) and multiple (or cascading) excitation (phantom : fluorescein, cryptocyanine, ethanol) were compared (Fig. 6 ).
We observe that there is good correlation between experimental and simulated fluorescence spectra given in Fig. 6 with a mean error of 5.3%. For the phantom and for the shortest distance between illumination and acquisition Fibers , standard deviations for the experimental and simulated spectra (respectively , ) are:, standard deviations are:
The emission intensity spectra of fluorescein between 480 and (peak at ) and of cryptocyanine at wavelengths above (emission peak at ) can be clearly noticed, and their proportional contributions are well reproduced by the simulation algorithm including cascading fluorescence. With reference to the absorption spectrum of cryptocyanine, we can even check in Fig. 6 that the amplitude of the fluorescein emission spectrum is clearly decreased in the wavelength band compared to the one of phantom in Fig. 5, indicating that absorption by cryptocyanine occurs and therefore that the cascade effect does exist.
Validation with Multilayered Phantoms
Experimental and simulated fluorescence spectra for multilayered phantoms (India ink/fluorescein), (intralipids/fluorescein), (India ink/intralipids/fluorescein) and (fluorescein) are given in Fig. 7 . The latter also highlights the way in which the fluorescence spectra are modified by the presence of absorbing and diffusing layers. Experimental spectra were normalized with reference to the maximum peak value of the experimental value of phantom , and simulated spectra with reference to the maximum peak value of the simulated value of phantom .
Looking comparatively at the spectra of and phantoms, it can be noticed that the addition of a purely absorbing top layer (India ink) leads to a global amplitude attenuation, while the overall shape of the spectrum remains similar. By comparing and spectra, we may observe that the presence of a diffusing top layer (intralipids) causes strong attenuation of the collected signal and also a change in the shape of the curve. In the case of phantom , the absorbing and diffusing top layer leads to attenuation and shape modification of the reflected spectra.
For these four different configurations, the results obtained experimentally and by simulation agree well. The lowest mean percent error obtained between experimental and simulated spectra is 2.5% for phantom (fluorescein), 6.0% for (intralipids, fluorescein), 7.1% for (India ink, fluorescein), and the maximum mean percent error is 9.4% for (India ink, intralipids, fluorescein). It can be noticed that the mean percent error is slightly higher for the solid phantom (2.5%) than for the phantom in solution (2.2%).
Figure 8 represents the measured and simulated fluore-scence spectra for the bilayers phantom with two fluorophores (eosin on the top layer and fluorescein at the bottom layer) and for the trilayer phantom with the same two fluorophore layers and an upper absorbing layer (India ink/eosin/fluorescein). A few explanations can be given as to the shape changes of spectra in Fig. 8. As presented in Fig. 3a, India ink absorbs twice as much at as at , therefore altering the spectrum unequally over the wavelength band. Thus, the eosin fluorescence appears proportionally with more intensity than the fluorescein fluorescence. Consequently, the modification of the shape in Fig. 8 is due to the absorption spectra of India ink and to the fact that this was added in a supplementary top layer. The mean percent error calculated between experimental and simulated spectra is 9.3% for the bilayer phantom with two fluorophores (eosin, fluorescein) and 10.6% for the trilayer phantom (India ink, eosin, fluorescein). It should be noted that the mean percent error increases with the addition of an absorbing layer, and with the number of layers (increased thickness). This trend is confirmed for the 12 other distances between illumination and acquisition fibers . For the interfiber distance, the average error varies from 3.1% for to 10.9% for . We can observe that the mean percent error values increase with the interfiber distances and with the medium complexity, i.e., the number of collected photons. Given the number of absorbing, fluore-scing, or diffusing layers, the increases in average error are probably due to the relative increase of statistical noise inherent to Monte Carlo simulation, because we used the same number of photons for all simulations (40 million photons per spectrum) whatever the phantom complexity.
Execution times vary with the composition of the phantoms, and more particularly with the optical parameters used in the simulation. Using a pentium 4, and 512-Mo RAM computer, execution times (20 million photons per spectra) range from for (fluorescein, agarose) to for three-layered phantom with two fluorophores (India ink, eosin, fluorescein). Simulation time is doubled between the monolayer monofluorophore phantom and bilayer bifluorophore .
Anatomical, physiological, and biochemical changes result in modifications of the absorption and light scattering in pathological tissues compared to healthy ones. These local changes of the optical properties affect the characteristics of reflectance diffuse spectra. The measurement of endogenous and exogenous fluorescence from tissue has become a very useful research and clinical diagnostic tool. Many authors have shown that fluorescence can indicate the presence of mali-gnant tissue, but one of the major limitations is that the fluorescence signal is mixed with absorption and scattering, which is intrinsic to the tissue components. As a consequence, Monte Carlo methods are widely used to model the fluore-scence emitted by the tissues. Part of the difficulty lies in the coupled simulation of absorption, diffusion, and multiple fluorescence in a multilayer medium, over a broad spectral band.
In this work, a Monte Carlo algorithm is developed for simulating the spectral absorption, diffusion, and multiple fluorescence in a multilayer medium, over a broad spectral band. As a primary validation, the correct operation of our program, including our fluorescence algorithm, is first verified in relation to other Monte Carlo simulation results referenced in literature for absorption and diffuse reflectance without fluorescence (same parametrization). The main experimental validation consists in comparing the fluorescence spectra measured on phantoms of various complexities to Monte Carlo simulated spectra of their corresponding models (mono- and multilayered). Four phantoms in solution with three different fluorophores and six mono- and multilayered solid phantoms with one or two fluorophores are used. The good correlation between simulated and experimentally measured results is quantified on the basis of mean error calculations. For phantoms in solution, the mean percent errors vary from 1.9% (one fluorophore) to 5.3% (two fluorophores). For multilayered phantoms, these errors range from 6.0% (two layers), to 9.4% for the trilayer phantom. The error values vary according to the complexity of the phantom, i.e., the model (for a fixed number of photons launched). Finally, the proper operation of “secondary fluorescence” (fluorescence emission due to the absorption of another already emitted fluorescent photon) is demonstrated.
This model proves appropriate for applications to biological tissues, in which the contribution of absorbers, scatterers, and fluorophores vary. More specifically, this work will be used for in-vivo tissue diagnosis such as precancer and cancer detection, where the identification of proper values of the optical parameters for these tissues takes into account the associated phenomenon of multiple fluorescence.
The authors would like to thank the Region Lorraine and the Ligue Contre le Cancer (CD 52, 54) for their financial support, and Céline Frochot (DCPR, Nancy) for complementary spectrophotometric measurements.