Medical imaging: Simulations of reflectance and fluorescence of human tissue

Methods

Optical Processes in Skin Tissue

When shining light into the skin, there are multiple light-matter interaction processes happening simultaneously, including reflection, absorption, scattering, and fluorescence. Reflection refers to the back-scatter of photon at the interface between different materials when the index of refraction changes. Since the main composition of human skin is water, the variation of index of refraction is minor therefore there is negligible amount of reflection inside human skin. Absorption is an intrinsic property of light-matter interaction where the amplitude of electromagnetic field attenuates as propagating. In human skin tissue, the total absorption can be decomposed into several different components, including absorption from water, fat, melanin, and blood content, which can be cast by an empirical equation:

$nu$

Monte Carlo Simulation

Due to the randomness nature of the scattering process and absorption process related to skin tissue, the optimal approach to simulate the skin behavior is to use Monte Carlo simuilation. Monte Carlo algorithm is a computationally algorithm that

Theoretical Model Extraction

If considering the activities in the epidermis layer, when illuminated by blue light pump beam, absorption will happen in both the tissue and the FAD, which is described by a single quatity, the absorptivity $\mu _{\alpha }$. In addition, FAD fluorescence is excited, and the fluorescence intensity is determined by the amount of blue photon absorbed by the FAD molecules. However, under regular material definition inside the simulation, it is impossible to separate the absorption by FAD from the absorption from the surrounding tissue, such that the actual fluorescence intensity will be inaccurate.

The multilayer model for fluorescence Monte Carlo simulation.

Here, the problem is solved by a multilayer model where the FAD molecules are collected into thin slabs distributed inside the medium. In this case, two different materials can be assigned and the absorption of FAD and surrounding tissue can be separated. However, by reshaping the geometry, the material absorpsivity will be correspondingly modified. Below, the equivalent absorptivity is derived to correct for this change.

The total optical absorption can be described as:

Absorption $\propto \int dV(-J)\cdot E=\int dV(-i\omega \varepsilon _{0}\chi |E|^{2})=\int dV\omega \varepsilon _{0}\chi ''|E|^{2}$

where $J$ is the current density, $E$ is the electric field, $\omega$ is the angular frequency of light, $\varepsilon _{0}$ is the free space electric permittivity, $\chi$ is the material susceptibility, and $\chi ''$ is the imaginary part of $\chi$ representing the optical absorption of a certain material. Under the representation of Beer-Lambert law, the electric field penetrating in the skin can be described as an exponential decay: $E=E_{0}\exp(-\mu _{\alpha }z)$, where $\mu _{\alpha }$ is the material absorpsivity in the unit of [cm]$^{-1}$. Substituting the exponentially decaying electric field into the absorption expression and applying our one-dimensional multilayer geometry:

Absorption $\cong A\int dz\exp(-2\mu _{\alpha }z)=A{\frac {1-\exp(-2\mu _{\alpha }L_{z})}{2\mu _{\alpha }}}$,

where $L_{z}$ is the total length in z direction. To have the same amount of absorption in the fluorophore when using the multilayer model, we need to equivalently adjust the optical property $\mu _{\alpha }$ of the material to balance the effect of shrinking the volume occupied by the fluorophore. Gathering the fluorophore into a thin slab centered at $L_{z}/2$, the absorption here is:

