Translator Disclaimer
1 September 2008 Confidence intervals on fit parameters derived from optical reflectance spectroscopy measurements
Author Affiliations +
Abstract
We validate a simple method for determining the confidence intervals on fitted parameters derived from modeling optical reflectance spectroscopy measurements using synthetic datasets. The method estimates the parameter confidence intervals as the square roots of the diagonal elements of the covariance matrix, obtained by multiplying the inverse of the second derivative matrix of Χ2 with respect to its free parameters by Χ2/v, with v the number of degrees of freedom. We show that this method yields correct confidence intervals as long as the model used to describe the data is correct. Imperfections in the fitting model introduces a bias in the fitted parameters that greatly exceeds the estimated confidence intervals. We investigate the use of various methods to identify and subsequently minimize the bias in the fitted parameters associated with incorrect modeling.

1.

Introduction

Fitting in vivo optical spectroscopic data with mathematical models that describe the underlying physical “principles” is common practice in biomedical optics. Examples of such analyses include the fitting of multicomponent fluorescence spectra in tissue (that can include tissue autofluorescence as well as exogenous fluorophores), the fitting of multicomponent Raman spectra to determine the tissue biochemistry, and the fitting of multicomponent extinction spectra of tissue chromophores using reflectance spectroscopic measurements. In vivo fluorescence, Raman, and reflectance spectroscopic measurements are all affected by the tissue optical properties. The effect of the optical properties on the measured spectra must first be accounted for using physical light-tissue interaction models to obtain the “intrinsic” fluorescence, Raman, and extinction spectra. These undistorted spectra can then be analyzed with linear component analysis techniques such as singular value decomposition or other least squares minimization techniques. However, if the physical models used in the “preprocessing step” are inaccurate, an error in the fitted component values will be introduced. Furthermore, even if the algorithms used to extract the intrinsic optical spectra from the measured spectra are completely accurate, a lack of knowledge of either the in vivo spectral shape or even the physical presence of some of the components may introduce another source of error in the analysis of in vivo optical spectroscopic measurements. Finally, even when the intrinsic spectra are determined correctly and full knowledge about the presence and spectral shapes of all the components is available, another source of error (uncertainty) in the analysis of optical spectroscopic measurements is introduced by a limited quality of the data (signal-to-noise ratio) that can prevent accurate determination of some or all of the components in the component analysis. The latter type of error is called a statistical error, while the first two types of errors are considered as systematic errors. These errors will introduce an uncertainty in the fitted parameter values extracted from the component analysis. Currently, it is not common practice in the field of biomedical optics to estimate the accuracy of the fitted parameter values from single spectroscopic measurements. However, we believe that it is useful to at least determine the statistical error in each of the fitted parameter values. For example, in the case of averaging multiple in vivo measurements, calculating the weighted average of the parameters using the statistical errors as weight factors would be more appropriate than simply calculating the unweighted average and ignoring the statistical uncertainties associated with the quality of the fits of the different spectra. Moreover, calculating the statistical errors in the fit parameters facilitates objective assessment of the quality of spectra and fits with the possibility to reject poor quality spectra or fit values from an in vivo dataset.

In principle, it is possible to calculate the parameter uncertainty by repeating the measurements and analyzing the standard deviation of parameters extracted from each of these multiple measurements. However, this is not practical for in vivo applications such as medical diagnosis, in which typically only one or a few measurements can be acquired. Scepanovic 1 recently published an error analysis for single Raman spectroscopic measurements. They presented an analytical formula for estimating uncertainty expressed as a function of measurement noise, signal strength, and spectral overlap. Their analysis characterizes uncertainty for linear fitting functions such as Raman and fluorescence spectroscopy, but has not been tested for nonlinear fitting functions encountered in reflectance measurements. In this work we validate a general method for determining the statistical error or confidence interval (CI) on fitted parameters derived from modeling single optical reflectance spectroscopic measurements using synthetic datasets. The method estimates the parameter CIs as the square roots of the diagonal elements of the covariance matrix, which is obtained by multiplying the inverse of the second derivative matrix of χ2 with respect to its free parameters by χ2v , with v being the number of degrees of freedom. While the method is in principle valid for all spectroscopic techniques (fluorescence, Raman, and reflectance), we demonstrate the validity of the approach using reflectance spectroscopy as an example. Reflectance spectroscopy is sensitive to the absorption and scattering properties of tissue. The absorption coefficient is proportional to the concentration of the tissue chromophores, and the scattering coefficient reflects the size and density of scattering structures in the tissue. Several methods have been employed to extract physically relevant information (such as blood volume fraction and blood saturation) from reflectance spectra, e.g., using an analytical approximation of the transport equation such as the diffusion equation,2, 3, 4, 5, 6 Monte-Carlo-based models,7, 8 or empirical methods.9, 10, 11 Irrespective of the approach, a model function minimization or fitting routine must always be used to find the optimal values for the free parameters of the models. In general, a least squares (χ2) optimization routine is used where the (weighted) square of the difference between the data points and the fitted curve is minimized. Most commercial software packages that employ least squares fitting routines automatically generate the covariance (also called variance-covariance or error) matrix. In principle, the diagonal elements of this matrix represent the variances of the parameters if 1. the fitting function is linear in the parameters and 2. the measurement uncertainties of the individual data points are known exactly. However, neither of these points is fulfilled for a single reflectance spectroscopy measurement. Nevertheless, we demonstrate that our simple method to estimate the statistical error of the fit parameters from single reflectance spectroscopic measurements yields correct CIs as long as the model used to describe the data is correct. Furthermore, we analyze the effect of imperfections in the fitting model (systematic errors) and show that even minor imperfections in the model may introduce a bias in the fitted parameters that greatly exceeds the estimated statistical CIs. Finally, we investigate the use of methods to identify and subsequently minimize the bias in the fitted parameters associated with incorrect modeling.

2.

Theory

In this work we use the empirical model used to analyze a particular type of reflectance spectroscopy that we have been investigating over the past 5years , differential path-length spectroscopy (DPS),11, 12, 13, 14 as an example. The model to which the differential path-length spectra are fitted can be written as:

Eq. 1

R(λ)=[a1(λλ0)a2+a3(λλ0)4]exp[0.4μatotal(λ)].
The scattering function (in square brackets) is modeled by a combination of Mie and Rayleigh scattering, given by power law functions with amplitudes a1 and a3 , and wavelength dependencies (λλ0)a2 and (λλ0)4 , respectively. Here λ0 is a normalization wavelength, which we usually set to 800nm . The differential reflectance signal is attenuated due to the presence of absorbers following Lambert-Beer’s law with a path length equal to the fiber diameter,15, 16 which in this exercise is assumed to be 0.4mm . The absorption coefficient μatotal(λ) is the sum of the absorption coefficients of the chromophores present in the interrogation volume, which in this example are assumed to be betacarotene and blood:

Eq. 2

μatotal(λ)={a4μaBcar(λ) +a5[a6μaHbO2(λ)+(1a6)μaHb(λ)] (1exp{a7[a6μaHbO2(λ)+(1a6)μaHb(λ)]}a7[a6μaHbO2(λ)+(1a6)μaHb(λ)])}.
Parameter a4 represents the concentration of betacarotene, a5 is the blood volume fraction, a6 is the microvascular blood oxygenation, and a7 is the average vessel diameter. Input spectrum μaBcar(λ) is the specific absorption coefficient of betacarotene,17 μaHbO2(λ) is the absorption coefficient of fully oxygenated whole blood, and μaHb(λ) is the absorption coefficient of fully deoxygenated whole blood.18

