1. Introduction

Optical energy has been used for decades in dermatology both for diagnosis and for treatment, where the whole electromagnetic spectrum can be exploited including UV-, visible- and infrared radiation (near, mid and far) to address various skin conditions and diseases [1–5]. Application of visible and NIR light for therapeutic purposes covers a whole spectrum of cutaneous interventions across both cosmetic and medical domains. Examples include removal of vascular and pigmentary lesions [6], unwanted hair [7], and tattoos [8] to wound healing [9], scar resurfacing [10], skin rejuvenation [11, 12], and stimulation of hair growth [13–15], treatment of psoriasis and eczema [16,17] and more [18]. Effective therapeutic modalities, which rely on all five types of light-tissue interactions (plasma ablation, photodisruption, photoablation, photothermolysis and photochemical reactions), and where the impact of photons on tissue depends on the wavelength, optical power density and exposure time [19], have already been very successfully embraced by a range of professional and home-use devices [20,21].

However, with ongoing quest to develop new light-based diagnostics and therapies and to further enhance the parameters space of the existing ones towards even higher predictive power, efficacy and safety, it is essential to assess the amount of light reaching the multiple different targets inside the tissue.

Several tools exist to predict the propagation of light in tissue. They include the solutions of the radiative transfer equation under strict approximation, the analytical application of the random walk theory to biological tissue [22], and the Monte Carlo methods [23].

Monte Carlo methods have gained most popularity and confidence. They are highly appreciated when there is a need to model of light propagation in inhomogeneous tissues with complex geometry, such as human skin. The latter is associated with a complex compartmentalized structure, with optical properties that are rapidly varying along spatial dimensions. And so Monte Carlo methods have become an indispensable tool for dermatological applications.

Briefly, the method relies on the random sampling of the trajectories of photon packets propagating through medium. The calculations of the trajectories are based on the random occurrence of scattering and absorbing events. An essential requirement is thus the knowledge of the optical properties (absorption and scattering) associated with tissue. The skin is generally approximated as a pile of flat infinite layers such as the epidermis, the dermis, and the subcutaneous fat with homogeneous properties, though blood capillaries and other structures are sometimes added [24]. The successful functioning of the Monte Carlo method in delivering value for clinical applications in dermatology (e.g. when selecting optical treatment parameters for in vivo studies on humans based on extrapolation of settings shown effective using in vitro studies on cells) thus requires knowledge of the optical properties of the skin layers that most realistically represent the in vivo case.

The literature contains numerous articles with data on the optical properties of skin layers and its appendages. However, the associated published literature is often inconsistent reflecting a rather hap-hazard approach of by individual research investigators. Disparity between the quantitative values of the optical properties employed will create much variability in the quantitative predictions of photon density in tissue, even while keeping all other parameters fixed [25]. This explanation is straightforward as in the model the absorption of photon is, at the first order, proportional to the absorption coefficient of the tissue (skin is a turbid media with absorption much lower than scattering at visible and NIR wavelengths, μ_A << μ_S) [23,26]. Any significant variability in the absorption coefficient of any particular skin layer will largely affect the calculated photon density. Therefore, the wide variation in optical properties of skin layers reported in the literature (up to 100 fold [25]) not only poses a problem when selecting treatment parameters based on modeling but also raises a question about the validity of the cited references reporting the values.

Our proposed approach towards resolving existing parameter controversy in the literature was to get a better grasp on the origin of the spread in the skin layers optical properties, and then to select a narrower optical window, i.e., a set of the optical properties more closely resembling ’true’ biological values. This would result in stronger predicting power of the Monte Carlo model, eventually delivering more reliable outputs for clinical applications.

This proposed approach was guided by a simple yet very powerful observation [27] that the quantitative value of absorption and scattering of any tissue, and particularly of any homogeneous skin layer, are rooted in its biochemical composition and thus should reflect it. The molecular constituents of the epidermis, dermis and subcutaneous fat layer are known [27,28] and these have known absorption spectra with prominent peaks [29]. This knowledge should then serve us as a critical ’tool’, based on objective arguments, facilitating selection of more robust dataset(s) rather than embarking on a journey towards the quantification of the optical values de novo.

The authors are fully aware that the complexity about the thicknesses of the skin layers (epidermis, dermis and subcutaneous fat layer) and their variability between healthy individuals and under certain dermatological conditions; all of which impact the outcome of the Monte Carlo predictions. As an instructive example, in 2006 Gambichler et al. published a study showing the strong dependence of epidermal thickness on age, gender and anatomic site [30]. However, the thickness of the skin layers is measurable in vivo by non-invasive techniques such as optical coherence tomography (OCT), which have already proven their reliable measurement, 1.5 mm depth with 15 μm resolution [31]. By contrast, the measurements of optical properties of human skin layers extend in questionable large ranges [25]. Therefore, we elected to keep skin layers thicknesses fixed throughout calculations and to focus instead on the variability of the skin layers optical properties.