Absorption $\cong A\int _{L_{z}/2-\Delta {z}/2}^{L_{z}/2+\Delta {z}/2}dz\exp(-2\mu _{\alpha }'z)$

$=A{\frac {\exp(-\mu _{\alpha }'L_{z})(\exp(\mu _{\alpha }'\Delta {z})-\exp(-\mu _{\alpha }'\Delta {z}))}{2\mu _{\alpha }'}}\cong A\exp(-\mu _{\alpha }'L_{z})\Delta {z}$.

Equaling the absorption and assuming the absorpsivity of FAD is small ($\mu _{\alpha }L_{z}\ll 1$), one can extract an approximate expression for the equivalent absorpsivity:

$\mu _{\alpha }'={\frac {1}{L_{z}}}\log({\frac {L_{z}}{\Delta {z}}}).$

Using the equivalent absorpsivity correction, the total absorption is the same. In addition, the blue photon absorbed by FAD will result in the fluorescence emission, such that the correct FAD fluorescence activity can be simulated.

Results

The properties of the dysplastic region directly determines the fluorescence signal. Here, the properties refer to the geometry size, shape, the fluorophore concentration, and the location inside the skin tissue. Below, we will carry out a few computational and numerical experiments to qualitatively describe the dependence between the fluorescence signal strength and these factors. The simulation is conducted based on MCMatlab Fluorescence Monte Carlo simulation package.

Qualitative Dependence on the Dysplastic Width

(a) Schematic of the geometry, which is defined as a one-dimensional variation in FAD concentration with a low concentration at the central part representing the dysplastic region. (b) The actual fluorescence signal measured as a function of x position. The contrast signal is extracted as the ratio between the magnitude of the dip and the maximum signal level. (c) The dependence betweewn the width of the dysplastic region and the fluorescence contrast signal. (d) The measured fluorescence signal for each width value with an increasing width from bottom to the top row.

For the first study, the simulated geometry is defined as a one-dimensional concentration variation in "x" as shown in the figure above with a lower FAD concentration at the central region representing the dysplastic region. The width of the dysplastic region is swept and the output fluorescence contrast signal is measured. It is shown that the contrast signal first goes up, then saturates at w~0.35 mm. At small width values, the signals are mixed up from the two edges where the concentration changes sharply, so that one cannot resolve the full contrast signal. When the two edges are separated sufficiently far, the full signal can be resolved, giving a relative fluorescence contrast as 0.3.

Qualitative Dependence on the Dysplastic Depth

Next, the qualitative relation is studied between the depth of the dysplastic region inside the skin and the fluorescence contrast signal. In this case, the geometry is defined as a rectangular wire structure as depicted in the figure below. In the simulation, the width of the rectangular wire is fixed to 0.4 mm and the thickness to 0.1 mm. The depth is set as the swept parameter, and it is defined as the distance between the top surface of the epidermis layer and the top surface of the rectangular wire.

(a) Schematic of the geometry, where the FAD concentration inside the rectangular rod is lower, representing the dysplastic region. The depth of the rectangular wire is swept. (b) The dependence betweewn the depth of the dysplastic region and the fluorescence contrast signal. The red dashed line illustrates an exponential fitting of the simulation results. (c) The measured fluorescence signal for each depth value with an increasing depth from bottom to the top row.

From the simulation, it can be observed that the fluorescence contrast signal decreases as the depth increases, which agrees with intuition. As for a deeper location inside the skin, it takes farther distance for the fluorescent light to travel out, therefore the distortion introduced by the scattering events as the light propagates through the skin tissue is stronger. In addition, the attenuation from both scattering and absorption is stronger as well, which further reduces the SNR. For depth greater than 0.1 mm, the contrast signal strength is nearly vanishing.

Qualitative Dependence on the FAD Concentration Difference

Similarly, the strength of the contrast signal should also be a function of the concentration difference of FAD between the dysplastic and the healthy tissue. The geometry is again defined as a rectangular wire structure as depicted in the figure below. In the simulation, the width of the rectangular wire is fixed to 0.4 mm, the thickness to 0.1 mm, and the depth is set to 0.04 mm. The volumetric concentration of FAD in the healthy region is fixed to 0.2 while the FAD concentration in the dysplastic region is varied from 0 to 0.2.

(a,c) Schematic of the geometry, (a) has a shallower dysplastic region while it is deeper in (c). The FAD concentration inside the rectangular rod is swept. (b,d) The dependence betweewn the concentration difference and the fluorescence contrast signal. (c) The measured fluorescence signal for each concentration value with a decreasing concentration contrast from bottom to the top row.

From the simulation, it can be observed that for 0.04 mm as the depth, the fluorescence contrast signal is linearly proportional to the FAD concentration difference, which again agrees with intuition. The concentration difference determines the distribution in the number of fluorescence emitter, which is equivalent to the source and directly proportional to the fluorescence intensity. When increasing the depth of the rectangular wire to 0.1 mm, again the linear dependence is shown however with an overall lower amplitude, representing the decrease in the number of fluorescence emitter.

Appendix

You can write math equations as follows: $y=x^{2}+5$

You can include images as follows (you will need to upload the image first using the toolbox on the left bar, using the "Upload file" link).