With our experimental setup,11, 12, 13, 14, 15, 16 the data acquired in each reflectance spectrum consist of N=2048pixel values Rdata(j) corresponding to the differential reflectance at wavelengths ranging from 350to1000nm . Since the transmission efficiency of the setup is different for different wavelengths, the signal-to-noise ratio (SNR) is wavelength dependent. Clearly it is desirable to put less weight on pixels with a low SNR in the fitting routine. To accomplish this, the spectra are first smoothed by binning p pixel values and using the Np=n average values Rdata(i) and corresponding standard deviations s(i) in the fitting routine. In our simulations, we have typically used a bin width of 15pixels , which corresponds to an optical resolution of 6nm , but we have also tested the effect of smaller and larger bin widths (10 and 30pixels , respectively). The data are then fitted to Eq. 1 using a least squares Levenberg-Marquardt fitting routine using the standard deviations of the binned data points as weight factors. Chi-square (the quantity that is minimized) is given by

Eq. 3

χ2=i=1i=n[Rfit(i)Rdata(i)s(i)]2,
where Rfit(i) are the fitted reflectance values according to Eq. 1, Rdata(i) are the corresponding measured bin values, and s(i) are the standard deviations of the bins. Note that in this approach, it is critical to bin the input spectra μaBcar(λ) , μaHbO2(λ) , and μaHb(λ) in exactly the same way as the measurement spectra. Furthermore, it is important to keep the bin width small enough so that the approximation
Rdata(i)=1pjj+pRdata(j)1pjj+pexp[0.4μatotal(j)]exp[0.41pjj+pμatotal(j)],
holds. Particularly for wavelength regions where the absorption coefficient varies rapidly, e.g., around the Soret absorption band of hemoglobin, this approximation will become invalid for large bin widths.

The minimum value of χ2 corresponds to values for the seven fit parameters (a1 to a7) that best describe the data. In principle, the CI of each of these parameters is given by the square root of the diagonal elements of the covariance matrix. The covariance matrix is the inverse of the second derivative matrix of χ2 with respect to its free parameters and is exact, if and only if; 1. the fitting function is linear in its parameters, and 2. the weight factors used in the minimization routine are exact.19, 20 In the case of our nonlinear model Eq. 1 with estimated weight factors s(i) , we estimate the covariance matrix by multiplying the inverse of the second derivative matrix of χ2 with respect to its free parameters by χ2v , with v being the number of degrees of freedom.20 The CI of each of the parameters is given by the square root of the diagonal elements of the estimated covariance matrix under the assumption that the model [Eq. 1] is correct. We validate that the estimated covariance matrix is reasonably accurate for our nonlinear fitting problem in Eq. 1.

The calculated CI represents the statistical error only and does not fully account for systematic errors due to inaccuracies in the physical model used to fit the data. For example, if an additional absorber was present in the sampling volume but not in the model fit function, or if the wavelength dependence of the scattering function is not modeled correctly, the fitted parameters will be biased with respect to their true values. We investigate the magnitude of the bias in the fitted parameters for the case of 1. a missing absorber in the model and 2. an error in the modeled scattering function, and investigate the use of methods to identify and subsequently minimize the bias in the fitted parameters associated with incorrect modeling.

3.

Simulations

We have generated 315 synthetic datasets, each consisting of 2048pixel reflectance values Rsyn(j) using Eqs. 1, 2 with the parameters as described in Table 1 .

Table 1

Parameter values used for the synthetic datasets.

ParameterParameter descriptionValues
a1 (−)Mie scattering amplitude1
a2 (−)Mie scattering slope 1
a3 (−)Rayleigh scattering amplitude0; 0.01; 0.1
a4 (μM) Betacarotene concentration0; 5; 30
a5 (%)Blood volume fraction0; 0.2; 1.0; 2.0; 5.0; 10; 20
a6 (%)Blood saturation0; 25; 50; 75; 100
a7 (μm) Average vessel diameter10

For each combination of parameter values, 100 different fits are performed by generating noise on the synthetic spectrum. The noise was generated by adding a random number between 1 and 1, multiplied by a specified noise amplitude to the pixels: Rdata(j)=Rsyn(j)+noiseamprandom[1,1] . Note that the noise generated in this way is independent of wavelength. Since the data Rsyn(j) are wavelength dependent, the relative noise amplitude is in fact wavelength dependent, which is typically the case for optical spectroscopic measurements. For our experimental setup, typical signal-to-noise ratios at λ0=800nm are better than 100, and we have used noiseamp=0.01 for most of the simulations, but also checked the validity of our approach for larger noise amplitudes (0.1 and 0.5). By definition, if the calculated confidence intervals are correct, the fitted parameter values should be within 2 CIs in 95 out of the 100 fits.

4.

Results

4.1.

Fitting with the Correct Model

We first analyze whether the CIs calculated by taking the square root of the diagonal elements of the estimated covariance matrix are correct, as a function of bin width and noise amplitude. The fits are performed using the correct model described in Eqs. 1, 2 with seven free parameters. Furthermore, we calculate the magnitude of the CIs as a function of blood volume fraction for this case.

4.1.1.

Influence of bin width

Table 2 shows the number of times each of the fitted parameter values is within 2 CIs of their true values, out of the 100 fits and averaged for the 315 synthetic datasets (standard deviations are indicated in parentheses) for different bin widths. If this number is close to 95, the CI was estimated correctly. If this number is much smaller, the CI was underestimated. The noise amplitude was set to 0.01.

Table 2

Average number of times (out of 100) each of the fitted parameter values is within 2 CIs of their true values when the correct model is used, for bin widths of 10, 15, and 30pixels . Standard deviations are indicated in parentheses. Noise amplitude=0.01 .

101530
a1 94 (4)96 (4)96 (7)
a2 91 (4)94 (4)94 (5)
a3 96 (3)97 (4)95 (8)
a4 94 (3)95 (5)84 (25)
a5 94 (3)95 (3)95 (7)
a6 94 (3)95 (5)85 (22)
a7 93 (4)95 (3)94 (9)

Table 2 shows that on average, close to 95% of the fitted parameter values is within 2 CIs, except for the largest bin width for parameters a4 (betacarotene concentration) and a6 (blood saturation). For this large bin width, the optical resolution of the binned data (12nm) becomes larger than the full width at half maximum (FWHM) of the narrowest betacarotene feature (8nm) , causing an underestimation in the CI of the betacarotene concentration a4 . Since betacarotene absorbs in the same wavelength region as oxyhemoglobin and deoxyhemoglobin, the accuracy of the blood saturation a6 is affected as well. From Table 2 it is clear that the best results were obtained for a bin width of 15, which is used in the remainder of this work.

4.1.2.

Influence of noise amplitude

Table 3 shows the number of times each of the fitted parameter values is within 2 CIs of their true values, out of the 100 fits and averaged for the 315 synthetic datasets (standard deviations are indicated in parentheses) for different noise amplitudes and a bin width of 15pixels . On average, close to 95% of the fitted parameter values is within 2 CIs, even for very high noise amplitudes. Table 4 shows the average CI for each of the parameters as a function of blood volume fraction a5 for the three different noise amplitudes. The CI increases with increasing noise amplitude, as expected. Thus for large noise amplitudes (0.1 and 0.5, corresponding to a signal-to-noise ratio at λ=800nm of 10 and 2, respectively), the fitted parameters are also correctly within 2 CIs in not significantly less than 95% of the simulations, because the CIs are dramatically increased compared to the low noise amplitude (0.01, or signal to noise of 100). Note that for the low noise amplitude, which is typical for our in vivo data,11, 12, 13, 14 the CI of the fitted blood saturation a6 is smaller than 1% for blood volume fractions a5> 1% . Even for a blood volume fraction as low as 0.2%, the CI of a6 is only 5% when the correct model is used in the fit. These numbers indicate that our instrument can extract the in vivo microvascular saturation with high accuracy, under the assumption that the model used to fit the data is correct.

