1. Introduction

The clinical value of Optical Coherence Tomography (OCT) [1] depends on high imaging speed to provide real time in vivo imaging [2], high spatial resolution to resolve small tissue structures [3], and sufficient contrast to discriminate between those structures. Contrast in OCT images originates from differences in reflectivity of different tissues, which are caused by their variation in (complex) refractive index n. Unfortunately, contrast is limited because for most tissues n only ranges from 1.3 to 1.4. Localized measurement of the attenuation coefficient μ_t can provide additional information, and may increase the clinical potential of OCT by allowing quantitative discrimination between different tissue types.

The attenuation coefficient can be measured from the OCT signal by fitting a model relation to this signal from a region of interest in an OCT image. Currently, two models are available. Widely used [4,5,6] is the so-called single scattering model, which assumes that only light that has been backscattered once contributes to the OCT signal. A model taking into account multiple scattering was introduced by Thrane et al [7] and has recently been used to extract optical properties of atherosclerotic lesions [8] and human skin [9].

Focusing optics in the sample arm suppress the detection of light scattered from outside the focal volume, similar to confocal microscopy. In clinically used probes and catheters, the optical components of the sample arm are fixed. Therefore, for quantitative extraction of μ_t, the confocal properties of the OCT system have to be taken into account, i.e., the change of the OCT signal with increasing distance between the probed location in the tissue and location of the focus [10,11]. We have recently derived a general expression for the confocal axial point spread function (PSF) for single mode fiber (SMF) based OCT systems [12]. The major advantage of this PSF is that it is described by one parameter only, the Rayleigh length, which can easily be determined experimentally.

In this paper we investigate the steps necessary to extract μ_t from OCT images of weakly scattering non-absorbing samples. This method provides a template that can be applied to other ranges of μ_t. In Section 2, we discuss the general principles of non-linear least squares fitting and introduce test statistics to judge the significance of the best fit values. In Section 3 we establish criteria for choosing an appropriate model for the OCT signal, and proceed to choose a model for weakly scattering media using calibrated samples. Section 4 investigates the range of validity of our PSF in scattering media. Section 5 combines these results to extract the attenuation coefficient from calibrated samples, in the clinically more relevant situation of a fixed focus. Section 6 discusses implications and limitations of this study.

2. Curve fitting

Discrimination of different tissues based on differences in their attenuation coefficient μ_t requires its accurate measurement from OCT data. This is done by defining a functional relationship between the OCT signal as a function of depth and μ_t, and then fitting this model to a region of interest in an OCT image. The curve fitting algorithms and statistics used throughout this paper can be found in textbooks [13,14]. Suppose we fit a model f with M adjustable parameters a_j to N data points (x_i ,y_i ±Δy_i). The maximum likelihood estimate of the model parameters a_j is found by minimizing the quantity χ ² given by:

i.e., by minimizing the sum of squared, weighted, residuals. It is appropriate to use the measurement errors for weighting, i.e., σ_i = Δy_i, for experimental data. For non-linear models the χ ²-minimization is an iterative process, implemented by the Levenberg-Marquardt method. The number of degrees of freedom (dof) of the fit is defined as N–M. To judge the significance of the best fit values of the parameters a_j , uncertainty estimates of these values plus some goodness-of-fit statistics have to be calculated. Note that often uniform weighting is used, i.e., each data point is assigned equal weight in the curve fitting (σ_i = 1 in equation 1).

Next to the best estimates, the standard error of each fitted parameter is calculated. If the standard error is small and the parameter is changed a little, the curve will fit much worse (i.e., higher χ ²). The magnitude of the weights σ_i therefore influences the ‘elbowroom’ of the parameter. Consequently, when the standard error is to be used as a reliable estimate of the uncertainty of the fit parameter, weighting with the measurement errors is essential. From the standard errors, 95% confidence intervals (c.i.) are calculated which are more insightful as uncertainty estimates: if the fitting is repeated on another data set from the same sample, the best fit value of the parameter is expected to fall within this c.i. 95 out of a 100 times. So-called parameter dependencies (between 0 and 1) are also calculated: a value (very) close to 1 indicates that the fit does not depend heavily on the parameter, and may point to over-parameterization (i.e. a change in the parameter can be compensated for by changing the other parameters).