This approach was implemented, in the first instance, by evaluating each of the reported datasets in terms of the presence of absorption bands of the strongest chromophores over the visible to NIR range, specific for each of the skin layers: i.e. melanin - in the epidermis, hemoglobin - in the dermis, and lipids - in the subcutaneous fat. Subsequently, sample preparation and handling procedures during optical measurements were also investigated in an attempt to link the differences in the reported literature values with variations in experimental lab protocols. Moreover, the coherence or clustering of data between the datasets originating from different research investigations was analyzed.

Secondly, after the sets of the optical properties were narrowed down, Monte Carlo simulations of light propagation in the skin layers were performed and the impact of the remaining spread in the input parameters on the resulting photon densities was evaluated.

As a final verification step towards recommending specific sets of the optical properties, the diffuse reflectances estimated based on the simulations were compared with an independent dataset of in vivo measurement of the diffuse reflectance performed in human subjects.

2. Materials and methods

2.1. Monte Carlo optical model

An in-house Monte Carlo computer routine code of light propagation in tissue, based on the original algorithm published by Wang [23], was used to calculate the map of photon density inside human skin. The geometry of the skin sample was represented by three horizontal layers of defined thicknesses, with optical properties uniform within the layers. The thicknesses of each of the modeled skin layers were fixed in accordance to anatomic references: epidermal thickness [30] was set to 100 μm; dermal thickness [32] to 1.4 mm and subcutaneous fat layer [33] to 2 mm. No variation of the skin thicknesses was permitted.

The input beam was collimated, homogeneous and square with a diameter of 3000 μm. The irradiation beam was positioned centrally across the skin and perpendicular to its surface. The skin sample had an area of 2 by 2 cm² and thickness of 3.5 cm, while the grid resolution was 20 μm in X and Y (radial) axes and 16.5 μm in the Z (propagation) axis. The number of photon packets sent per simulation was 100e⁶.

2.2. Literature review of the optical properties of the skin layers

As the goal of our study was, in the first instance, to investigate variability in reported values of the optical properties of the human skin compartments throughout the visible to NIR spectral window, we reviewed existing literature data reporting on the optical properties of the epidermis, dermis and subcutaneous fat layer over 400–1000 nm range. For simplicity the model of the skin was limited to a three-layer skin structure.

Following that, the quantitative values of the optical properties of the skin layers, including the absorption and scattering coefficients (μ_A, μ_S), the anisotropy factor (g) and the refractive index (n), were extracted from 10 literature sources reporting optical properties of the skin layers published between 1981 to 2013 [6,24,27,34–40] which correspond to all the studies reporting optical properties of human skin layers to our knowledge.

A challenge we immediately encountered was that most of the datasets originating from one particular individual reference alone were not sufficient to perform Monte Carlo simulations over the visible to NIR spectrum. This was because not all the relevant optical properties (including μ_A, μ_S, g) of all the skin layers (epidermis, dermis, subcutaneous fat) were measured over the wide optical range. To solve this, we decided that the optical values from one or more references had to be combined to create the ’complete’ datasets.

Additionally, the values of some particular optical properties were reported in one publication only and/or its value was only present for one particular wavelength. In these cases, the value of that variable had to be fixed for all the constructed datasets of the optical properties and had to be assumed invariant with respect to light wavelength.

For example, the anisotropy factor of the subcutaneous fat layer and it refractive index were taken from Meglinski et al. [24] equal to 0.75 and to 1.44, respectively, for all the datasets of the optical properties that were constructed here, and was assumed to be independent of the wavelength. The index of refraction of both epidermis and dermis was extracted from Ding et al. [41].

Some references reported the optical properties of a more complex skin structure including for example the differences between blood-free and blood-containing dermis [24]. In such cases, when the structure was not a three-layer skin model, we chose to take the average of the values for each sub-compartment. Thus, for example, the average of the coefficients of the blood-free and blood-containing dermis was taken as the value representing the entire dermis.

3. Results

3.1. Optical properties of the skin layers

In general, approaches to quantify the optical properties of the biological tissue are categorized in two groups: mathematical models based on biochemical composition of a particular layer and experimental measurements (Tables, 1, 2 and 3).

Table 1. Factors accounted by the mathematical models for estimating the absorption coefficient of the skin layers in the bibliographic references included in the study.

View Table | View all tables in this article

Table 2. Factors explained by the mathematical models for estimating the scattering coefficient of the skin layers in the bibliographic references included in the study.

View Table | View all tables in this article

Table 3. Sample origin and methods used in the experimental measurements of the optical properties of the skin layers in the bibliographic references included in the study. N_sample represents the number of samples used for the measurement.

View Table | View all tables in this article

3.1.1. Mathematical models

In our literature search of the optical properties of the skin compartments we found 5 mathematical models [6,24,27,40,42] reporting analytical expressions to derive the sought after optical variable. Chronologically, the model of Svaasand [6] was first followed by Douven [42], Jacques [27], Meglinski [24] and Altshuler [40] each making a step forward, bringing further modifications and improvements.