Table 3

Average number of times (out of 100) each of the fitted parameter values is within 2 CIs of their true values when the correct model is used, for noise amplitudes of 0.01, 0.1, and 0.5. Bin width=15pixels .

0.010.10.5
a1 96 (4)93 (3)93 (3)
a2 94 (4)93 (3)93 (3)
a3 97 (4)93 (3)93 (3)
a4 95 (5)93 (3)92 (3)
a5 95 (3)93 (4)91 (6)
a6 95 (5)94 (2)94 (3)
a7 95 (3)92 (3)92 (3)

Table 4

Average CI for each of the parameters as a function of blood volume fraction a5 for the three different noise amplitudes.

Noiseamplitudes0%0.2%1%2%5%10%20%
a1 [−]0.010.100.500.30.52.50.00030.0022.40.00030.0020.010.00030.0020.010.00040.0020.010.00060.0030.020.0010.010.05
a2 [−]0.010.100.500.0020.010.070.0010.010.070.0010.010.060.0010.010.060.0020.010.070.0020.020.090.0040.030.2
a3 [−]0.010.100.500.00050.0020.080.00020.0010.020.00020.0010.0060.00020.0010.0060.00030.0020.0080.00050.0030.010.0010.010.05
a4 [μM] 0.010.100.500.150.94.50.10.94.60.11.04.90.21.15.30.31.57.00.42.3110.96.030
a5 [%]0.010.100.50 104 104 105 0.0080.07 104 0.0090.080.40.010.080.40.010.10.50.020.140.70.040.31.4
a6 [%]0.010.100.50 106 108 107 5.034 107 1.07370.64200.42.4110.42.0100.42.110
a7 [μm] 0.010.100.50 106 107 107 0.86 103 0.21.370.10.740.050.420.040.31.70.080.73.4

Figures 1a and 1b show a typical fit and the residuals (data-fit) , respectively, for the synthetic spectrum with a3=0.01 , a4=5μM , a5=2% , and a6=75% for a noise amplitude of 0.01 and a bin width of 15. The fitted values for the parameters and the CIs were a1=1.0003±0.0002 , a2=1.000±0.001 , a3=0.0099±0.0002 , a4=5.0±0.2μM , a5=2.00±0.01% , a6=74.9±0.7% , and a7=10.04±0.05μm . All fit values are unbiased and within 2 CIs of the true values. The residual spectrum shows only noise and is featureless.

Fig. 1

(a) Typical fit and (b) the residuals for the synthetic spectrum with a3=0.01 , a4=5μM , a5=2% , and a6=75% for a noise amplitude of 0.01 and a bin width of 15pixels .

054044_1_034805jbo1.jpg

In conclusion, when the correct model is used to fit the data, the fitted parameter values are unbiased, and the CIs are estimated correctly using the diagonal elements of estimated covariance matrix obtained by multiplying the inverse of the second derivative matrix of χ2 with respect to its free parameters by χ2v , with v being the number of degrees of freedom.

4.2.

Fitting with the Wrong Model: Missing Absorber

In this section, we investigate the magnitude of the bias in the fitted parameters introduced when the fitting model is missing an absorber. Two different situations are investigated: the absorption band of the missing absorber has 1. little or 2. significant spectral overlap with the other absorbers in the component fit. In Sec. 4.2.1 the missing absorber is an artificial Gaussian centered at 700nm with a width of 25nm . Such an absorber has only little overlap with the main absorption bands of oxyhemoglobin and deoxyhemoglobin. In Sec. 4.2.2 the missing absorber is betacarotene, which does overlap with the absorption bands of hemoglobin. Figure 2 shows the normalized absorption coefficients of the four absorbers used in the simulations. The noise amplitude is set to 0.01 and the bin width is 15 in this section.

Fig. 2

Normalized absorption coefficients of the four absorbers used in the simulations.

054044_1_034805jbo2.jpg

4.2.1.

Nonoverlapping Absorber: Gaussian Centered at 700nm

An additional 210 synthetic datasets were generated using the Gaussian absorber instead of betacarotene in Eq. 2. Fitting was performed with only six free parameters by excluding a4μaBcar(λ) from Eq. 2. Figures 3a and 3b show a typical fit and the residuals, respectively, for the synthetic spectrum with a3=0.01 , a4=5 , a5=2% , and a6=75% . This “concentration” (a4=5) of Gaussian absorber results in a hardly visible maximum signal decrease of 2.5% at λ=700nm . The fitted values for the parameters and the CIs were a1=0.9977±0.0007 , a2=0.992±0.004 , a3=0.0098± 0.0005 , a5=1.88±0.03% , a6=74±2% , and a7=9.6± 0.2μm .

Fig. 3

(a) and (c) typical fit and (b) and (d) the residuals for the synthetic spectrum with a3=0.01 , a4=5 and 30, a5=2% , and a6=75% , fitted without the Gaussian absorber in the fitting model.

054044_1_034805jbo3.jpg

Figures 3c and 3d show a typical fit and the residuals, respectively, for the synthetic spectrum with a3=0.01 , a4=30 , a5=2% , and a6=75% . This “concentration” (a4=30) of Gaussian absorber results in a clearly visible maximum signal decrease of 14% at λ=700nm . The fitted values for the parameters and the CIs were a1=0.982±0.004 , a2=0.93 ±0.02 , a3=0.011±0.003 , a5=1.3±0.2% , a6=62±11% , and a7=7.4±0.8μm . Even though the statistical CIs of these fits are much larger than the statistical CIs of the parameters using the correct model (Sec. 4.1.2), only the Rayleigh amplitude (a3) and saturation (a6) fit values are within 2 CIs of the true values. The residual spectra very clearly show the features of hemoglobin absorption as well as the missing Gaussian absorber.

Table 5 shows the number of times each of the fitted parameter values is within 2 CIs of their true values, out of the 100 fits and averaged for the 105 synthetic datasets for different concentrations of Gaussian absorber.

Table 5

Average number of times (out of 100) each of the fitted parameter values is within 2 CIs of their true values for different concentrations of missing Gaussian absorber.

0530
a1 96 (4)6 (21)3 (16)
a2 94 (3)16 (36)28 (44)
a3 97 (4)86 (31)78 (40)
a4
a5 95 (3)2 (7)3 (15)
a6 95 (5)80 (36)71 (44)
a7 95 (3)1 (3)1 (3)

Table 6 shows the average bias ( absolute value of absolute difference between true and fitted parameter values) and CI for each of the parameters for different concentrations of Gaussian absorber. Blood volume fractions a5<1% were excluded from the calculation of the averages due to the extremely large bias and CI of the blood parameters ( a5 , a6 , and a7 ) for small blood volume fractions (see for example Table 4). It is observed that, except for parameters Rayleigh amplitude (a3) and saturation (a6) , the average bias in the parameters is much larger than twice the average statistical CI of the parameters.

Table 6

Average bias and CI for each of the parameters for different concentrations of missing Gaussian absorber.

Gauss absorberconcentrationAverage biasAverage CI
a1 [−]5300.0050.030.0020.008
a2 [−]5300.020.130.0070.04
a3 [−]5300.00160.0090.00150.007
a4
a5 [%]5300.21.00.060.3
a6 [%]5301.391.59
a7 [μm] 5300.94.60.31.5