Note that a small c.i. can also be found when a best fit does not follow the data very well, for example due to an inappropriate model. The correlation coefficient R² (between 0 and 1) is calculated for each fit. A R² close to 1 indicates the best fit comes close to the data points. Arbitrarily, a R²>0.8 is assumed as reasonable. However, a fit using physically unrealistic parameter values can also have a high R². The χ ²-minimization ‘assumes’ that the weighted residuals have a Gaussian distribution and have the same standard deviation along the best fit curve; and that the deviation of a point from the curve is not correlated to the deviation from the next or previous point. These assumptions are tested using a Shapiro-Francia test and a runs-test. The Shapiro-Francia test computes the correlation coefficient W’ between a normal distribution and the distribution of the weighted residuals, and it is interpreted in the same way as R². The runs-test compares the observed number of runs (a series of consecutive data points either above or below the curve) to the expected number of runs. Again arbitrary, W’ > 0.8 and p-value of the runs test p_runs > 0.05 are considered acceptable. The assumption of independent scatter of the weighted residuals is a priori violated due to low pass filtering (either in hardware or software) in the OCT data acquisition, or possible speckle averaging [6] prior to fitting. When necessary, this is dealt with by using each n’th point of the data set in the fitting, where n corresponds to the number of points in the duration of the averaging time of the filters. If all points are used in the fit, p_runs is disregarded. All curve fitting algorithms and statistical analysis algorithms are implemented in LabVIEW 7. The algorithms were verified against commercially available curve fitting software (Microcal Origin 6.0, Microcal Software, Inc.).

3. Choosing a model for the OCT signal.

The main question in choosing a model for the OCT signal is: can multiple scattering effects be ignored? Since imaging depths generally do not exceed ≈ 1 mm, this may well be justified for weakly scattering media (μ_t < 6 mm^-1 for the samples used in this paper). We first compare the single scattering model [4,5,6] to a multiple scattering model [7] using the same calibrated scattering samples with μ_t ranging from 2 mm^-1 to 6 mm^-1 described in [15] and used in a similar analysis in [8]. The experiments are performed under the condition of dynamic focusing (i.e., the focal plane coincides with the probing location) such that the influence of the confocal properties during the depth scan is constant.

In OCT, interference between the light returning from the sample arm and the reference arm takes place only when the path length difference between both arms is matched within the coherence length of the light source (coherence gating). In the following, the OCT signal i(z) refers to the amplitude of the interference signal, and we define z=0 at the sample interface. Ideally, z is the probing depth in the sample but the term ‘location of the coherence gate in the sample’ is more accurate. In the single scattering model, only light that has been backscattered once contributes to the OCT signal and the OCT signal is given by Beer’s law:

The factor 2 accounts for round trip attenuation, the square root appears because the detector current is proportional to the sample field rather than intensity.

The contribution of multiple scattering to the OCT signal has been described by Thrane et al [7]. Following their terminology, the OCT signal for dynamic focusing is expressed as the root mean square heterodyne signal current:

The last two terms under the square root describe the contribution due to multiple scattering. Here, $w_s^2$/$w_h^2$=(1+[2w ₀/ρ ₀(z)]²) where w_h and w_s are the 1/e intensity radii of the probe beam with and without scattering, respectively, w ₀ is the 1/e intensity radius at the focusing lens and ρ ₀ is the lateral coherence length, given by: ${\rho }_0\left(z\right)=\sqrt{\frac{3}{{\mu }_{s}z}}\left(\frac{{\lambda }_0}{\pi {\theta }_{\mathit{rms}}}\right)\left(\frac{\mathit{nf}}{z}\right),$, where λ ₀ is the center wavelength of the light source,f is the focal length of the objective lens and θ_rms is the root-mean-square scattering angle. Note that μ_s in equation 3 is the scattering coefficient. For non absorbing samples as used in this paper μ_s=μ_t. For the epoxy samples A1-E1, the attenuation coefficient μ_t increases from ≈2 mm^-1 to ≈7 mm^-1; n=1.55. The scattering anisotropy g = 0.75 which means that Eq. (3) is slightly outside of the range of validity of small-angle forward scattering. However, in [8] good agreement between Eq. (3) and experiment was found using the exact same samples, at 1300 nm. θ_rms in Eq. (2) is approximately given by √(2(1-g)) = 0.71.