All mathematical models considered here are built on a similar concept. In particular, quantification of the absorption/scattering coefficient of each of the skin layers is derived from the baseline absorption/scattering (i.e. the absorption/scattering of a skin layer without its major absorber(s)/scatterer(s)), as well as the contribution of the major absorber(s)/scatterer(s) of the layer. The latter ones, in their turn, can be measured or modeled separately (Tables 1 and 2).

We observed different degrees of complexity between the models. More specifically, the earlier models, such as the ones of Jacques [27] or Svaasand [6], have a lower number of factors taken into account as compared to the models of Altshuler [40] and Meglinski [24]. Indeed, neither the Jacques [27] nor Svaasand [6] models take into account the contribution of water to the absorption and scattering of the layers. In this current study, we omitted the model of Douven [42] simply because it represents an intermediate state between the two more extreme models, the one by Svaasand [6] and the one by Altshuler [40]. Thus, the optical properties obtained using the remaining four mathematical models were investigated.

Logically and as expected, there are natural links between the models as they represent a kind of ’evolutionary continuum’. In particular, the measurements of the baseline absorption and scattering of the skin proposed by Jacques [27] are re-used in the model of Meglinski [24], while the ones of Svaasand [6] are used in the model of Altshuler [40]. Some of the mathematical models also choose to have (partly) a quantification of the absorption and/or scattering coefficient, which is dependant on an experimental measurement such as Svaasand [6] or Altshuler [40], and where the variation with wavelength is approximated to a simple mathematical law. For example, the definition of the scattering coefficient (in epidermis and dermis) of Svaasand [6] is scaled on an empirical value of the scattering coefficient measured by Anderson [36] and Wan [35] at 577 nm and the variation with wavelength is assumed to be 1/λ law. Similar approaches are found in the model of Altshuler [40]. The similarity of the approach to quantify the optical properties of the skin layers between the models is directly reflected by the consistency of the resulting values of the optical properties (Fig. 1(dashed-lines)), where in particular, the absorption coefficients of the epidermal and dermal compartments are nearly agreeing between the studies selected here for analysis.

 

Fig. 1 Absorption and scattering coefficients versus wavelength from bibliographic sources for epidermis (A, D), dermis (B, E) and subcutaneous fat layer (C, F). Solid lines represent data extracted from the experimental measurements, dashed-lines - from the mathematical models.

Download Full Size | PDF

3.1.2. Experimental measurements

In total 6 unique studies [34–39] that report measurements of the absorption and scattering coefficients of the skin layers were found (Table 3). The skin type used in most studies was classified as Caucasian or fair skin, with the exception of the studies of Anderson [36] and Salomatina [37] where the skin type of the sample was not reported. The tissue samples were obtained either after surgery or post-mortem, where measurements were performed within a relatively large time frame, ranging from 1 hour post-surgery to 5 days after excision (Table 3). The body location of the samples was also highly variable and included face, abdomen, groin and breast. The handling of the skin before and during the measurement was also not consistent between the studies; the separation of layers was either performed via mechanical means e.g., using a razor or thermally. Moreover the preparation of the sample included various steps such as rinsing or not with saline or the addition of mechanical elements (such as microscopic slide, coverslip, plastic holders etc.) to stabilize the sample.

All these differences are expected to have not only an impact on the actual biological state of the sample such as hydration, swelling, amount of blood, degree of cellular necrosis, but also introduce uncertainties about surface roughness and geometry of a tissue slab, parameters of crucial importance for optical measurements. As a consequence, the results of the optical measurements (such as e.g., diffuse transmittance and reflectance and collimated transmittance) performed on these tissue samples, prepared under different conditions, will differ. Moreover, so too will the output of the back calculations strongly impact the estimated optical properties of the skin layers (due to the exponential power law between a bulk parameter, such as e.g. the diffuse transmittance and a microscopic parameter, such as e.g. the scattering coefficient).

What was consistent in all bibliographic references selected for study was the estimation of the optical properties based on the measurement of diffuse reflectance and transmittance of the sample. By contrast, the back-calculation methods, performed to quantify the absorption and scattering coefficients, differed from one reference study to the other i.e., the solution of the one dimensional diffusion approximation, the Kubelka-Munk model and the inverse Monte Carlo method.

These differences between the methods applied to estimate the optical properties of the skin layers, are reflected by the reported quantitative values (Fig. 1(solid-lines)), where there is little or, in a more strict sense, no quantitative agreement can be seen between the results (Fig. 1).

3.2. Origin of the spread in the optical properties

Our review of the reported optical properties of the skin compartments revealed a very large spread, both for the values of the absorption and of the scattering coefficients (Fig. 1). The range for the same variable, reported, by different research groups, extended to 100-fold.

3.2.1. Absorption coefficient of the epidermis