4.2.2.

Overlapping Absorber: Betacarotene

Similar to Sec. 4.2.1, fitting was performed with only six free parameters by excluding a4μaBcar(λ) from Eq. 2. Figures 4a and 4b show a typical fit and the residuals, respectively, for the synthetic spectrum with a3=0.01 , a4=5μM , a5=2% , and a6=75% . This concentration (a4=5μM) of betacarotene results in a hardly visible maximum signal decrease of 2.5% at λ=450nm . The fitted values for the parameters and the CIs were a1=0.9987±0.0006 , a2=0.984±0.003 , a3=0.0103±0.0004 , a5=2.02±0.03% , a6=91±1% , and a7=9.9±0.1μm . Even though the statistical CIs of this fit are larger than the statistical CIs of the parameters using the correct model (Sec. 4.1.2), only half of the fit values ( a3 , a5 , and a7 ) are within 2 CIs of the true values. The fitted saturation (a6) is very biased with respect to the true value (91 versus 75%), while the CI is increased only marginally compared to the case where betacarotene was included in the fit. The residual spectrum clearly shows the features of hemoglobin absorption as well as betacarotene absorption.

Fig. 4

(a) and (c) typical fit and (b) and (d) the residuals for the synthetic spectrum with a3=0.01 , a4=5 and 30μM , a5=2% , and a6=75% , fitted without betacarotene in the fitting model.

054044_1_034805jbo4.jpg

Figures 4c and 4d show a typical fit and the residuals, respectively, for the synthetic spectrum with a3=0.01 , a4=30μM , a5=2% , and a6=75% . This concentration (a4=30μM) of betacarotene results in a clearly visible maximum signal decrease of 14% at λ=450nm . The fitted values for the parameters and the CIs were a1=0.995±0.003 , a2=0.91±0.02 , a3=0.012±0.002 , a5=2.4±0.1% , a6=155 ±7% , and a7=10.7±0.6μm . Even though the statistical CIs of this fit are much larger than the statistical CIs of the parameters using the correct model (Sec. 4.1.2), only half of the fit values ( a1 , a3 , and a7 ) are within 2 CIs of the true values. The fitted saturation (a6) is very biased with respect to the true value (155 versus 75%), while the CI is increased only marginally compared to the case where betacarotene was included in the fit. Note that the fitted blood volume fraction a5 has a smaller bias than expected from the figure; this is due to the fitted negative amount of deoxyhemoglobin (which also resulted in the physically impossible 155% saturation). The residual spectrum very clearly shows the features of hemoglobin absorption as well as betacarotene absorption.

Table 7 shows the number of times each of the fitted parameter values is within 2 CIs of their true values, out of the 100 fits and averaged for the 105 synthetic datasets for different concentrations of betacarotene.

Table 7

Average number of times (out of 100) each of the fitted parameter values is within 2 CIs of their true values for different concentrations of missing betacarotene.

0530
a1 96 (4)24 (37)23 (40)
a2 94 (3)48 (44)42 (44)
a3 97 (4)43 (42)38 (44)
a4
a5 95 (3)58 (39)46 (42)
a6 95 (5)10 (29)22 (40)
a7 95 (3)45 (42)39 (41)

Table 8 shows the average bias and CI for each of the parameters for different concentrations of betacarotene, for blood volume fractions a5> 1% .

Table 8

Average bias and CI for each of the parameters for different concentrations of missing betacarotene.

BetacaroteneconcentrationAverage biasAverage CI
a1 [−]5300.0030.020.0010.004
a2 [−]5300.0120.070.0040.02
a3 [−]5300.0030.020.00070.003
a4
a5 [%]5300.060.380.030.14
a6 [%]53084416
a7 [μm] 5300.52.10.21.2

It is observed that the average bias in the parameters is larger than twice the average statistical CI of the parameters, especially for the saturation a6 .

4.3.

Fitting with the Wrong Model: Incorrect Scattering Model

In this section, we investigate the magnitude of the bias in the fitted parameters that is introduced when a wrong scattering function is used in the fitting routine. This is established by fitting with only six free parameters by excluding the Rayleigh scattering term a3(λλ0)4 in Eq. 1 from the fit function. Figures 5a and 5b show a typical fit and the residuals, respectively, for the synthetic spectrum with a3=0.01 , a4=5μM , a5=2% , and a6=75% . The fitted values for the parameters and the CIs were a1=1.0092±0.0009 , a2=1.052 ±0.006 , a4=6.0±0.9μM , a5=1.75±0.05% , a6=73±4% , and a7=9.6±0.3μm .

Fig. 5

(a) and (c) typical fit and (b) and (d) the residuals for the synthetic spectrum with a3=0.01 and 0.1, a4=5μM , a5=2% , and a6=75% , fitted without Rayleigh scattering in the fitting model.

054044_1_034805jbo5.jpg

Figures 5c and 5d show a typical fit and the residuals, respectively, for the synthetic spectrum with a3=0.1 , a4=5μM , a5=2% , and a6=75% . The fitted values for the parameters and the CIs were a1=1.105±0.003 , a2=1.35 ±0.02 , a4=10±3μM , a5=1.2±0.2% , a6=65±21% , and a7=10±2μm . Even though the statistical CIs of these fits are much larger than the statistical CIs of the parameters using the correct model (Sec. 4.1.2), only half the fit values are within 2 CIs of the true values. The blood volume fraction is particularly underestimated with respect to the true value. The residual spectrum very clearly shows the features of hemoglobin absorption as well as a strong peak at the start of the spectrum associated with an incorrect scattering function.

Table 9 shows the number of times each of the fitted parameter values is within 2 CIs of their true values, out of the 100 fits and averaged for 70 synthetic datasets (only the extreme betacarotene concentrations 0 and 30μM were considered) for different Rayleigh amplitudes.

Table 9

Average number of times (out of 100) each of the fitted parameter values is within 2 CIs of their true values for different Rayleigh amplitudes with the Rayleigh term omitted from the fit function.

00.010.1
a1 91 (3)0 (0)0 (0)
a2 91 (3)0 (0)0 (0)
a3
a4 94 (3)62 (43)9 (27)
a5 94 (3)11 (25)1 (4)
a6 95 (3)86 (25)80 (34)
a7 94 (3)62 (42)69 (43)

Table 10 shows the average bias and CI for each of the parameters for different Rayleigh amplitudes. Again, blood volume fractions a5<1% were excluded from the calculation of the averages due to the extremely large bias and CI of the blood parameters ( a5 , a6 , and a7 ) for small blood volume fractions. For the largest Rayleigh amplitude, blood volume fractions a5<2% were excluded from the calculation of the averages due to the large bias and CI of the blood parameters ( a5 a6 , and a7 ) for these blood volume fractions as well.

Table 10

Average bias and CI for each of the parameters for different Rayleigh amplitudes.

RayleighamplitudeAverage biasAverage CI
a1 [−]0.010.10.00980.1030.00060.002
a2 [−]0.010.10.0420.300.0040.02
a3 [−]
a4 [μM] 0.010.121914
a5 [%]0.010.10.231.40.050.2
a6 [%]0.010.131039
a7 [μm] 0.010.10.71.80.41.3

It is observed that the scattering parameters ( a1 and a2 ) as well as the absorber concentrations ( a4 and a5 ) are very biased compared to the calculated statistical CIs, whereas the blood saturation and vessel diameter ( a6 and a7 , respectively) are not affected as much by imperfections in the scattering function.