The SMF (Fibercore, SM750, mode field diameter 5.3 μm.) based OCT setup used in the experiments includes a Ti:Sapphire laser (Femtolasers, λ₀ = 800 nm, Δλ= 125 nm FWHM). Reference mirror and the focusing lens in the sample arm are mounted on two voice-coil translators (QuickScan V-102.2L, Physik Instrumente). Scan speed was 1 A-scan/s. Dynamic range was 111 dB. The detector current is demodulated using a lock-in amplifier and low-pass filtered in software prior to storage. All data acquisition software is written in LabVIEW 6. The collimating lens and focusing lens in the sample arm are both Edmund Optics Achromats P45-793, f=25 mm, NA=0.08. Chromatic aberration expressed as max.-min. effective focal length is 10 μm in the bandwidth of our light source. Depth of focus in air is 126±6 μm (corresponding to 2 × Rayleigh length measured in air). The lateral resolution (determined by the spot size of the focused sample beam) is appropimately 7 μm. We recorded OCT images of epoxy samples A1-E1. OCT images contained 500 A-scans of 4096 points (0.24 μm axial, 20 μm lateral increment). After imaging, a 6x6 pixel moving average filter (approximately corresponding to 1 coherence length) was applied to the data to reduce speckle. Note that this averaging reduced the standard deviation of the data by a factor √6. A larger averaging kernel would lead to further reduction but also to undesirable loss of resolution. The average and standard deviation of 50 A-scans was calculated for use in the curve fitting. Although averaging more A-scans would yield a smoother data set for the fitting, for clinical OCT images we expect to be able to use only the limited number of 50–100 A-scans for a specific tissue region.

Three different models were fit independently to the average A-scans using the standard deviation as weights: (I) the single scattering model of Eq. (2) with an added offset i₀ and multiplier A. For each fit, i₀ was fixed to the average noise level of the data set; A and μ_t were the free running parameters. Model (II) is based on Eq. (3), with added offset i₀ and multiplier A. θ_rms was fixed to its theoretical value and i₀ was fixed to the average noise level of the data set, but was allowed to vary if this improved the fit. A and μ_t were the free running parameters. Model (III) is the same as II, with all parameters i₀, A, μ_t and θ_rms allowed to vary. For all models, convergence of the algorithm is checked by using different initial guess values for the model parameters. In all fits, every 20’th data point was used, corresponding to the response time of the software low-pass filter (see Section 2).

 

Fig. 1. comparison of attenuation coefficients extracted from epoxy samples A1-E1 using integrating sphere measurements (from ref [15]), and curve fitting using models I,II and III. The error bars represent 95% confidence intervals of the extracted fit parameter.

Download Full Size | PDF

The results of the best fit values of μt, using models I, II and III are summarized in Fig. 1. The error bars represent the 95% confidence intervals of the fitted μ_t. The first criterion in choosing a model is whether or not the best fit values of the parameters, and their c.i. and dependency, are physically reasonable. Both model I and II give physically acceptable values for the attenuation coefficient and corresponding confidence intervals. In the fits of model II to samples C1-E1 the offset i₀ was allowed to vary, because it reduced χ ² compared to fitting with fixed i₀. However, this causes larger confidence intervals of μ_t compared to model I because the fit is now not as ‘tight’: within the limits given by the measurement errors, a change in μ_t can be compensated for to a larger extent by a change in A and i₀. The fits of model III did not yield physically possible values for θ_rms (which tended to unrealistically large values, effectively annihilating the multiple scattering contributions). Both the offset i₀ and the multiple scattering contributions in Eq. (3) result in the OCT signal approaching a constant value, with increasing depth. Therefore, variation of i₀ and θ_rms has the same effect. Fixing i₀ to the noise level of the data set did not resolve this problem. We conclude that model III does not model the OCT signal for this range of scattering coefficients, relatively low scattering anisotropy, and our measurement setup very well.