One can observe the existence of two clusters of reported absorption coefficients of the epidermis. One is formed by the values estimated using mathematical models of Jacques [27], Meglinski [24] and Altshuler [40] together with those resulting from the empirical measurement of Wan [35]. The other one, having in general a lower absorption, is populated by the values originating from the model of Svaasand [6] and the measurements by Marchesini [34] and Salomatina [37] (Fig. 1(A)).

In an attempt to apply our rational-based approach towards evaluating the optical properties of human skin, where we strongly rely on the fact that an absorption spectrum of any skin compartment is rooted into its biochemical composition, we consider first the epidermis layer.

The absorption of the human skin epidermis throughout the visible and NIR spectral range is strongly defined by that of melanin contained in melanosomes of the keratinocytes and basal melanocytes. The human Caucasian epidermis (light-pigmented) will contain on average 5 % of melanosomes [29], with large variabilities within individuals. The melanosome absorption spectrum were previously measured [49], and shown to have absorption at least higher than 100 cm⁻¹ over the whole visible to NIR range. However, the density of melanosome does not directly translate in density of melanin as many factors are involved such as skin type and sun exposure. The melanin density in human epidermis of Caucasian skin has also been estimated and is higher than 2.5 mg.mL⁻¹ for most skin types [50]. The extinction coefficient of melanin was also estimated and found higher than a few cm⁻¹.(mg.mL⁻¹)⁻¹ [51]. Using either the melanosome density and the melanosome absorption coefficient or the melanin concentration and the melanin extinction coeffcient, one could logically state that the value of the absorption coefficient of the epidermis should be at least equal to 5 cm⁻¹, where additional absorbers could increase this further.

When looking at the reported literature data, one can observe, however, that the epidermal absorption coefficients of the datasets of Salomatina [37], Svaasand [6] and Marchesini [34] are significantly close or lower (Fig. 1(A)) than this lower limit.

Additionally, the epidermis absorption coefficient given by the model of Svaasand [6] appeared quantitatively lower as compared to the other mathematical models, despite showing a consistent trend for the variation of the coefficient with wavelength (Fig. 1(A-dashed lines)).

An explanation of this shift towards a lower absorption level as identified. First of all, Svaasand [6] selected a single value 0.25 cm⁻¹ to model the background absorption of the epidermis (absorption without pigments present). It corresponds to the absorption of low-absorbing tissue (uterine tissue and eye) in the wavelength range 600–900 nm and it is assumed to be independent of wavelength (Table 1). In contrast to that assumption, Jacques [27] and Meglinski [24] relied on the wavelength-dependent absorption of ex vivo neonatal skin (which is nearly free of melanin and in optical sense is very much translucent, yet containing water and other chromophores inherently contained in human keratinocytes and interstitial fluid, e.g., flavoproteins), which showed a much higher number as compared to a single value used by Svaasand, in particular in the low wavelength range (an order of magnitude higher at 450 nm).

Secondly, the definition of the absorption of epidermis due to melanin is not consistent between inspected mathematical models. While Jacques [27] and Meglinski [24] share a similar mathematical formula, Svaasand [6] and Altshuler [40] methods were guided by an empirical definition. The definition of Svaasand [6] of the absorption of epidermis due to melanin yields quantitative values lower than any of the other models, i.e. almost 10 times lower than that of Altshuler [40] and 4 times lower than that of Jacques [27] at 450 nm.

When looking at the experimental sets (Wan [35], Salomatina [37] and Marchesini [34]) one can see that the data of Marchesini [34] stand out because of their low absorption level and the steep variation with wavelength (Fig. 1(A)). This steep decrease in the NIR range contrasts with a slowly decreasing absorption spectrum of melanin over the visible to NIR band [29,49].

The sample handling reported by all three studies (Wan [35], Salomatina [37] and Marchesini [34]) differ in one important aspect, the addition of saline to the sample before the measurement. This step was present in the methods of Salomatina [37] and Marchesini [34] (Table 3). It might be possible that the saline buffer has increased the volume of the epidermis sample and therefore reduced the absorption per unit of length as well as affecting tissue scattering and surface roughness, therefore reducing the global absorption level of the coefficients of Salomatina [37] and Marchesini [34].

Putting all together, it appears that only the mathematical models of Jacques [27], Meglinski [24] and Altshuler [40] and the experimental dataset of Wan [35] are coherent in terms of both basic spectroscopic arguments and order of magnitude that one can expect from the absorption of the epidermis.

3.2.2. Absorption coefficient of the dermis

The main absorbing components of the dermis are blood (strongly dominating across the visible spectrum) and water (which becomes prominent in the NIR region). Therefore, one should expect that the estimated absorption spectrum of human blood-containing dermis should show the characteristic absorption peaks of the hemoglobin, specifically at 420 nm and 540 nm in the short wavelength region of the visible spectrum [52], and that of water in the NIR area with a specific absorption band at 970 nm [46].