4.4.

Identification and Minimization of Bias

In the previous sections we have shown that when the fitting model is correct, the fitted parameter values are unbiased and the CIs can be estimated correctly using the reduced chi-square normalized covariance matrix. When the fitting model is incorrect, however, the fitted parameters are biased and the calculated statistical CIs cannot correctly account for this bias. This is expected, since the systematic errors associated with using the wrong model are not fully accounted for by our purely statistical analysis. It is practically impossible to account for systematic errors using statistical methods only, and hence it should be realized that the calculated CIs represent the lower limits of the true CIs of the parameters. In this section we present a method to identify and minimize systematic errors in the fitting routine.

Identification of systematic errors is easily performed based on the shape of the residual spectrum. When the fitting model is correct, the residual spectrum is featureless and shows random noise only [Fig. 1b]. When the fitting model is incorrect, the residual spectrum will always show features of the chromophore absorption coefficients, and in case the wrong scattering function is used, an additional (slow) background variation over the whole wavelength range can be observed [Figs. 3b, 3d, 5b, 5d]. The major difference in the residual spectra of fits with missing absorbers [Figs. 3b, 3d, 4b, 4d] and fits with an incorrect scattering function [Figs. 5b and 5d] is that in the first case, the residuals have a (strong) negative peak near the location of the absorption maximum of the missing absorber, while in the second case all peaks in the residuals are near the absorption maxima of the chromophores used in the fit. In Sec. 4.4.1, we discuss an approach to minimize the systematic errors associated with missing absorbers, and in Sec. 4.4.2, we discuss a method to minimize the systematic errors associated with an incorrect scattering function.

4.4.1.

Minimization of Bias: Missing Absorber

A missing absorber in the model can be identified by a strong negative peak in the residuals (if defined as data minus fit) at a location away from absorption maxima of the chromophores used in the fit. After identification of a missing absorber, one can reduce the systematic error by refitting the spectrum, but now including an additional absorber that peaks near the wavelength of the negative peak in the residuals. Since a Gaussian function is a good first approximation for the shape of an absorption band, we propose to use a Gaussian function for this unidentified absorber, following Eq. 4:

Eq. 4

μaunknown(λ)=a8exp[0.5(λλpeaka9a10)2].
Here λpeak is the location of the negative peak in the residual spectrum and a8 , a9 , and a10 are additional free parameters representing the “concentration,” the peak shift, and the width of the unknown absorber, respectively. The peak shift parameter a9 accounts for the fact that the location of the negative peak in the residual spectrum is not necessarily the exact location of the absorption maximum of the absorber, but will be close to it. Therefore, this parameter should be constrained to small values to force the location of the absorption maximum of the Gaussian to be near the location of the negative peak in the residual spectrum. We have used a9<10nm in our analysis, but did not investigate in detail how the results depend on this constraint.

When this strategy is employed for the case of the missing absorbers in the fitting routine (Sec. 4.2), dramatic improvement in the fit quality is achieved. Table 11 shows the number of times each of the fitted parameter values is within 2 CIs of their true values, out of the 100 fits and averaged for the 105 synthetic datasets for different concentrations of Gaussian absorber. It is observed that on average, nearly 95% of the fitted parameter values is within 2 CIs of their true values, which indicates that the parameters and their CIs are estimated correctly.

Table 11

Average number of times (out of 100) each of the fitted parameter values is within 2 CIs of their true values for different concentrations of Gaussian absorber with an additional absorber in the fit function.

0530
a1 96 (4)93 (4)92 (6)
a2 94 (3)94 (6)94 (7)
a3 97 (4)95 (8)94 (10)
a4
a5 95 (3)91 (8)89 (12)
a6 95 (5)91 (10)89 (11)
a7 95 (3)93 (5)93 (6)

The location of the peak λpeak was determined by the mimimum of the residual spectra, and the absorption maximum of the unknown absorber is λpeaka9 . For the small Gaussian absorber concentration a4=5 , the average fitted values of the missing absorber were λpeaka9=701nm and width a10=24nm . For the larger Gaussian absorber concentration a4=30 , the average fitted values of the missing absorber were λpeaka9=699nm and width a10=27nm . Thus in both cases the absorption maxima and widths of the unknown absorber reproduce the true absorption spectrum very accurately.

Table 12 shows the number of times each of the fitted parameter values is within 2 CIs of their true values, out of the 100 fits and averaged for the 105 synthetic datasets for different concentrations of betacarotene.

Table 12

Average number of times (out of 100) each of the fitted parameter values is within 2 CIs of their true values for different concentrations of betacarotene with an additional absorber in the fit function.

0530
a1 96 (4)72 (16)39 (35)
a2 94 (3)53 (22)11 (20)
a3 97 (4)38 (24)9 (19)
a4
a5 95 (3)42 (21)17 (23)
a6 95 (5)60 (33)55 (33)
a7 95 (3)44 (26)25 (32)

The location of the peak λpeak was determined by the mimimum of the residual spectra, and the absorption maximum of the unknown absorber is λpeaka9 . For the small betacarotene concentration a4=5μM , the average fitted values of the missing absorber were λpeaka9=453nm and width a10=40nm . For the larger betacarotene concentration a4=30μM , the average fitted values of the missing absorber were λpeaka9=455nm and width a10=37nm . These values were determined from the fits with blood volume fractions a5<10% only; for larger blood volume fractions, the minimum of the residual spectra was colocated with the blood absorption peaks. As a result, no missing absorber could be identified for blood volume fractions a510% . Figure 6 shows the normalized betacarotene absorption spectrum and the average extracted absorption spectra of the unknown Gaussian absorber determined from the fits. Good agreement in both the peak position and the width is found between the true betacarotene absorption coefficient and the absorption coefficients of the unknown Gaussians determined from the fits. Despite these promising results, Table 12 shows that on average, much less than 95% of the fitted parameter values is within 2 CIs of their true values. This indicates that there are still systematic errors and that the resulting fit values are biased. The disappointing results are also due to the fact that the method fails to identify a missing absorber for large blood volume fractions. Table 13 shows the average bias ( absolute value of the absolute difference between true and fitted parameter values) and CI for each of the parameters for different concentrations of betacarotene, with an unknown Gaussian absorber in the fit that has a maximum absorption at a wavelength no more than 10nm away from the minimum of the residuals, as explained before. Blood volume fractions a5<1% were excluded from the calculation of the averages, as discussed previously. It is observed that the average bias in the parameters is mostly larger than twice the average statistical CI of the parameters, but much smaller than before (Table 8) when no additional Gaussian absorber was included in the fit. In particular, the bias on the saturation is decreased by a large amount.

Fig. 6

Normalized betacarotene absorption spectrum and the average extracted absorption spectra of the unknown Gaussian absorber determined from the fits. Solid line: true betacarotene spectrum. Dotted line: Gaussian approximation from fits with large betacarotene concentration. Dashed line: Gaussian approximation from fits with small betacarotene concentration.

054044_1_034805jbo6.jpg

Table 13

Average bias and CI for each of the parameters for different concentrations of betacarotene, with an unknown Gaussian absorber in the fit, which has a maximum absorption at a wavelength no more than 10nm away from the minimum of the residuals.

BetacaroteneconcentrationAverage biasAverage CI
a1 [−]5300.0020.0050.0080.001
a2 [−]5300.0040.0110.0040.005
a3 [−]5300.0010.0060.0010.001
a4
a5 [%]5300.040.240.020.05
a6 [%]5301.53.50.61.5
a7 [μm] 5300.31.60.10.3