Both model I and II appear appropriate for describing the OCT signal. The goodness-of-fit of both models is judged by their R² value which was 0.8 for sample A1 and larger than 0.95 for samples B1-E1, for both models. Violations of the assumptions of non-linear least squares fitting are checked by judging W’ and the p-value of the runs-test. There is no evidence these models are inappropriate since W’ > 0.86 and p_runs > 0.05 for all samples. Consequently, the model yielding the smallest value of χ ² should be chosen. By definition, variations in χ ² < 1 are not significant [13]. A fit using a model with more parameters (less dof) likely has smaller χ ² and an F-test can assess whether this reduction in χ ² is worth the cost of having less dof. A comparison per sample shows that χ ² is not reduced between model I and II; whereas model II has an additional fit parameter for the larger scattering coefficients. Therefore the simpler model I is chosen.

We conclude that there is no evidence against using the single scattering model I for describing the OCT signal with depth, for μ_t < 6 mm^-1 and our measurement setup, and we will use this model in the remainder of this paper. For other ranges of μ_t the analysis outlined in this section should be repeated to establish the appropriate model.

4. Modeling of the confocal PSF

In clinically used probes and catheters, dynamic focusing is in general not possible, and the influence of the confocal point spread function on the OCT signal has to be accounted for to quantitatively extract attenuation coefficients. We have recently introduced a general expression for the PSF of single mode fiber based OCT systems. In this section, we investigate the range of validity of this expression in scattering media.

The axial confocal PSF for these OCT systems h(z) is given by:

Here, z _cf is the position of the confocal gate and z_R is the ‘apparent’ Rayleigh length used to characterize the PSF. The Rayleigh length z ₀ of a Gaussian beam is given by z ₀ = πnω²/λ₀ with ω the beam waist at the focus and λ0 the center wavelength of the light source. The apparent Rayleigh length is related to z ₀ through z_R = αz ₀ where a is used to distinguish specular reflection (α=1) from diffuse reflection (α=2) [12]. This distinction is based on theoretical grounds assuming single backscattering (or more generally, assuming the beam is not distorted prior to and after backscattering). The Rayleigh length of our system measured on a mirror in air is z ₀ = 63 ± 3 μm.

 

Fig. 2. Schematic illustration of the OCT measurement method. A: OCT images of the scattering samples are taken for different positions z_cf of the focus in the sample (indexed i,j,k). B: From each image, the average A-scan is calculated. C: All average A-scans are combined into a data set, shown as a gray scale image, where the horizontal axis corresponds the focus position in the sample (i.e. to the location of the confocal gate), and the vertical axis corresponds to depth (i.e. to the location of the coherence gate or position of a ‘reflector’).

Download Full Size | PDF

The PSF can be measured by moving a reflector through the focus of the sample beam and recording the detector output. Equivalently, the reflector can be held at a fixed position z while moving the focusing lens. In OCT, we can use the coherence gate to select a ‘reflector’ inside a sample. To systematically evaluate the PSF for specular and diffuse reflections inside scattering media, OCT images of the samples described below were recorded for 100 different positions of the confocal gate z _cf as illustrated in Fig. 2(a). From each image the average A-scan was calculated (Fig. 2(b); red curves represent fits as discussed in Section 5 below) and the average A-scans were combined to a data set shown as a gray scale image in Fig. 2(c), where the horizontal axis corresponds to the position of the confocal gate z _cf and the vertical axis to the position of the coherence gate (or ‘reflector’) z. A similar data set containing the standard deviation of the OCT signal was also constructed for use as weights in the curve fitting. Each row of this data set (constant z) is thus proportional to the PSF at fixed position of a ‘reflector’ inside the sample. Figure 3(a) shows the data set for epoxy sample A1 and the dotted line corresponds to a ‘diffuse reflector’ at depth z = 0.3 mm. The apparent Rayleigh length z _R is then extracted at different depths by fitting eq. 4 to rows of the data set (see Fig 3(b)).