All the mathematical models assessed here quantitatively agree and show the expected absorption peaks in the visible and NIR part of the spectrum (Fig. 1(B-dashed-lines)). Also the absorption coefficient of Simpson et al. [38] derived from indirect optical measurements in the range 600–1000 nm shows the expected water peak at 970 nm and quantitatively agrees with the mathematical models (Fig. 1(B)). On the contrary, experimentally obtained datasets of Salomatina [37] and Anderson [36] do not reveal any of these expected spectral properties (Fig. 1(B)) and show higher global absorption level than most datasets. As far as sample preparation is concerned and was reported in literature, we hypothesize that rinsing the tissue with PBS and the rehydration with saline may have altered the content of the dermis and therefore the measured data of Salomatina [37] (Table 3).

3.2.3. Absorption coefficient of the subcutaneous fat layer

The main components of the subcutaneous fat layer are lipids and water [53]. The absorption spectrum of purified pig fat was reported to show a characteristic absorption peak at 930 nm [54]. No distinct clustering of the reported optical properties of subcutaneous fat layer was found in the analyzed datasets (Fig. 1(C)), except for a quantitative agreement of the model of Meglinski [24] and the experimental measurement of Simpson [38] in the spectral range 600–1000 nm (Fig. 1(C)).

The dataset of Simpson [38] is the only one showing the fat-related absorption peak at 930 nm (Fig. 1(C)). The model of Meglinski [24] does not show a specific fat-related absorption feature but otherwise is quantitatively in agreement with the measurements of Simpson [38]. In the design of their skin model, Meglinski et al. [24] surprisingly did not consider the fat as an absorbing component of the subcutaneous layer, and only considered blood and water (composition is 5 % of blood and 70 % of water and no lipids). This perhaps explains why in that dataset the fat-related absorption peak was not present. Given the molecular composition of the subcutis, one should genuinely expect a significant contribution of fat absorption [54]. Following the general equation for the absorption of a layer reported by Jacques [27], a possible revision of the model of Meglinski [24] when applied to the subcutaneous fat layer can be made by adding the purified fat absorption spectrum from van Veen [54].

3.2.4. Scattering coefficient of epidermis and subcutaneous fat layer

The values of the scattering coefficients of the epidermis and subcutaneous fat layer reported in the selected literature did not cluster for any skin layer (Fig. 1(D,F)). All studies tended to show a relatively weak agreement regarding the variation of their coefficient with wavelength. It might be argued that the values of the scattering coefficient of the subcutaneous fat as reported by Salomatina et al [37] are at least 2-fold higher than those of Bashkatov [39] and Simpson [38] (Fig. 1(F)). We can only speculate that similar to the case of epidermal absorption reported by Salomatina [37], sample preparation where PBS-based rinsing and a subsequent rehydration with saline may have played a role, resulting in much different values than those reported by Bashkatov [39] and Simpson [38].

3.2.5. Scattering coefficient of dermis

Finally, when looking at the reduced scattering coefficient of the dermis, one can distinguish two clusters based on the data from Altshuler/Svaasand and Jacques/Simpson/Salomatina (Fig. 1(E)).

Altschuler et al. [40] used a formula from Svaasand [6] to model the scattering of the blood. This was then translated to estimate the values of the scattering coefficient of the dermis (Fig. 1(E)).

A more simplistic approach was undertaken by Jacques et al. [27] where the definition of the scattering coefficient of the human dermis did not take into account the scattering of the blood. This could explain the apparent lower value of the scattering coefficient reported by Jacques [27] as compared to Altshuler [40] and Svaasand [6]. The experimental measurement of Salomatina [37] and Simpson [38] are in agreement with the dataset of Jacques [27]. This agreement might be related to a reduction of blood content of the human dermis during sample preparation (washing) for the optical measurement and thus be less representative for in vivo case.

3.3. Selection of the datasets of optical properties

As stated earlier, the objective of the study was to investigate the high variability of the reported values of the optical properties of the human skin compartments and to select subsets of the optical properties of the skin layers, presumably more closely representing realistic values using a rational-based approach. Guided by this strategy, only the most coherent datasets reflecting contribution of the major chromophores of the skin compartments were retained (Table 4). They were combined in 4 datasets containing the optical properties of each of the 3 main human skin layers (Table 5). We note that for the reason of limited available data, only one dataset of the subcutaneous fat layer was created combining the individual datasets of Simpson [38], Meglinski [24] and Bashkatov [39].

Table 4. Summary of the rational for the inclusion or exclusion of bibliographic references from analysis

View Table | View all tables in this article

Table 5. Selected datasets used in simulations together with the bibliographic references for the absorption and scattering coefficients (absorption / scattering) of epidermis, dermis and subcutaneous fat layer

View Table | View all tables in this article

3.4. Impact of the variability of skin layer optical properties on the predicted photon density distribution in depth and in section

After having reviewed available data on the optical properties of the skin and after having narrowed down a very wide parameter window to four sets of the values to be more coherent in resembling (or deviating least from) the in vivo optical properties, we applied Monte Carlo model to predict light propagation in the skin.