4.4.2.

Minimization of Bias: Wrong Scattering Model

When a wrong scattering model is used, the residuals show features of the chromophore absorption coefficients on top of an additional background variation over the whole wavelength range. After identification of a wrong scattering model, one can reduce the systematic error by refitting the spectrum, but now including an additional scattering term. We propose to use a second-order polynomial function as an additional scattering term, since it provides a good first-order approximation to any nonperiodic continuous function, such as usually encountered in tissue scattering:

Eq. 5

R(λ)=[a1(λλ0)a2+a3(λλ0)2+a8(λλ0)+a9] exp[0.4μatotal(λ)].
Thus the scattering function is now modeled by Mie scattering plus a second-order polynomial. Whereas the Mie scattering model was based on physical principles where the parameters a1 and a2 are associated with the density and size distribution of scatterers,21 respectively, the new scattering model is primarily designed to describe the scattering function mathematically without associating physical meaning to the parameters. In other words, the scattering function is modified with the aim to reduce the bias in the concentration estimates of the chromophores, at the cost of modeling with a physically uninterpretable scattering function. Table 14 shows the number of times the fitted microvascular parameter values are within 2 CIs of their true values, out of the 100 fits and averaged for 35 synthetic datasets (betacarotene was exluded in this subanalysis for convenience) for different Rayleigh amplitudes. It is observed that on average, the in the case of a small Rayleigh component, nearly 95% of the fitted parameter values is within 2 CIs of their true values, which indicates that the parameters and their CIs are estimated correctly. For the larger Rayleigh amplitude, the accuracy of the results decreases but is still much better than before (Table 9) without the additional second-order polynomial scattering term. As for the Mie and Rayleigh scattering parameters, it should not be expected that they are accurately reproduced, since the mathematical scattering function used in the fit is very different from the scattering function of the synthetic datasets. Finally, Fig. 7a shows the synthetic scattering function
Ssyn(λ)=1(λλ0)1+0.1(λλ0)4,
and a fit of this function using the scattering function of Eq. 5. Figure 7b shows the residuals. Clearly the scattering function of Eq. 5 is a good but not perfect approximation for Ssyn(λ) , leading to the results of Table 14.

Fig. 7

(a) The synthetic scattering function Ssyn(λ) and a fit of this function using the scattering function of Eq. 5, and (b) the residuals.

054044_1_034805jbo7.jpg

Table 14

Average number of times (out of 100) the fitted microvascular parameter values are within 2 CIs of their true values for different Rayleigh amplitudes fitted with an additional second order polynomial scattering function.

00.010.1
a5 96 (3)92 (5)81 (20)
a6 95 (3)92 (3)80 (10)
a7 95 (3)89 (6)75 (24)

5.

Discussion and Conclusion

We have validated a simple method for determining the confidence intervals on fitted parameters derived from modeling optical reflectance spectroscopic measurements using synthetic datasets. The method estimates the parameter confidence intervals as the square roots of the diagonal elements of the covariance matrix, which is obtained by multiplying the inverse of the second derivative matrix of χ2 with respect to its free parameters by χ2v , with v being the degrees of freedom. Most commercial fitting software packages, employing least squares analysis, automatically generate the covariance matrix. CIs are easily extracted after multiplication of the diagonal elements with χ2v to account for the uncertainty in the estimation of the statistical fluctuations and corresponding standard deviations s(i) of the data points. We have shown in Sec. 4.1 that this method yields correct CIs for our nonlinear fitting function, as long as the model used to describe the data is correct. The statistical CIs do not depend on the number of pixels used to bin the data, as long as the optical resolution associated with this number does not exceed the optical resolution of the biological features. Furthermore, the CI of the fitted parameters correctly increases with decreasing signal-to-noise ratio, while in all cases the fitted parameter values are unbiased with respect to their true values. The CI calculated in this way represents the statistical error on the fit parameters and is a lower boundary of the true fit parameter uncertainty. We have also investigated the sensitivity of our parameters to systematic errors associated with incorrect models in results Secs. 4.2, 4.3. Even when there are only small amounts of missing absorber (Sec. 4.2), such that the fit does not deviate from the data by more than 2 to 3% over the entire wavelength range [Figs. 3a, 3b, 4a, 4b], the fitted parameters may be biased by large amounts compared to the calculated statistical CIs. Figures 4a and 4b are particularly instructive: the quality of the fit looks excellent, the residuals are smaller than 2.5% over the entire wavelength range, and yet the bias in the saturation is more than 15%. Thus imperfections in the fitting model introduce a bias in the fitted parameters that greatly exceeds the estimated statistical CIs. Since multicomponent spectral fitting of in vivo optical spectroscopic data is particularly sensitive to these problems, we have also introduced methods to identify and subsequently minimize the bias in the fitted parameters associated with incorrect modeling (Sec. 4.4). Identification of systematic errors is based on the shape of the residual spectrum. When the fitting model is correct, the residual spectrum is featureless and shows random noise only [Fig. 1b]. When the fitting model is incorrect, the residual spectrum always shows features of the included chromophore absorption coefficients, and in case the wrong scattering function is used, an additional background variation over the whole wavelength range can be observed [Figs. 3b, 3d, 5b, 5d]. In the case of a missing absorber, a strong negative peak in the residuals at a location away from absorption maxima of the chromophores used in the fit is observed. The systematic error associated with a missing absorber was reduced by refitting the spectrum, but now including an additional absorber that peaks near the wavelength of the minimum (maximum negative) in the residuals. Dramatic improvement in the fit quality was achieved when using a Gaussian absorber with three free parameters: amplitude, peak shift with respect to minimum of residuals, and peak width. In the case where the missing absorber truly was a Gaussian, the fits were perfect, and nearly 95% of the fitted parameter values was within 2 CIs of their true values, which indicates that the fitted parameters and their CIs were estimated correctly. Furthermore, the absorption maxima and widths of the unknown Gaussian absorber reproduced the true Gaussian absorption spectrum very accurately. In the case where the missing absorber was betacarotene, which has a non-Gaussian absorption spectrum overlapping the absorption spectra of oxy- and deoxyhemoglobin, the results were not as good. For large blood volume fractions, the minimum of the residual spectra was colocated with the blood absorption peaks, and no missing absorber could be identified for large (a510%) blood volume fractions. For smaller blood volume fractions, the extracted absorption spectra of the unknown Gaussian absorber determined from the fits showed good agreement in both the peak position and the width with the true betacarotene absorption coefficient. Despite this good agreement, much less than 95% of the fitted parameter values was within 2 CIs of their true values for both high and low betacarotene concentrations. This indicates that there are still systematic errors and that the resulting fit values are biased. However, the bias in the parameters was shown to be much smaller than without an additional Gaussian absorber included in the fit; the bias on the saturation was decreased by a large amount (from 8 to 1.5% for small betacarotene concentration, and from 44 to 3.5% for large betacarotene concentration).

Palmer 22 have previously used a similar approach in a clinical dataset of malignant and nonmalignant breast tissue. They included a Gaussian absorber with fixed width and peak position in their fits based on the shape of their previously obtained residuals in a large clinical dataset, under the assumption that the spectral shape of the missing absorber is the same across the whole set of tissue samples measured. Interestingly, they found that although the quality of the fits was improved with the addition of the Gaussian absorber, the conclusions regarding the extracted absorber and scattering parameters from the fits were not significantly affected by its addition, and so all subsequent analysis was carried out without the inclusion of this Gaussian absorber. Based on our present research, this result is surprising. Possibly this statement was made with respect to the average of the fitted parameters over the entire dataset. In that case, the biases in the fitted parameters for the individual spectra may be canceled on averaging. Furthermore, the biological variation in the parameters for multiple samples is probably larger than the bias associated with the systematic error introduced by the incomplete fitting model, as we describe next.