To measure z _R for specular reflection, a mirror was placed at z = 3 mm inside diluted Intralipid solutions with varying μ_t and the PSF at the mirror was extracted from the data sets. The samples were prepared from a stock solution (n = 1.35, changes in n due to dilution are neglected). The attenuation coefficient of the stock solution was determined at μ_t = 4.6 ± 0.2 mm^-1 by dynamic focusing OCT and a fit using model I. Images were recorded for 100 different, fixed positions of the confocal gate (-1.0 to 1.0 mm with respect to sample boundary, measured in air) and 20 A-scans of 2048 points (axial increment 0.73 μm, lateral increment 10 μm). Signal from the mirror is integrated by a 50×1 moving average filter. Equation (4) is expanded with an offset i0 and multiplier A. In the fitting, i₀ was fixed at the average noise level; A, z_R and z_cf were free running parameters. The calculated s.d. of the OCT data is used for weighting.

 

Fig. 3. (a) Gray scale image of recorded data set of sample A1; (b) confocal PSF (dots) along dashed line in A, and best fit (solid); (c) blue circles: z_R for specular reflection and 95% c.i; red squares: average z_R and s.d. for diffuse reflection in different (mfp) intervals. Dashed lines: expected z_R for both cases.

Download Full Size | PDF

Figure 3(c), blue dots, shows the best fit values for the apparent Rayleigh length z_R and 95% c.i. as a function of the expected number of mean free paths 〈mfp〉 = μ_t×z. The value of z_R at 〈mfp〉 = 0 corresponds to a measurement in water. We expect z_R for specular reflection in Intralipid to be n×z₀ = 85 + 4 μm which corresponds well to the experimental data. In all fits, R² > 0.93 and W’ > 0.87. However, values of p_runs were all < 0.05. This is due to small asymmetry of the PSF.

To measure z _R for diffuse reflection, data sets were collected for the calibrated samples A1-E1. OCT images were recorded for 100 different positions of the confocal gate (-0.33 to 1 mm; with respect to sample boundary) and contained 100 A-scans of 4096 samples (axial increment 0.73 μm, lateral increment 40 μm). Equation (4) was fitted to rows of the data sets, i.e. at different depths inside the samples, with an offset i₀ and multiplier A. All parameters are free running in the fitting. To prevent loss of spatial resolution, no speckle averaging was performed. Consequently, the s.d. of the OCT data was too large to serve as weighting factors and uniform weighting was used. As a result, the calculated 95% c.i. are unreliable estimates of the accuracy of the fitting parameter (see Section 2). Fits with unrealistic best fit values for the fit parameters, low R² or W’ were discarded. Consequently, z _R could be measured up to 3.5 〈mfp〉. Data from A1-E1 were combined, and the average z _R and its s.d. for different 〈mfp〉 intervals was calculated. The results are shown in Fig. 3(c), red squares. In all fits, R² > 0.80 and W’ > 0.80. Most values of p_runs were < 0.05 due to small asymmetry of the PSF. Measured values of z_R are larger than expected based on theory (n×2×z₀ = 195 + 9 μm). On average, α= 2.6 + 0.8. About 10% of the difference can be accounted for by the samples being placed under an angle with respect to the probe beam, to avoid high excess noise levels. Oblique incidence of the beam leads to a broader waist in the sample [16], and consequently a larger z _R. Furthermore, in the derivation of Eq. (4) for diffuse reflection, it is assumed that the beam is not distorted by the sample prior to, and after being reflected. Most likely this assumption is not true. The total signal can be split into a contribution due to single (back)scattering and multiple scattering (as is done in [7]). Theoretically, the first contribution can then be described using equation 4 and α=2. This description also suffices if the beam is not distorted due to scattering, i.e., when multiple purely forward scattering occurs. To model any other multiple scattering contribution, an expression for the multiple scattered beam and its waist position should be derived, and the analysis given in appendix A of [12] should be repeated. Such an expression will inherently include the phase function of the scattering particles. From the data of the previous section it seems that in the focus the single scattering assumption is valid. Outside the focus the contribution of multiple scattered light to the OCT signal may be larger, which would lead to broadening of the PSF, i.e., to our larger observed apparent Rayleigh length.