The resulting predictions of light propagation along the axial axis, i.e., in depth, now fell within a 6-fold range, in absolute terms. The greatest difference was observed with the ratio between the maximum photon density in datasets 1 (Jacques [27]) and 4 (Altshuler [40]) at 655 nm (Table 6).

Table 6. Ratios between the maximum of photon density reached in the skin with the datasets 1 (Jacques [27]), 2 (Wan [35]/Simpson [38]) and 3 (Meglinski [24]) over the maximum obtained with the dataset 4 (Altshuler [40])

View Table | View all tables in this article

Overall the ratios between the estimated photon densities obtained using the four selected datasets was stable over the 4 tested wavelengths. In particular, dataset 1 was systematically predicting the lowest maximal photon density, reaching around 25% of the magnitude obtained using the dataset 4. The values estimated using dataset 2 (Wan [35]/Simpson [38]) and the dataset 3 (Meglinski [24]) were reaching 50% and 75% of the maximum when using the dataset 4 (Table 6).

When looking in greater depth at the relative predictions of photon density level using the four datasets, a relatively good agreement was also observed. In particular, a maximum 2-fold difference was reached at a 2 mm depth between the photon density estimated using the datasets 1 and 4 at 450 nm (Fig. 2). Furthermore, estimated axial beam profiles expressed in terms of photon density in the upper dermis (700 μm below the skin surface) varied by maximum a factor of 5, which was observed again between the datasets 1 and 4 at 450 nm (Fig. 3). Within 2 mm distance from the irradiation point, the relative quantitative difference between the predictions using each of the datasets were agreeing with each other within a factor 2 only.

 

Fig. 2 Relative photon density versus depth obtained from the Monte Carlo predictions of optical transport in a three-layer human skin model using the selected optical properties datasets, shown in semi-logarithmic scale. The photon density presented was extracted from a rectangular cylinder of sizes 400 μm by 400 μm centered around the propagation axis.

Download Full Size | PDF

 

Fig. 3 Beam profile in the dermal layer (photon density) versus radial axis obtained from the Monte Carlo predictions of optical transport in a three-layer human skin model using the selected optical properties datasets, shown in semi-logarithmic scale. The photon density presented is the beam profile measured at 700 μm under the skin surface and averaged over 164 μm in the direction perpendicular to the propagation axis

Download Full Size | PDF

At 530 nm, and 850 nm, a similar trend towards a decrease of the photon density level in axial and lateral dimensions was observed using all datasets (Fig. 2 and 3), where there was a quantitative agreement between them to within a factor 2, at all spatial locations.

At 450 nm, however, dataset 4 stood out, where the photon density associated with it had a steeper decrease with depth in the epidermal layer, a slower decrease with depth in the dermal layer and a slower decrease with radius compared to the other three datasets (Fig. 2 and 3). As such, this resulted in a lower relative photon density in the epidermis and a higher photon density in the dermis (Fig. 2 and 3(450 nm)). This trend is due to a relatively higher epidermal absorption and scattering as well as due to a lower dermal absorption and scattering coefficients of the dataset of Altshuler [40] (Fig. 4). Similar difference were also observable at 530 nm (Fig. 2 and 3), however, to a lower extent as the absorption and scattering coefficients comprising the dataset 4 are quantitatively closer to the ones originating from the other three datasets (Fig. 4).

 

Fig. 4 Absorption versus scattering coefficients corresponding to datasets 1, 2, 3 and 4 for the epidermis (left) and dermis (right). The symbol indicates the dataset: 1 (○), 2 (+), 3 (x) and 4 (□). The colour of the dots indicates the wavelength, blue 450 nm, green 530 nm, red 655 nm and black 850 nm.

Download Full Size | PDF

3.5. Comparison of skin diffuse reflectance estimated using the Monte Carlo simulations of each selected optical properties dataset and an independent source of in vivo measurement

On our path to form a recommendation for a more robust source optical properties of human skin layers, the datasets reported in literature were consolidated to include only those believed to more closely approximate the in vivo case.

The final step towards forming such a recommendation was to compare the skin reflectance estimated using the Monte Carlo model and the selected optical properties to an independently obtained in vivo measurements on human skin. Human skin reflectance measured on 28 individuals with Caucasian skin type was published by Cooksey et al., data from which is available at the National Institute of Standards and Technology (NIST) [55]. The variability of the diffuse reflectance originating from individual human subjects as reported by the NIST reaches about 30% (Fig. 5(dotted-black lines)).

 

Fig. 5 Diffuse reflectance of human skin computed based on the Monte Carlo predictions of the propagation of light in a three-layer skin model using the selected optical properties datasets. Also shown is an empirical measurement of the diffuse reflectance of human Caucasian skin as measured by the NIST [55]

Download Full Size | PDF

Remarkably, we were positively surprised to observe that the diffuse reflectance estimated using a Monte Carlo model and all four selected datasets fell within the range of diffuse reflectance measured by the NIST, specifically for the visible wavelengths (Fig. 5(450, 530, 655 nm)).