The systematic error associated with a wrong scattering model was reduced by refitting the spectrum, but now including a second-order polynomial as an additional scattering term. The new scattering model is primarily designed to describe the scattering function mathematically without associating physical meaning to the scattering parameters, with the aim to reduce the bias in the concentration estimates of the chromophores. The implementation of this polynomial scattering function greatly reduced the bias in the microvascular parameters, although the approach has its limits for large deviations in the model scattering function from the true scattering function. The mathematical description of the “true” scattering function in case of biological (in vivo) tissue actually depends on the measurement geometry and on the tissue itself. In the case of diffuse measurements, the power law behavior of scattering seems to be an adequate description of the scattering function. In the case of nondiffuse measurements, more complicated functions (including periodic elements23, 24, 25) may be necessary to fully describe the scattering function.

Other potential sources of systematic error that are not investigated in detail in this work are: 1. incorrect wavelength calibration of the spectrometer, and 2 incorrect or unknown specific absorption coefficients of the chromophores. Regarding the first point, we observed that a 1-nm offset in the spectrometer calibration typically resulted in a 3% bias in the fitted saturation a6 (data not shown). Since it should always be possible to calibrate the wavelength axis of a spectrometer to within 1nm , we consider this source of error to be minor. Regarding the second point, the specific absorption coefficients of chromophores are usually measured when dissolved in water, ethanol, hexane, or chloroform. Apart from water, these environments are dramatically different from the in vivo environment, and it is not easy to predict how the environment changes the specific absorption coefficient of the chromophore. Spectral shifting as well as spectral broadening may occur, which would also result in biased fit parameters, the amount of bias depending on the amounts of shifting and broadening. Therefore, it is imperitive to use specific chromophore absorption coefficients in the fitting routine that were measured in an evironment that comes as close to the in vivo environment as possible. This is particularly important for component fits of in vivo fluorescence and Raman spectra, whose shapes can be altered even more dramatically than extinction spectra due changes in the microenvironment. Alternatively, any features in the residual spectra obtained from fitting in vivo spectra using ex vivo component spectra could be used to modify these component spectra and determine the true in vivo fluorescence, Raman, or extinction spectra. Moreover, environmental changes (such as pH variations) that are reflected in altered spectral shapes of the components will lead to features in the residual spectra. This provides an opportunity to characterize and monitor the microenvironment by analyzing the shape of the residual spectra.

As can be seen in Table 4, for a noise amplitude typical for our in vivo data, the CI of the fitted parameters is quite small (e.g., the CI for blood saturation a6 is smaller than 1% for blood volume fractions a5> 1% ). This CI represents the statistical error associated with an individual spectrum. It is critical that this is interpreted differently from a standard deviation calculated from averaging multiple measurements. The latter standard deviation represents the biological variation in the parameters, which for most applications will be larger than the statistical uncertainty in the parameters related to the quality of the data. In fact, in our experience the biological variation in the parameters for multiple measurements is even larger than the bias associated with the systematic errors discussed in results Secs. 4.2, 4.3, and for most applications it will be sufficient to specify these biological standard deviations only. However, completely disregarding the statistical CIs can also cause problems. For example, if the blood volume fraction is very low, the statistical uncertainty of the blood related fit values becomes very large (Table 4). Therefore, in the case of averaging multiple measurements, we feel that calculating the weighted average of the parameters, with the statistical CIs as weight factors, is more appropriate than calculating the unweighted average and ignoring the statistical CIs. Furthermore, calculating the statistical errors in the fit parameters facilitates objective assessment of the quality of spectra, and fits with the possibility to reject poor quality spectra or fit values from an in vivo dataset, thereby avoiding unnecessary pollution of the dataset by poor quality spectra.

Table 8 shows that a missing absorber that overlaps with the spectra of oxy- and deoxyhemoglobin can cause a large bias in the extracted saturation. This is particularly true for small concentrations of blood, where the calculated statistical uncertainty of the saturation is much smaller than the bias in the saturation introduced by a missing spectrally overlapping absorber (data not shown). Therefore, to be safe, we have in previous publications restricted the calculations of average blood saturation to spectra with blood volume fractions larger than 1%.11, 12, 13, 14 It is important to note that the saturation CI was calculated for our specific DPS measurement geometry with 400-μm fibers, featuring a path length of 0.4mm . Since the saturation CI will depend on signal attenuation rather than on blood volume fraction, the saturation CI (and potential bias) will be smaller for similar blood volume fractions in the case of DPS with larger fiber diameters, or any other measurement geometry with a longer path length.

In this work, the parameters were not constrained, i.e., all free parameters were allowed to range from minus infinity to plus infinity. However, in reality the parameters are physically constrained to a10 , 0a24 , a30 , a40μM , a50% , 0%a6100% , and a75μm (the diameter of a red blood cell). We did not implement these boundary conditions here to fully explore the sources and magnitudes of the statistical and systematic errors. When fitting an in vivo spectrum, however, it would be appropriate to constrain the parameters by implementing these physical boundary conditions to exclude nonphysical parameter values. However, in that case, care should be taken that the final value of any parameter is not at one of these artificially imposed limits. Boundary conditions on the individual parameters may also be artificially imposed by transformations of the variables. For example, by making the transformation a612+[arctan(a6)]π , the transformed parameter a6 is restricted to 0<a6<1 for <a6< . In this way the fit is performed without boundary conditions on the fit parameters while still imposing physical boundary conditions on the transformed parameters. If the transformation of any fit parameter a is given by af(a) , the statistical error of the transformed parameter is given by σ=σ[df(a)da] .20 For our example, a612+[arctan(a6)]π we find σa6=σa6π(1+a62) , where σa6 is the square root of the sixth diagonal element of the reduced chi-square normalized covariance matrix obtained from the fit.

In conclusion, we show that we can correctly estimate the statistical error on parameters extracted from fits to reflectance spectroscopic data. We demonstrate the use of various methods to identify and subsequently minimize the bias in the fitted parameters associated with systematic errors. Detailed analysis of the residual spectrum is potentially useful not only for minimization of bias, but also for monitoring changes in the in vivo microenvironment.

Acknowledgment

This research has been supported by the Dutch Technology Foundation STW, applied science division of NWO, and the Technology Program of the Ministry of Economic Affairs.

References

1. 

O. R. Scepanovic, K. L. Bechtel, A. S. Haka, W. C. Shih, T. W. Koo, A. J. Berger, and M. S. Feld, “Determination of uncertainty in parameters extracted from single spectroscopic measurements,” J. Biomed. Opt., 12 064012 (2007). https://doi.org/10.1117/1.2815692 1083-3668 Google Scholar

2. 

T. J. Farrell, M. S. Patterson, and B. C. Wilson, “A diffusion theory model of spatially resolved, steady-state diffuse reflectance for the non-invasive determination of tissue optical properties,” Med. Phys., 19 879 –888 (1992). https://doi.org/10.1118/1.596777 0094-2405 Google Scholar

3. 

R. M. P. Doornbos, R. Lang, M. C. Aalders, F. W. Cross, and H. J. C. M. Sterenborg, “The determination of in vivo human tissue optical properties and absolute chromophore concentrations using spatially resolved steady-state diffuse reflectance spectroscopy,” Phys. Med. Biol., 44 967 –981 (1999). https://doi.org/10.1088/0031-9155/44/4/012 0031-9155 Google Scholar