For further analysis, we extracted intensity vs. depth profiles from the datasets corresponding to a fixed distance between the coherence and confocal gate (Fig. 4(a)), and subsequently fitted the single scattering model I to the data. In all fits, R² > 0.8. Prior to fitting a 4×4 moving average filter was used to reduce speckle. The offset in the fitting was fixed to the average noise level. Figure 4(b) shows the fitted attenuation coefficients vs. offset from the focus (all 95% confidence intervals < 0.1 mm^-1; not plotted for clarity). A negative offset from the focus indicates that the coherence gate is located between the lens and the location of the confocal gate. From the samples with higher μ_t, e.g., C1 - E1, we see the fitted attenuation coefficient is indeed less for the out-of-focus situation, indicating a larger contribution of multiply scattered light. For the other samples this effect is less distinct. The asymmetry of the curves of Fig. 4(b) around the zero-offset point explains the slight asymmetry found in the PSF mentioned above. These measurements indicate an increased contribution of multiple scattering in the tails of the PSF and a resulting increase of the apparent Rayleigh length. In any case, the PSF is distinctly broadened for diffuse reflection compared to specular reflection. From the curve fits we find no evidence that the PSF of Eq. (4) should be inappropriate to describe the PSF in scattering media even though α is larger then theoretically predicted. Furthermore, the off-focus fitted attenuation coefficients are within 1 mm^-1 from the corresponding dynamic focusing value for all samples, which suggest that the single scattering model may still yield results of reasonable accuracy.

 

Fig. 4. (a) extraction of data with a fixed distance between coherence and confocal gates. The dashed lines represent zero (1) distance i.e. corresponding to dynamic focusing, and non-zero (negative) distance (2). (b) fitted attenuation coefficient vs distance between confocal and coherence gate for all 5 epoxy samples.

Download Full Size | PDF

We conclude that the combination of confocal and coherence gating to a large extent suppresses the multiple scattered part of the specularly reflected beam at up to 7 〈mfp〉 and of the diffuse reflected beam, e.g., at tissue interfaces, at least up to 3.5 〈mfp〉 and expect that the combination of our ‘single scattering’ point spread function and the single scattering model for the OCT signal can describe the OCT signal recorded in fixed-focus geometry with reasonable accuracy.

5. Extraction of μ_t from fixed-focus OCT data

From the results of Sections 3 and 4, we expect that a model based on the PSF as given by Eq. (4), using α=2, in combination with the single scattering model of Eq. (2) will allow extraction of attenuation coefficients of weakly scattering samples in the fixed-focus geometry with sufficient accuracy. The OCT signal as a function of depth i(z) is then proportional to the product of the PSF and single scattering model:

To verify Eq. (5), it was fitted to columns of the recorded data sets of the epoxy samples A1-E1, described in Section 4, which are the average A-scans at a specific position of the confocal gate z _cf inside the sample (see Fig. 2). Prior to constructing the data set for the present analysis, speckle in each of the OCT images was reduced by a 4×4 moving average filter (corresponding to approximately 1 coherence length). In the curve fitting, Eq. (4) was expanded with an offset i₀ and multiplier A; i₀ was fixed at the average noise level, z _cf was fixed at its pre-set position and z _R was fixed at n×2×z₀ = 195 μm (α=2). A and μ_t are free running parameters in the fit. Standard deviation of the OCT data was used for weighting. We expect this method to be close to what will be possible for clinical images. Figure 5 shows an average A-scan of sample B1 and the best fit to the data (red curve). From this fit, μ_t = 3.92 ± 0.15 mm^-1 which corresponds well to the attenuation coefficient measured using dynamic focusing (μ_t = 3.69 ± 0.37 mm^-1). To further illustrate the distinction between specular and diffuse reflection in h(z), the blue line shows a fit using z _R fixed at 98 μm (α=1) which gives μ_t = 3.09 ± 0.17 mm^-1. As expected, α=2 provides the better fit. Other examples are shown in Fig. 2(b).

 