One can easily recognize here two clusters in the calculated values of the diffuse reflectance, particularly at 655 nm: one cluster originating from the dataset of Altshuler [40] and Wan [35]/Simpson [38] and the other - based on that of Jacques [27] and Meglinski [24] (Fig. 5(655 nm)). These two clusters originate in the differences of absorption coefficients of the epidermis and dermis between these two pairs of datasets, where Altshuler [40] and Wan [35] showed higher epidermal absorption and lower dermal absorption in particular for wavelength longer than 600 nm (Fig. 1(A,B)).

The diffuse reflectance associated with the dataset of Meglinski [24] particularly is lower than that associated with the rest of the datasets (Fig. 5) at 450 nm and 530 nm. This differences can be explained by the combined effect of a higher dermal absorption, as well as a lower epidermal absorption as discussed earlier (Fig. 4).

4. Discussion

The objective of this study was to investigate the origin of the variability of the optical properties of human skin compartments reported in literature, to underpin the underlying rational behind this variability and to propose a much more robust dataset of values that better represent the in vivo case. Having achieved this, one would be equipped with a more realistic approach to predict light propagation through cutaneous compartments.

Our review of the literature revealed 4 mathematical models and 6 experimental measurement results reporting on the optical properties of one or several of the skin layers. Strikingly we found a very large spread between the reported values of optical properties of the skin layers.

More specifically, published values of the absorption and scattering coefficients of the epidermis, dermis and subcutaneous fat layer were found to extend over a very large range, where up to 100-fold differences were present for a selected coefficient of one of the skin layers (Fig. 1).

This existing stark variation poses a major problem for the application of the Monte Carlo optical model for the prediction of the light propagation in human skin. It also complicates the further tuning of treatment parameters and the interpretation of values for skin diagnostics, as the calculated photon density will be directly and strongly impacted by the choice of absorption and scattering coefficients for the skin components. Furthermore, such a large spread in the previously reported values raises a question about the validity of the studies reporting on the optical properties of the skin layers.

Therefore, as part of the overall aim of this work, we aimed to select subsets of optical properties of the skin layers, presumably more closely representing realistic values using a fact-based approach from already existing datasets. In an attempt to achieve this, a rational-guided method was adapted, relying on a simple yet powerful concept [27]. Light absorption and scattering of any tissue, and particularly any homogeneous skin layer, should simply be influenced by the biochemical composition of that structure.

This idea was implemented, in the first instance by evaluating each of the reported datasets for the presence of the absorption bands of the strongest chromophores over visible to NIR range, specific for each of the skin layers: namely melanin - in the epidermis, hemoglobin - in the dermis, and lipids - in the subcutaneous fat. Secondly, sample preparation and handling prior to and during the optical measurements were studied in an attempt to link differences in the reported values with variations in the associated experimental lab protocol. And finally, the coherence or clustering between the datasets originating from different research investigations was analyzed.

As a result of these approaches, it was possible to exclude from further considerations three sets of the optical properties estimated based on experimental measurements (Salomatina [37], Anderson [36] and Marchesini [34]) and one - based on a mathematical model (Svaasand [6]). All these datasets were discarded as they did not contain the spectroscopic features that one should expect to be present based on what we already know of the biochemical content of the skin layers. In particular, the mathematical model of Svaasand [6] and the experimental measurements of Salomatina [37], Anderson [36] and Marchesini [34] did not show the expected spectroscopic features in the quantitative values of their respective absorption coefficient (Table 4). Finally, as a result of such a screening, guided by our selective strategy, 6 published studies were retained (Wan [35], Simpson [38], Jacques [29], Meglinski [24], Altshuler [40] and Bashkatov [39]) as these reported what we believe to be the most coherent with each other, and were characterized by realistic optical properties reflecting contribution of the major chromophores of the skin compartments (Table 4).

When applying these datasets the variability of the absorption and scattering coefficients decreased from a factor 100 to under a factor 3 at fixed wavelength (Fig. 4). Looking at the sample preparation and handling, we propose that any experimental measurement of the optical properties of skin layers should try to avoid as much as possible altering the biochemical content of the skin layers. The importance and impact of the sample handling routine should not be underestimated. For example, the quantitative values of the optical properties measured by Salomatina [37] showed signs that the chemical composition of the layers was severely affected during measurement, potentially during the washing/rinsing phase. Critically, the absorption peaks of blood, water and fat were not visible in the corresponding layers (dermis and subcutaneous fat layer). Besides, the absorption coefficients of all layers (epidermis, dermis and subcutaneous fat layer) were quantitatively too similar to reflect their key anatomic distinctiveness: a slow decrease from around 10 cm⁻¹ at 400 nm to 0.5 cm⁻¹ at 1000 nm [28]. Measurements of the optical properties should ensure that the constituent pigments and other optically-active compounds are preserved in their original state.

Likewise, for the application of mathematical models to calculate optical properties of specific skin layers based on baseline properties and additional chromophores, we recommend the inclusion of all the contributions of the relevant pigments to the absorption of the layer. For example, the addition of fat to the absorption coefficient of the subcutaneous fat layer in the model of Meglinski might improve the validity of the model [24].