4. 

G. Zonios, L. T. Perelman, V. Backman, R. Manoharan, M. Fitzmaurice, J. Van-Dam, and M. S. Feld, “Diffuse reflectance spectroscopy of human adenomatous colon polyps in vivo,” Appl. Opt., 38 6628 –6637 (1999). https://doi.org/10.1364/AO.38.006628 0003-6935 Google Scholar

5. 

N. Ghosh, S. K. Mohanty, S. K. Majumder, and P. K. Gupta, “Measurement of optical transport properties of normal and malignant human breast tissue,” Appl. Opt., 40 176 –184 (2001). https://doi.org/10.1364/AO.40.000176 0003-6935 Google Scholar

6. 

J. C. Finlay and T. H. Foster, “Hemoglobin oxygen saturations in phantoms and in vivo from measurements of steady-state diffuse reflectance at a single, short source–detector separation,” Med. Phys., 31 1949 –1959 (2004). https://doi.org/10.1118/1.1760188 0094-2405 Google Scholar

7. 

P. Thueler, I. Charvet, F. Bevilacqua, M. St. Ghislain, G. Ory, P. Marquet, P. Meda, B. Vermeulen, and C. Depeursinge, “In vivo endoscopic tissue diagnostics based on spectroscopic absorption, scattering, and phase function properties,” J. Biomed. Opt., 8 (4), 495 –503 (2003). https://doi.org/10.1117/1.1578494 1083-3668 Google Scholar

8. 

G. M. Palmer and N. Ramanujam, “Monte Carlo-based inverse model for calculating tissue optical properties. Part I: theory and validation on synthetic phantoms,” Appl. Opt., 45 1062 –1071 (2006). https://doi.org/10.1364/AO.45.001062 0003-6935 Google Scholar

9. 

T. J. Pfefer, L. S. Matchette, C. L. Bennett, J. A. Gall, J. N. Wilke, A. J. Durkin, and M. N. Ediger, “Reflectance-based determination of optical properties in highly attenuating tissue,” J. Biomed. Opt., 8 (2), 206 –215 (2003). https://doi.org/10.1117/1.1559487 1083-3668 Google Scholar

10. 

R. Reif, O. A’Amar, and I. J. Bigio, “Analytical model of light reflectance for extraction of the optical properties in small volumes of turbid media,” Appl. Opt., 46 7317 –7328 (2007). https://doi.org/10.1364/AO.46.007317 0003-6935 Google Scholar

11. 

A. Amelink, M. P. Bard, S. A. Burgers, and H. J. Sterenborg, “In vivo measurement of the local optical properties of tissue by use of differential path-length spectroscopy,” Opt. Lett., 29 1087 –1089 (2004). https://doi.org/10.1364/OL.29.001087 0146-9592 Google Scholar

12. 

M. P. Bard, A. Amelink, V. Noordhoek Hegt, W. J. Graveland, H. J. Sterenborg, H. C. Hoogsteden, and J. G. Aerts, “Measurement of hypoxia-related parameters in bronchial mucosa by use of optical spectroscopy,” Am. J. Respir. Crit. Care Med., 171 1178 –1184 (2005). 1073-449X Google Scholar

13. 

R. L. van Veen, A. Amelink, M. Menke-Pluymers, C. van der Pol, and H. J. Sterenborg, “Optical biopsy of breast tissue using differential path-length spectroscopy,” Phys. Med. Biol., 50 2573 –2581 (2005). https://doi.org/10.1088/0031-9155/50/11/009 0031-9155 Google Scholar

14. 

A. Amelink, O. P. Kaspers, H. J. Sterenborg, J. E. van der Wal, J. L. Roodenburg, and M. J. Witjes, “Non-invasive measurement of the morphology and physiology of oral mucosa by use of optical spectroscopy,” Oral Oncol., 44 65 –71 (2008). 0964-1955 Google Scholar

15. 

A. Amelink and H. J. C. M. Sterenborg, “Measurement of the local optical properties of turbid media using differential pathlength spectroscopy,” Appl. Opt., 43 3048 –3054 (2004). https://doi.org/10.1364/AO.43.003048 0003-6935 Google Scholar

16. 

O. P. Kaspers, H. J. C. M. Sterenborg, and A. Amelink, “Controlling the optical path length in turbid media using differential path-length spectroscopy: fiber diameter dependence,” Appl. Opt., 47 365 –371 (2008). https://doi.org/10.1364/AO.47.000365 0003-6935 Google Scholar

17. 

S. W. van de Poll, “Raman spectroscopy of atherosclerosis,” 123 Univ. of Leiden, (2003). Google Scholar

19. 

F. James, Minuit—function minimization and error analysis, 1998). Google Scholar

20. 

P. R. Bevington and D. K. Robinson, Data Reduction and Error Analysis for the Physical Sciences, McGraw-Hill, New York (1969). Google Scholar

21. 

J. R. Mourant, T. Fuselier, J. Boyer, T. M. Johnson, and I. J. Bigio, “Prediction and measurements of scattering and absorption over broad wavelength ranges in tissue phantoms,” Appl. Opt., 36 949 –957 (1997). https://doi.org/10.1364/AO.36.000949 0003-6935 Google Scholar

22. 

G. M. Palmer, C. Zhu, T. M. Breslin, F. Xu, K. W. Gilchrist, and N. Ramanujam, “Monte Carlo-based inverse model for calculating tissue optical properties. Part II: Application to breast cancer diagnosis,” Appl. Opt., 45 1072 –1078 (2006). https://doi.org/10.1364/AO.45.001072 0003-6935 Google Scholar

23. 

L. T. Perelman, V. Backman, M. Wallace, G. Zonios, R. Manoharan, A. Nusrat, S. Shields, M. Seiler, C. Lima, T. Hamano, I. Itzkan, J. Van Dam, J. M. Crawford, and M. S. Feld, “Observation of periodic fine structure in reflectance from biological tissue: a new technique for measuring nuclear size distribution,” Phys. Rev. Lett., 80 627 –630 (1998). https://doi.org/10.1103/PhysRevLett.80.627 0031-9007 Google Scholar

24. 

M. B. Wallace, L. T. Perelman, V. Backman, J. M. Crawford, M. Fitzmaurice, M. Seiler, K. Badizadegan, S. J. Shields, I. Itzkan, R. R. Dasari, J. Van Dam, and M. S. Feld, “Endoscopic detection of dysplasia in patients with barrett’s esophagus using light-scattering spectroscopy,” Gastroenterology, 119 677 –682 (2000). https://doi.org/10.1053/gast.2000.16511 0016-5085 Google Scholar

25. 

V. Backman, R. Gurjar, K. Badizadegan, I. Itzkan, R. R. Dasari, L. T. Perelman, and M. S. Feld, “Polarized light scattering spectroscopy for quantitative measurement of epithelial cellular structures in situ,” IEEE J. Sel. Top. Quantum Electron., 5 1019 –1026 (1999). https://doi.org/10.1109/2944.796325 1077-260X Google Scholar
©(2008) Society of Photo-Optical Instrumentation Engineers (SPIE)
Arjen Amelink, Dominic J. Robinson, and Henricus J. C. M. Sterenborg "Confidence intervals on fit parameters derived from optical reflectance spectroscopy measurements," Journal of Biomedical Optics 13(5), 054044 (1 September 2008). https://doi.org/10.1117/1.2982523
Published: 1 September 2008
JOURNAL ARTICLE
14 PAGES


SHARE
Advertisement
Advertisement
Back to Top