Fig. 5. Average A-scan of sample B1, with the focus fixed at z=0.3 mm (black line); best fits to the data with equation 5 using α=2 (red line) and α=1 (blue line).

Download Full Size | PDF

The 95% confidence intervals of μ_t are smaller than the c.i. obtained for dynamic focusing in Section 3. In the present data, 50 additional A-scans were used for averaging, which would maximally reduce the s.d. of the OCT data by √50 ≈ 7 for fully developed speckles. Recall that the size of the measurement errors used for weighting determines the sensitivity of χ ² to small changes in the fit parameter, which is expressed in the magnitude of the confidence intervals. Nevertheless, the 95% c.i. may not represent a reliable estimate of the accuracy with which the attenuation coefficient can be determined, because of possible dependence of the fitted μ_t on the focus position.

Figure 6, left panel, shows the best fit values of the attenuation coefficient for the fixed focus geometry vs. the location of the focus with respect to the sample boundary. In all fits, α=2 is used. For all samples, the attenuation coefficient is underestimated when the focus is located near the sample boundary, and overestimated when the focus is located inside the sample. This effect is larger for larger mt. For each fit, R² > 0.8 and W’>0.8 (except for sample A1 which had R² > 0.7). Since all points of an average A-scan were used in the curve fitting, the runs-test did not yield valid results (Section 2), therefore the ‘reasonability’ of the best fit values, R² and W’ were used as the only goodness-of-fit criteria.

 

Fig. 6. Left panel: ‘fixed focus’ attenuation coefficient vs. position of the focus in the sample (with respect to sample boundary). Right panel: ‘fixed focus’ attenuation coefficient, averaged over all focus positions vs. attenuation coefficient determined by dynamic focusing (DF). The line y=x is drawn as a guide to the eye. In all fits, α=2 is used.

Download Full Size | PDF

The accuracy of the μ_t measurements regardless of the focus position is calculated by averaging the fixed focus μ_t over z_cf; the result is shown in the right panel of Fig. 6. For all samples, the s.d. in μ_t was < 0.8 mm^-1. The correspondence with the dynamic focusing μ_t’s shown in Fig. 6 is excellent. Moreover, t-tests showed that the means of the measured distributions of μt were all significantly different from each other (p<0.002). This accuracy is expected to be sufficient for discrimination between tissue structures, e.g., for calcified and lipid lesions of atherosclerotic tissue [0]. The systematic underestimation of μ_t when the focus is located near the sample boundary can partly be caused by the influence of the sample boundary itself, or may point to under-parameterization of our model. However, from the fits no statistical evidence was found that shows Eq. (5) to be inappropriate for describing the OCT signal in the fixed focus geometry, for this range of scattering coefficients and our setup.

We conclude that we can accurately retrieve the attenuation coefficients of weakly scattering samples using the single scattering model in combination with the PSF of Eq. (5).

6. Discussion and conclusion

In this study we investigate the steps necessary for extracting the attenuation coefficient from weakly scattering samples (μ_t < 6 mm^-1) by fitting a model describing the OCT signal to a region of interest in the OCT image. For this range of attenuation coefficients, we compare a single scattering model to a model taking into account multiple scattering using OCT measurements on calibrated scattering samples. The scattering anisotropy of this samples (g=0.75) is lower than that generally found in biological tissues. The multiple scattering model is valid only for small-angle forward scattering, which may explain why the full multiple scattering model III performed worse then, e.g., model I taking into account only single scattering. On the other hand, in Ref. [8] exactly the same samples were used in a similar analysis at 1300 nm; here valid results were indeed obtained. No statistical evidence is found opposing the use of a single scattering model for this μ_t range. For larger μ_t the analysis should be repeated, and this will be the subject of further research. It is important to note that choices in the fitting procedure (e.g., size of speckle reduction kernel, number of A-scans for averaging) influence the confidence intervals of the best fit values. We choose these parameters in accordance with what we expect to be relevant for clinical images; increasing kernel size or number of A-scans would lead to unacceptable loss of resolution or image speed.