We suggest that one would benefit from preferentially using the values reported by Wan [35], Simpson [38], Jacques [29], Meglinski [24] and Altshuler [40] as their datasets are coherent with what can be expected based on the biochemical content of the skin layers.

Prior to performing Monte Carlo calculations of light propagation in the skin using the selected sources of literature, 4 complete datasets were constructed by merging individual datasets from 6 sources [24, 27, 35, 38–40] (Table. 5). This was necessary as individual publications did not report the complete set of the optical parameters over the total visible to NIR range for all three optical properties (absorption and scattering coefficients and anisotropy of scattering) for each of the three skin layers (epidermis, dermis, subcutaneous fat).

As the next step in our approach we calculated the light propagation in the skin using a Monte Carlo model and the selected sets of optical properties and evaluated the resulting differences in photon densities in terms of dependency on the wavelength, and axial- and radial dimensions. More specifically, the resulting impact of the variability of the optical properties on the output photon densities, and thus on the uncertainty in Monte Carlo prediction, was reported.

Logically, with the reduction in the variability of the optical properties of the skin layers at the input of the model, our simulations of light propagation revealed a much smaller quantitative disagreement between the photon densities calculated using different datasets, varying up to a factor of just 6 (Fig. 2).

We hypothesize that this variability could be reduced by further tuning the mathematical models and reducing the variability of the optical properties of the datasets. Indeed, most mathematical models have many factors that are only estimated and so might not be set to the most realistic values (melanin content of the epidermis, blood content of the dermis, list of chromophores of the layer, etc.). However, it will never be possible to find exactly defined single value for those due to the variability of the skin properties existing between individuals, body location, skin type, etc [30,56] and also the difficulty of measuring these quantities. For example, the melanosome density of the epidermis was shown to vary between 0.5 to more than 5 % in Caucasian skin adults [27]. In another example, a study on 83 human subjects found that epidermal thickness was measured with a standard deviation higher than 10% in adults. Additionally, hair follicle size can be highly variable in skin of different body sites within the same individual and between individuals [56].

Nevertheless, variability in the optical properties of the skin layers will always exist, and so this remaining variability between the selected datasets might be partly due to a natural variability existing in large group of individuals. We believe that the quantitative disagreement of the calculated photon density, due to the variability of the optical properties of the skin layers, might only be partly solved with a ’better’ measurement or model of the optical properties of the skin layers.

The last step towards narrowing down the parameter window was a comparison of the estimated values of human skin diffuse diffuse reflectance using our Monte Carlo model to an independent in vivo measurement [55]. While the variability of the diffuse reflectance originating from individual human subjects as reported by the NIST was about 30% (Fig. 5(dotted-black lines)), the diffuse reflectance calculated using Monte Carlo model (and associated with all of the four selected datasets) fell within the range reported by the NIST, specifically for the visible wavelengths (Fig. 5(450, 530, 655 nm)). This verification step supports the proposed selection of the datasets of the optical properties.

In addition to the variability of the optical properties, the predictions of the optical propagation in the skin by the Monte Carlo model are also impaired by an (over-)simplified representation of the complex geometry of the skin, deeply rooted in its physiology. The case of the absorption coefficient of the human epidermis is a perfect example. Melanin is the main pigment of the epidermis absorbing throughout the whole visible part of the spectrum. Presence, amount and compartmentalization of melanin has a complex fingerprint in the epidermis which depends on skin type and sun exposure [57, 58]. Melanin exists in two forms, eu- and pheo-melanins, with different optical properties [59]. Importantly, it is never present as totally ’freely-dispersed’ pigment in the skin, which conflicts with most representation of melanin assumed in the mathematical models. Originally it is produced by the melanocytes where it resides in melanosomes. Once it is transmitted to the keratinocytes, the melanosomes experience some degree of degradation by the keratinocyte phago-lysosomal machinery and form melanin granules [58]. Overall, about 70–80% of all melanin in the epidermis remains associated with the the basal layer of the epidermis, which will influence how light is absorbed locally in human epidermis i.e., likely leads to more heat generation in the basal layer. Most mathematical models reduce this complexity to an homogeneous distribution of a single pigment assumed to represent melanin in all its forms. This over-simplification will definitely lead to errors in the accurate calculation of the optical propagation in the human epidermis. This example can be extended further to the blood spatial distribution in the dermis or the roughness of the skin surface and the epidermal-dermal junction. Accurate representation of the skin geometry and physiology in the mathematical models will also result in better prediction of propagation of light in the skin.

We conclude with a recommendation to use the identified more robust datasets when estimating light propagation in human skin using Monte Carlo simulations or otherwise follow our proposed strategy to screen any new datasets to determine their biological relevance. The quantitative variability of the optical properties of the skin layers, might however only be partly reduced using a ’better’ measurement or model, as there always be the remaining variability due to a natural diversity existing between individuals of different skin type, ages and skin condition.