A library based fitting method for visual reflectance spectroscopy of human skin

The diffuse reflectance spectrum of human skin in the visible region (400–800 nm) contains information on the concentrations of chromophores such as melanin and haemoglobin. This information may be extracted by fitting the reflectance spectrum with an optical diffusion based analytical expression applied to a layered skin model. With the use of the analytical expression, it is assumed that light transport is dominated by scattering. For port wine stain (PWS) and highly pigmented human skin, however, this assumption may not be valid resulting in a potentially large error in visual reflectance spectroscopy (VRS). Monte Carlo based techniques can overcome this problem but are currently too computationally intensive to be combined with previously used fitting procedures. The fitting procedure presented herein is based on a library search which enables the use of accurate reflectance spectra based on forward Monte Carlo simulations or diffusion theory. This allows for accurate VRS to characterize chromophore concentrations in PWS and highly pigmented human skin. The method is demonstrated using both simulated and measured reflectance spectra. An additional advantage of the method is that the fitting procedure is very fast.


Introduction
Visual reflectance spectroscopy (VRS) of human skin involves measurement of diffuse reflectance in the 400-800 nm wavelength region. The measured spectrum is compared with a modelled reflectance spectrum based on a layered skin model and optical diffusion theory. Fitting parameters in the model are usually important skin chromophores such as epidermal melanin and haemoglobin (Dawson et al 1980, Feather et al 1989, Zonios et al 2001, Matcher 2002. VRS is appealing as a diagnostic method for skin characterization because it is quick, non-invasive and relatively inexpensive. To extract values for skin chromophores, a curve is fit to the measured spectrum by adapting the values of model skin parameters. A number of procedures have been employed for VRS of human skin, such as manual fitting (e.g., Svaasand et al (1995), Viator et al (2004), Zonios et al (2001)), a Nelder-Mead direct search simplex method (Douven and Lucassen 2000) and a genetic algorithm (Zhang et al 2004) to optimize the fit. Each of these procedures uses a forward calculation of the model reflectance spectrum based on an analytical solution of the optical diffusion equations. This approach enables very quick calculation of a spectrum (<1 s) but it is implicitly assumed that the diffusion approximation (i.e. light transport is dominated by scattering rather than absorption) is valid. However, for patients with port wine stain (PWS) and highly pigmented skin, this approximation tends not to be valid, therefore limiting the accuracy of VRS in these clinical situations.
In order to apply VRS for PWS and highly pigmented skin, one may use Monte Carlo simulations (Verkruysse et al 1999) to calculate diffuse reflectance spectra of skin models. However, these simulations are currently computationally too intensive and repeated forward calculations of reflectance spectra using one of the fitting procedures mentioned above may take several days and thus be impractical for clinical use.
Herein, we introduce and demonstrate the use of a novel procedure for VRS which largely avoids repeated forward calculations. Instead, we make use of a digital library which stores a large number of reflectance spectra and their corresponding skin parameters. The procedure then consists of a simple comparison of the measured spectrum with the library spectra. Thereafter, the library spectrum that fits best is used. This is a somewhat unconventional fitting method and, to our knowledge, has not been applied previously to VRS. The overall goal of this project is to assess feasibility of a library based procedure for VRS which would enable fitting with Monte Carlo based reflectance spectra. Herein, the library is not built with spectra from Monte Carlo simulations but with spectra from a diffusion theory based analytical expression. However, it will be demonstrated that a library based fitting procedure is not only feasible but fast and rigorous in identifying skin parameters from diffuse reflectance spectra as well. Results will be shown for simulated and measured (normal and PWS skin) reflectance spectra.

Skin reflectance model
To model reflectance of human skin, three components can be distinguished: (1) skin geometry, (2) optical properties of the structures in the skin model, and (3) method to calculate a corresponding reflectance spectrum.
First, we chose a simple two-layered model for skin geometry. The first and second layers represent the epidermis and dermis, respectively, and each layer is treated as optically homogeneous. This simplification allows for an analytical solution of the diffusion equations for light and, subsequent calculation of reflectance (Svaasand et al 1995, Farrell et al 1992, Douven and Lucassen 2000. Second, the optical properties (absorption coefficient µ a , scattering coefficient µ s and scattering anisotropy factor (g)) are treated as constants which then represent weighted averages of skin constituents in each layer. Several expressions for the composite optical properties of human skin layers have been published (Svaasand et al 1995, Meglinski and Matcher 2002, Douven and Lucassen 2000, Verkruysse et al 1993. In essence, they are all very similar in that they proportionally average the optical properties of the epidermis and melanin in the first layer, and those of dermal blood for all underlying layers. The spectral behaviour of chromophores such as melanin (in the epidermis) and oxy-and deoxy-haemoglobin (in the dermis) is usually expressed by convenient analytical or empirical expressions (Svaasand et al 1995, Douven andLucassen 2000). Examples for the average absorption by melanin in the epidermis µ M a , and epidermal and dermal scattering µ s , are given in equations (1) and (2), respectively. Throughout this paper, wavelength (λ) is given in nm: Values for µ a,m,694 , which is µ M a (λ = 694 nm), range from 0.2-2.5 mm −1 for light and dark skin, respectively (Norvang et al 1997). Values for µ s,577 , which is µ s (λ = 577 nm), range from 15-60 mm −1 .
Equations (1) and (2) were taken from Douven and Lucassen (2000) who derived them from the earlier work of Svaasand et al (1995). For this paper, we adopt all the expressions for composite optical properties from Douven and Lucassen with one exception. Those authors considered the µ a of skin without chromophores (µ a,T ) to be a free parameter during the fitting procedure although it was constrained to be the same for all visible wavelengths. Similarly, Svaasand et al (1995) used a constant value of 0.025 mm −1 for all visible wavelengths. We chose to adopt a wavelength-dependent expression (µ aT (λ)) proposed by Saidi (1992), which was used by Meglinski and Matcher (2002) to calculate composite skin optical properties: Saidi's original expression used a value of 7.84 × 10 7 for µ C . We will allow µ C to vary by ±10% which is reasonable given the large range (approximately one order of magnitude from 400-800 nm) covered by equation (3) and the values used by Svaasand et al (1995) and Douven and Lucassen (2000). Other model parameters involved in the calculation of composite optical skin properties are blood oxygenation and average blood vessel diameter (Svaasand et al 1995, Verkruysse et al 1997. For the purpose of this paper, which is to demonstrate the use of a library based fitting procedure, the selection of optical properties is not critical. Since the expressions of composite optical properties are exactly the same as those presented by Douven and Lucassen (except for µ aT , indicated above), we refer the interested reader to their paper.
A number of analytical expressions for skin reflectance have been published by various authors (Farrell et al 1992, Keijzer 1993, Svaasand et al 1995, Sinichkin et al 2002. For the present paper, we chose to use the Keijzer model for diffuse irradiance which is described in detail by Douven and Lucassen (2000). In this model, the deepest layer (the dermis in our case) is assumed to be infinitely thick, implicitly assuming subcutaneous tissue (fat and muscle) to be of negligible impact to the reflectance spectrum. This is a fair assumption for most visible wavelengths. It is important to realize that each of the analytical expressions is accurate only when the diffusion approximation is satisfied (Star 1995): Yoon et al (1989) determined that errors in diffusion approximation theory are relatively large for g > 0.6, which is the case for human skin. However, even for those wavelengths where the diffusion approximation is clearly violated, we ignore this violation in this study and use the analytical model. This violation may affect the accuracy of predicted tissue parameter values to some extent; however, it is not critical for testing the feasibility of the library based 10 6 µ C 7 × 10 7 -9 × 10 7 1 × 10 7 3 fitting method, the purpose of this study. Eventually, we intend to replace the library based on diffusion theory with that based on transport theory.

Building the library
The requirement for using a library based fitting method is to be able to calculate and store spectra for a large number of skin parameter combinations. Fitting a measured spectrum consists of a simple comparison with library spectra and selecting the best match. The size of the library is only limited by the size of the storage medium (e.g., a PC hard drive). We chose to start modestly and build a relatively small library of simulated spectra for a simple two-layered skin model (epidermis and dermis) while varying a total of nine skin parameters presented in table 1. The range of values for each parameter is collected from several sources in the literature (Svaasand et al 1995, Garden et al 1986, Saidi 1992, Meglinsky and Matcher 2002. With the sampling ranges and increments listed in table 1, the total number of library spectra is approximately nine million.

The fitting procedure
The spectrum that we are attempting to fit is defined as the input reflectance spectrum R inp (λ), whether it be measured or simulated. The library based method comprises two separate steps: a global and local search. The first step is a simple comparison of R inp (λ) with entries in the library of simulated spectra. The best matching library spectrum and its corresponding values for skin parameters are selected for further optimization. In the second step, we make use of interpolation and a simple linear optimization algorithm to further adapt the parameter values (from the library search) so that the modelled spectrum R model (λ) closely matches R inp (λ).
As a measure of the fit quality we consider the shape of the spectral residue (R model (λ) − R inp (λ)) as well as its Euclidean norm (equation (5)). We will refer to the latter as the average residue , where N is the total number of discrete wavelengths. In our calculations we select λ = 5 nm resulting in N = (800-450 nm)/ λ = 71 discrete wavelengths. Small or large values for  Figure 1. The wavelength ranges used to categorize the library spectra are shown in (A) along with a simulated input spectrum R inp (λ) and three selected spectra (i), (ii) and (iii) from the library. The average residue values for all the spectra found in the library for R inp (λ) are plotted versus their corresponding values for the dermal blood volume fraction B2 in (B). The circles indicate residues for the spectra shown in (A). Spectral residues for these spectra are shown in (C). Table 2. Wavelength ranges used to categorize the library spectra in files with similar spectra.
Wavelength range (nm) 450-550 550-575 650-800 450-800 Number of categories 100 100 100 100 per search key indicate a good or poor match, respectively. Following other authors (Svaasand et al 1995, Douven andLucassen 2000), we omit wavelengths 400-450 nm in our evaluation of . This omission is not critical for testing the feasibility of our library based fitting method. In this paper, we will plot spectra from 400-800 nm for illustration purposes only.

Global optimization
In the global library search for R model (λ) that matches R inp (λ), we seek those spectra that have a small value of for R inp (λ). Searching all nine million library spectra and evaluating for each would take an impractically long time. However, an essential characteristic of the library is that its components can be sorted by search keys, which enable a fast, targeted search. We chose to use four search keys to sort and categorize all spectra into smaller groups of similar spectra. Each search key is defined by a wavelength region (see table 2) over which the average reflectance value of a spectrum is evaluated. Figure 1(A) illustrates the wavelength regions 'b', 'y', 'r' and 'a' along with a simulated R inp (λ) with added Gaussian noise of amplitude 0.001. The skin parameter values for this example of R inp (λ) are presented in table 3. For a reflectance spectrum, the average values in the search key wavelength regions uniquely define a filename. For example, search key values for R inp (λ) were calculated as 0.18, 0.21, 0.42 and 0.28 for 'b', 'y', 'r' and 'a', respectively, defining the filename 'b18y21r42a28'. This filename immediately links R inp (λ) to all library spectra with the same average reflectance values in the search key wavelength regions.
Typically, a library file contains several hundred or thousand spectra. For R inp (λ) we found 2446 spectra 'b18y21r42a28'. For each of these spectra we calculated and plotted those values as a function of the corresponding value for dermal blood volume fraction of the second skin layer (B2) in figure 1(B). Ideally, one would be able to guess the value of B2 for R inp (λ) from this figure because the spectra that match R inp (λ) significantly better than others would have significantly smaller values for . However, we observe that no spectrum clearly stands out with a significant smaller residue than other spectra and that low values for are found for B2 values as low as 2% and as high as 12%.
As an illustration of the goodness of the fit between R inp (λ) and R model (λ) we selected three spectra that correspond to a very low (2%) or very high value (12%) for B2, and one arbitrary spectrum that corresponds to a relatively large average residue value. These three cases (i), (ii) and (iii) are plotted in figure 1(A) along with R inp (λ) and indicated by the circles in figure 1(B). Figures 1(A) and (B) show that with just the library search we cannot obtain a good prediction for the value of B2 for R inp (λ). However, the spectral residues (figure 1(C)) illustrate that the library spectra deviate from R inp (λ) by larger values than just the noise level of R inp (λ) which suggests that further optimization of the skin parameters may be achieved.

Local optimization
In the following, we use R lib (λ, P ) to refer to a library spectrum where P represents the set of parameters P i (epidermal thickness, epidermal melanin concentration, dermal blood volume fraction, etc). The group of spectra found in a library file for R inp (λ) will be referred to as {R lib (λ, P )}.
Each spectrum in {R lib (λ, P )} looks similar to R inp (λ). However, since they were grouped together because of similar average values, a considerable difference in specific wavelength regions for some spectra may still exist (see spectrum (iii) in figures 1(A) and (C)). Since the current library is relatively small, it would be a coincidence if it contained one spectrum that perfectly matched R inp (λ). If one spectrum R 0 (λ, P 0 ) is selected out of the group {R lib (λ, P )}, a better fitting spectrum can be found by slightly modifying the skin parameters P 0 .
For optimization of the parameters P 0 , we first calculate the partial differential spectra for the chosen library spectrum R 0 (λ, P 0 ). For each parameter P i in the set P 0 , we define a partial The spectral residue (noisy curve) is fitted with a linear combination of differential spectra (smooth curve) by optimizing scaling values k i . differential spectrum according to equation (6): We use the analytical model to calculate R 0 (λ, P i ) and R 0 (λ, P i + P i ) for very small values of P i to obtain a good estimate of the sensitivity of R(λ) for a change in each parameter. We chose the spectrum indicated by (ii) in figure 1(A) to serve as R 0 (λ, P 0 ) in the following discussion of local optimization. Differential spectra for eight skin parameters are shown in figure 2(A). The next step is to consider the spectral residue for R inp (λ) and library spectrum R 0 (λ, P 0 ) as shown in figure 1(C) (indicated by (ii)). The concept of the optimization here is to minimize the difference between a newly calculated model spectrum R 1 (λ, P 1 ) and R inp (λ) by finding an optimal linear combination of differential spectra ∂R(λ) ∂P i,0 . Mathematically, we aim to find an optimal set of scaling values k i , such that equation (7) is minimized: Table 4. Skin parameters from the library before (P i,0 ) and after optimization (P i,1 ) using the scaling values k i . Since the differential spectra are calculated with P i,0 equal to 1% of the value of P i,0 , a value of −51.5 for k i means that the corrected value P i, 1 is 51.5% lower than the original value of P i, 0 . Range limits for the skin parameters are also listed. Once the optimal values for k i (table 4) have been found, corrected skin parameters P i,1 are calculated according to equation (8), For this study, we used a standard algorithm in Matlab R to optimize k i , but in principle any optimization algorithm (e.g. simplex, genetic algorithm) could be used.
In figure 2(B), the spectral residue of figure 1(C) is shown along with the linear combination of partial differential spectra for scaling values k i obtained by minimizing equation (7). The spectral residue can be approximated nicely with the differential spectra of figure 2(A).
Constraints on the skin parameter values (see table 4) are implemented by restricting corresponding k values during the optimization. The spectrum R 1 (λ, P 1 ) for the new set of parameters P 1 fits R inp (λ) better than R 0 (λ, P 0 ); the average residue value is reduced from 0.0053 to 0.0035. This procedure can be repeated several times until the best fitting spectrum (smallest ) is determined. With our current small library, we found that up to ten repeated local optimization steps were needed. With a more extensive library this number can likely be reduced.
The above procedure of local optimization was done for a spectrum with a constant B2 value of 12% (the k-value for B2 was not allowed to vary). We performed similar local optimizations for library spectra with B2 values of 2, 3, . . . , 11%. These spectra are represented by the squares in figure 3(A) and are the same as the spectra in figure 1(B) with the lowest average residues. Each of these spectra is used as R 0 (λ, P 0 ) in further optimizing skin parameters by attempting to minimize the residue value . For each value of B2, it was possible to find a better match with R inp (λ) by letting the multi-linear optimization algorithm alter the seven skin parameters. The optimized spectra for each value of B2 are represented by the diamonds connected by the curve in figure 3(A). The curve now clearly shows a 'dip' for a B2 value of 5% which matches the B2 value of R inp (λ) (see table 3). Because the optical absorption and scattering spectra of skin are relatively smooth, we do not expect that the optimization method can fit the high-frequency features of the Gaussian noise. Therefore, the minimum residue value of approximately 0.001 corresponds to the level of noise that was added to the simulated spectrum.
The curve in figure 3(A) also shows that the residues for all values of B2 have been reduced significantly by the local optimization indicating that even very low or high B2 values can produce a reasonably nice fit to R inp (λ) when other skin parameters are adapted. To

VRS measurements
VRS measurements were performed with a commercial spectroscopy device (Ocean Optics Inc., Orlando, FL). Light from a halogen source (LS-1) is guided into an integrating sphere (model ISP-REF 38.1 mm diameter) to illuminate diffusely the skin through a 10 mm diameter aperture. The diffuse reflectance is collected through the same aperture. A switch on the integrating sphere was set to exclude the contribution of specular reflectance to the measured signal. Light from the integrating sphere is then guided to a spectrophotometer (model SD-2000) via an optical fibre. Before each measurement the spectroscopy device was calibrated using 2 and 99% reflectance standards (Labsphere, North Sutton, NH). In order for the spectroscopic device to thermally stablize, the halogen source was switched on at least 30 min before a measurement was taken. We always took care to place the integrated sphere on the skin gently and with minimal pressure, thus avoiding possible occlusion of microvasculature.
We will present measurements from two PWS patients whom we will refer to as patients I and II. Patient I is Caucasian and patient II is Asian. The PWS of patient I was located on the right cheek and was light pink in appearance. The PWS of patient II was on the upper left side of the face and was dark red in appearance. Measurements on PWS and adjacent normal skin were taken for each patient. Each measurement was repeated at least three times.

Reflectance measurements on normal and PWS skin
Figure 4(A) shows two examples of measurements (noisy curves) from normal skin for patients I and II plus the corresponding model curves that were found to be the best fits using our method. The general reflectance for patient II (Asian) is expected to be lower than that for patient I (Caucasian) because of the higher skin pigmentation in the former. The best fits correspond to . Measured spectra on normal skin of two PWS patients plus their corresponding best fits using our library based fitting method are shown in (A). The dermal blood volume fractions at the best fits are 5% and 12% for patients I and II, respectively, are shown in (B). Each curve in (B) represents the fitting results for one reflectance measurement. Only one reflectance spectrum for each patient is shown in (A). dermal blood volume fractions B2 of 5 and 12% for patients I and II respectively. These values were obtained from the minima in the residue/B2 curves shown in figure 4(B). Each curve represents the fitting results for one measured spectrum. The fact that the curves are similar reflects that the reflectance measurements of normal skin were well reproducible. Please note that only one measured reflectance spectrum for each patient is shown in figure 4(A).
Other skin parameter values that were found as the best fit are presented in table 5. The epidermal melanin content for patient II (Asian) is indeed found to be higher than that for patient I (Caucasian) as expected. Figures 5(A) and (B) show similar results for the PWS sites of patients I and II. From figures 5(A) and (B) one can see that the fits are not as good as those for normal skin. The minimum residue values found are 0.0088 and 0.0066 for patients I and II, respectively. The reflectance spectra for PWS on both patients were not well reproducible which explains that the curves in figure 5(B) are not as similar to one another as the curves in figure 4(B). Based on the curves in figure 5(B), dermal blood volume fractions are estimated to be 20% or higher and approximately 20% for patients I and II, respectively. The dermal blood volume fractions at the best fits are 20% or higher and 18% for patients I and II, respectively, as can be found from the residue plotted versus B2 in (B). Each curve in (B) represents the fitting results for one reflectance measurement. For each patient, we show only the best fitting reflectance spectrum in (A).

Reflectance measurements
Even though we think that the library found the best fits possible, we acknowledge that many of the corresponding skin parameter values are unrealistic. The blood volume fractions, for both normal and PWS skin, appear to be too high. Similarly, the oxygenation values, at 90% or higher, appear too high. For the PWS of patient II we find an epidermal melanin concentration almost twice as high as the adjacent normal skin. Since PWS is classified as a vascular and not a pigmented lesion, it is not likely that the pigmentation of the PWS skin is indeed higher than for normal skin on the same patient. An alternative fitting procedure which was developed in our laboratory is based on a genetic algorithm and has found similar unrealistic values for the same measured spectra (Zhang et al 2004). We propose several explanations for the unrealistic values predicted by the fitting method. First, the spectral behaviour of some of the skin constituents may not necessarily be accurate. For example, a range of slightly different curves, tabulated values or approximations for the absorption of oxy-and deoxy-haemoglobin (e.g. van Kampen and Zijlstra (1965), Takatani and Graham (1987), Wray et al (1988), Douven and Lucassen (2000)) are available. Similarly, the value of µ a,T is not well known and, consequently, the value used in different studies varies greatly. The power with which melanin absorption decreases with wavelength varies from −3.3 (Meglinski and Matcher 2002) to −4 in the literature (Svaasand et al 1995). Preliminary analysis shows that these inaccuracies or uncertainties have a significant effect on the modelled spectra and, therefore on the corresponding skin parameters.
Second, the high dermal blood volume fractions clearly invalidate the use of the diffusion approximation. The approximation error is strongly wavelength dependent, and, consequently so is the error in reflectance level prediction (see Douven and Lucassen (2000)). As a result, the shape of the spectra is somewhat distorted and the skin parameters that follow from a fit are distorted as well.
Third, the numerical aperture used in the study is limited to 10 mm, resulting in loss of a portion of the diffuse reflected light. The error introduced by this effect is greater when the scatterer path lengths in skin are large which is predominantly in the red part of the spectrum (e.g., Takiwaki et al (2002)). Similar to the violation of the diffusion approximation, this error affects the values for parameters that result from a fit. Monte Carlo simulations or improved analytical models of spectral reflectance may be required to account for the limited numerical aperture as well as the discussed violations of the diffusion approximation. Finally, it is likely that the objective in the optimization procedure (which is to minimize the residue) is not optimal for an accurate prediction of some of the skin parameters. In our definition of the residue (equation (5)) the difference between modelled and measured reflectance has equal weight at all considered wavelengths. If certain wavelength regions were to be given more weight in the calculation of the residue, the accuracy of predicted chromophore concentrations might improve. In fact, initial studies (Feather et al 1989, Dawson et al 1980 suggest that using only a few well-chosen wavelengths rather than the entire spectrum would improve the predictive value of VRS. The method presented here can be easily adapted to give certain wavelengths or wavelength regions a larger weight in the minimization by simply modifying equation (5). This work is ongoing.

Limiting the range of solutions by constraining parameter values
The results in figure 3(A) illustrate that combined global and local optimization can reproduce the original dermal blood volume fraction (B2 = 5%) of the simulated spectrum. However, it also shows that optimization at B2 values significantly different from 5% can still produce relatively small residues. This is because considerable freedom exists for other skin parameters to compensate for changes in the reflectance caused by a discrepancy in the B2 value. For example, in the local optimization for the results in figure 3, oxygenation levels were allowed to vary between 60-100%. Dermal scattering µ s,d 577 nm was allowed to vary from 15-50 mm −1 . Figure 3 suggests that a near perfect fit with an average residue value of 0.0025 would find a range for B2 values between 4 and 10%. In practice, it is very difficult to find fits of such high quality.
Intuitively, it seems plausible that when information about the range of skin parameters is available (e.g., from literature or additional measurements with other diagnostic methods), one could narrow the range of possible solutions. This can be demonstrated with the spectrum and corresponding library search data of figure 3. In the second optimization step, the oxygenation is allowed to vary between realistic values of 75 and 85%. To limit dermal scattering, we choose two ranges, µ s,d 577 nm = 25-30 mm −1 or 35-40 mm −1 . Figure 6 shows that limiting the ranges of the above two parameters indeed results in curves with a more distinct 'dip' and thus a smaller range for B2 values at a certain residue value. This range represents the measurement error or accuracy range of using VRS to 'measure' B2. Often, the accuracy of values, determined with VRS is not reported in the existing literature. This may be because previously used fitting procedures aim at finding a single solution rather than exploring the entire solution space. The presented library method automatically generates information on the range, and thus estimated accuracy of the determined parameter values.
When attempting to improve the accuracy of the value for B2 by limiting the range of other parameters, care should be taken because an incorrect restriction of dermal scattering results in a 'dip' for B2 = 8%, an overestimation of the dermal blood fraction.

Conclusions
Our method is capable of quickly and correctly finding skin parameters for a simulated spectrum. However, it did not do well at finding realistic values for all of the skin parameters on measured spectra. We have outlined possible explanations for the latter.
We conclude that VRS for PWS diagnostics presents challenges of different natures. Violation of the conditions for the diffusion approximation is the easiest challenge to overcome. An improved diffusion model or a Monte Carlo based library would address that problem. An immediate advantage of the fitting procedure, described here, however, is that it is objective and much faster than existing methods. Speed is important for real-time VRS diagnostics or hyper-spectral imaging in which many spectra are measured simultaneously.
Two other challenges in VRS are uncertainty in the spectral behaviour of optical properties and the ill-posedness of the inverse problem. Both are somewhat related to each other because an uncertainty in skin parameters increases the ill-posedness of the problem. More freedom to vary a skin parameter allows for more ways to compensate for discrepancies in the spectrum induced by skin parameters. This results in a larger range of values for each parameter that satisfactorily fits the spectrum, i.e., the ill-posedness of the problem is increased. We would like to emphasize that the ill-posedness of VRS has not been discussed elaborately in the literature and curves such as in figures 3 and 6 have not been shown before. This may be due to the fact that when a single reasonable fit is obtained (e.g. through manual fitting) it is not obvious that there may be many alternative solutions as well. Existing methods do not systematically investigate the solution space but work towards one particular solution. Therefore they do not indicate the robustness of a skin parameter value found from a fit. With our method, the systematic investigation of robustness is much easier and an impression of the robustness, or measurement error for a parameter value, is automatically obtained.
In future studies, we intend to enhance the diagnostic capability of VRS by using more accurate spectral properties of the skin constituents. We will further investigate the illposedness for parameters other than B2, and explore ways to reduce it by limiting the range of values for some of the parameters.