In clinical practice, dynamic focusing will not be possible and the influence of the focusing optics on the OCT signal has to be accounted for. We have recently introduced a description of the axial confocal point spread function of single mode fiber based OCT systems. The major advantage of this expression is that it is characterized by one single parameter, the Rayleigh length, which for existing probes and catheters can be easily measured. Therefore, implementation of this PSF in data-analysis algorithms for clinical systems will be straightforward. Our model predicts different apparent Rayleigh lengths, characterized by the parameter α, for specular (α=1) and diffuse (α=2) reflection. This difference is experimentally verified; however, the Rayleigh length for diffuse reflection was larger than expected (α=2.6±0.8). Part of the difference (about 10%) is caused by the fact that the sample is placed under an angle with respect to the beam. Analysis of out-of-focus data indicates an increased contribution of multiply scattered light to the tails of the PSF, which also leads to broadening. A more comprehensive model of the PSF would therefore have to take into account the optical properties of the sample, and more specifically the phase function of the scatterers. This will be subject of further study. We note here that the derivation of Eq. (5) is based on Gaussian optics which is only valid in the paraxial approximation (low NA) and on the assumption that the probe beam is not distorted by the tissue. Other descriptions for the confocal PSF for OCT systems have been given in literature and were discussed in [12]. In most reported PSF’s no distinction is made between specular and diffuse reflection. This distinction is made in the OCT signal analysis of Thrane et al [7], however contrary to our findings, the PSF for diffuse reflection is not broadened according to this theory; i.e., in terms of our model, α=1 for both diffuse and specular reflection.

Combining the single scattering model and this axial PSF allows extraction of the attenuation coefficient of calibrated samples at different fixed positions of the focus inside the scattering medium. The accuracy with which the attenuation coefficient could be determined, regardless of the focus position, is approximately 0.8 mm^-1. The precision of the individual measurements (i.e., defined as 95% confidence intervals in this paper) is much higher and in general less then 10%, comparable to for example [11] and [17].

The PSF and the OCT signal description are independent in our model. This allows easy determination of μt in multi-layered tissue; i.e., in small regions of interest in the OCT image. In models based on, e.g., Eq. (3) this extraction may be less straightforward. In clinical images, we expect to have about 50-100 A-scans (corresponding to a certain tissue) available for averaging. Therefore, in the experiments we use the same number for averaging.

Clinical implications

Preliminary studies indicate we can discriminate between lipid rich and calcified lesions in ex vivo images of human atherosclerotic tissue due to a 3-fold lower attenuation coefficient of the former. This identification was not possible based on gray levels and structural appearance in the OCT image alone [18]. Figure 7 shows a typical example of the samples used in that study. The thick grey line represents the average of 100 A-scans corresponding to the highlighted area in the OCT image in the inset. The arrow marks the known position of the focus. The attenuation coefficient of regions of interest (identified as intima and a lipid rich lesion in histology) was fitted using Eq. (4) and the procedures outlined in this paper. The extracted attenuation coefficients for this sample are indicated in the figure, including 95% confidence intervals of the fitted μ_t.

 

Fig. 7. Average OCT A-scan data (thin gray line) of the region depicted in the OCT image in the inset; and the fitted signal using equation 4 (thick black lines) and corresponding attenuation coefficients for the inimal region and a lipid rich region. The arrow shows the location of the focus in the sample.

Download Full Size | PDF

In conclusion, we have shown that a single scattering model accurately retrieves attenuation coefficients for dynamic focusing OCT data up to 7 mm^-1; we have verified our expression for the confocal PSF at specular and diffuse reflection inside scattering media; and have shown that for attenuation coefficients up to 6 mm^-1, our PSF with a single scattering model can extract μ_t with an accuracy of about 0.8 mm^-1 for the clinically relevant fixed focus geometry. This PSF can easily be implemented for existing clinical OCT systems which may allow quantitative discrimination between different tissues in vivo, and therefore increase the clinical potential of OCT.