Comparison of tissue dispersion measurement techniques based on optical coherence tomography

Abstract. The effects of dispersion on optical coherence tomography (OCT) images have long been documented. The imbalance of spectral broadening, caused by dispersion mismatches in the two arms of the OCT interferometer, can result in significant resolution degradation. Efforts to correct this phenomenon have resulted in improved image quality using various techniques. However, dispersion is also present and varies in tissues. As a result, group velocity dispersion (GVD) can be used to detect changes in tissues and provide useful information for diagnosis. Several methods can be utilized to measure the GVD from OCT images: (i) the degradation of the point spread function (PSF), (ii) the shift (walk-off) between images taken at different wavelengths, (iii) the changes in the second derivative of the spectral phase, as well as two new methods, which do not require a reflector and are applicable in intact tissues, i.e., using (iv) the speckle degradation, and (v) the speckle cross correlation. A systematic, experimental, evaluation of these methods is presented to elucidate the capabilities, the limitations, and the accuracy of each technique when attempting to estimate the GVD in scattering samples. The most precise values were obtained from the estimation of the PSF degradation, whereas using the phase derivative method was only applicable to minimally scattering samples. Speckle broadening appears to be the most robust method for tissue GVD measurements.


Introduction
Dispersion, a result of wavelength-dependent index of refraction variations, causes pulse-width broadening with detrimental effects in many pulsed-laser applications. It is also considered to be one of the major causes of resolution degradation in optical coherence tomography (OCT) imaging and, thus much effort has been expended over the years to develop techniques to counteract its effects. As a result, several dispersion compensation methods have been developed. Conventional methods rely on placing the right amount of dispersion balancing material in one interferometer arm of the OCT setup to counteract the effects of the dispersion in the other arm. Initially, dispersion balancing was achieved by a fused-silica prism pair is inserted with faces contacted and index matched to form a variable-thickness window in the reference arm. However, this method is usually only practical for second-order dispersion. [1][2][3] Gratingbased phase delay scanners can also be used for second-order dispersion compensation. Using the rapid-scanning optical delay, one can adjust dispersion by displacing the diffraction grating from the focal plane of the lens to achieve transformlimited interferogram profiles. 4 Dual optical fiber stretchers can also be used for dispersion compensation while affording some degree of tunability. Using two fiber stretchers made up of different fiber types, an all-fiber tunable dispersion compensator can be implemented with the delay and the dispersion in the two arms of the interferometer can be adjusted independently. However, these approaches require bulky equipment and high voltage so they are not easy to use. 5 Recently, a fiber-stretchingbased dispersion compensator has been combined with a grating-based, time domain OCT system to compensate for second-and third-order dispersion, but the system becomes increasingly complicated. 6 Numerical dispersion compensation, based on the use of the fractional Fourier transform (FT), is also possible providing new perspectives on the nature and role of group-velocity dispersion in Fourier domain OCT. 7 Even after an OCT interferometer is optimized, dispersion differences are still present due to the dispersive properties of the material of tissues that are imaged. Since this dispersion is specific to the sample that is causing the effect, the concept of extracting useful material information from the OCT signal has emerged. 8,9 The idea of using dispersion, as a source of contrast, is not new. For example, there have already been reports, in the literature, of using the dispersion of biomolecules to quantify their concentration. For example, the dispersion of hemoglobin was used to extract the concentration of hemoglobin in intact red blood cell. 10 The relation between dispersion and biochemical composition was further demonstrated using quantitative dispersion microscopy, which has confirmed that the dispersion of live HeLa cells agrees well with the dispersion measured for pure proteins solutions. 11 Variations in the dispersion of different types of normal skin have also been identified in vivo with coherent reflection measurements of different skin types. 12 Given the dramatic changes in cellular biochemistry caused by cancer, 13 which are discernible by other optical techniques such as Raman spectroscopy, 14 it is highly likely that dispersion can also be used as a contrast mechanism in OCT imaging of early cancer and result in more accurate disease diagnosis. However, if dispersion is to be used for quantitative measurements, it is very important to know the capabilities and limitations of the different methods of measuring dispersion with OCT. In this paper, we compare different techniques for estimating the group velocity dispersion (GVD) from OCT images, in order to evaluate their accuracy and applicability to highly scattering media such as tissues.

Theory and Methods
There are several methods described in the literature on how to measure GVD from the OCT signals or images. They depend on resolution degradation, image feature walk-off shift, and phase differences (Fig. 1). In addition, two recently developed methods for dispersion estimation in scattering tissues, without the need of a reflector or a strong distinct scatterer, were also included. 8,15 Each is described further in the following sections along with the mathematical framework for the GVD calculations for each case. Since the index of refraction of the material was needed for some of the calculations, this was measured using the method described by Tearney et al. 16

Resolution Degradation
Dispersion causes degradation of the image point spread function (PSF) and, therefore, of the image resolution. The resolution degradation is evident in the pulse width of the broadened Gaussian resolution envelope t d . If the width of the original envelope t τ and the sample thickness L are known, or can be measured, the GVD can be estimated. The original Gaussian resolution envelope width is related to the pulse FWHM by E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 1 ; 6 3 ; 4 6 8 t τ ¼ t FWHM 1.665 (1) and increases with sample thickness L according to 17,18 E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 2 ; 6 3 ; 4 1 9 t 2 Thus the GVD is given by E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 3 ; 6 3 ; 3 6 6 GVD ¼ Experimentally, the original envelope width t τ can be measured from a free-space portion of the reflector (i.e., not covered by the sample) and the degraded width t d from the reflector peak behind the tissue. The sample thickness is measured as the distance between the sample top surface and the plane of the freespace portion of the reflector.

Walk-Off Shift
When dispersion is present, beams at different wavelengths perceive different path lengths. The result is in an apparent shift in the OCT images taken at different center wavelengths, an effect called "walk-off." The mathematical relationship between the shift Δz and GVD is derived below starting with E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 4 ; 3 2 6 ; 7 5 2 where c is the speed of light, n is the index of refraction, and ω is the optical frequency. The differential walk-off Δz between two spectrally separated independent source spectra as a function of the index of refraction is related to the difference in index of refraction by E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 5 ; 3 2 6 ; 6 6 5 where L is the sample thickness. Thus the GVD can be written as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 6 ; 3 2 6 ; 6 0 4 GVD ¼ Δz cLΔω : where Δz is the differential walk-off in the images acquired from two sources, Δλ is the source bandwidth, λ 0 is the center wavelength, and L is the sample thickness. 18,19 For the experimental verification, a single broad source can be used. After acquisition of the interferometric signal, the spectrum is split into two by multiplication with two shifted Gaussian envelopes. Subsequently, two OCT images are formed, corresponding to the two different spectra, from which the walkoff is measured. This approach results in images with reduced resolution but it assures that the OCT images are acquired at exactly the same location and are, therefore, comparable.

Phase Difference
In spectral interferometry, dispersion can be estimated from the interference spectrum produced by two time-delayed beams. The interferometric signal, for a given time or phase delay, can be expressed as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 8 ; 3 2 6 ; 3 1 7 IðωÞ ¼ jE 0 ðωÞj 2 þ jEðωÞj 2 þ fðωÞ expðiωτÞ þ f Ã ðωÞ expð−iωτÞ: (8) The delay here is due to index of refraction variations and appears as a phase shift between different wavelengths. Hence, the GVD can be estimated from the phase changes of the spectrum of the OCT signals. Note that E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 9 ; 3 2 6 ; 2 2 4 includes the phase information on the spectral phase difference as ΔφðωÞ ¼ arg½fðωÞ. To derive the phase difference ΔφðωÞ, the inverse FT of IðωÞ is performed and the following relationship is obtained: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 0 ; 6 3 ; 7 3 0 Then an FT is applied to the component fðt − τÞ to transfer it back to the spectral domain and the complex amplitude becomes E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 1 ; 6 3 ; 6 4 7 fðωÞ ¼ jE 0 ðωÞjjEðωÞj exp½iΔφðωÞ þ ωτ: The phase of this complex amplitude minus the linear delay part ωτ yields to the spectral phase difference between the two beams. Finally, the phase frequency first derivative yields the group delay E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 2 ; 6 3 ; 5 7 2 GDðωÞ ¼ − ∂ðΔφÞ ∂ðωÞ (12) and, finally, the GVD is given by the second derivative of the spectrum phase as [20][21][22] E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 3 ; 6 3 ; 5 0 6 GVDðωÞ ¼ − Experimentally, the spectrum from a single reflector, positioned below the sample, can be extracted by isolating the signal of that particular peak from the real part of the FT of the interferogram (i.e., multiplying the real A-scan with a Gaussian envelope centered at the peak), and then applying an inverse FT to get the spectrum. The second derivative of that spectrum provides the GVD as a function of wavelength. To compare this solution with the results of the previous methods, the average GVD can be calculated.

Using the Speckle Degradation to Estimate the GVD
The methods described above are very difficult to apply in vivo and are limited only to samples where strong and distinct reflectors are present. A developed method based on speckle width degradation does not require such reflectors and can thus be applied in vivo. This method uses the dispersion induced change in the speckle size to estimate the image PSF broadening and then calculate the GVD. 8 However, speckle variations are difficult to estimate due to the randomness of the speckle signal. To do so, a portion of an OCT image that contains speckle, at depth z from the sample surface, denoted as i s ðzÞ, is related to a similar portion of the OCT image from the sample surface, i.e., i s ð0Þ, by a depth-dependent speckle-degrading impulse response sdfðzÞ such that E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 4 ; 3 2 6 ; 7 3 0 i s ðzÞ ¼ sdfðzÞ Ã i s ð0Þ; (14) where * is the convolution of the two terms and z is the depth, which takes the values of z ¼ d 0 , 2d 0 ; : : : ; L, d 0 being the system resolution. To estimate the impulse response sdfðzÞ, a Wiener-type minimization is used. 23 For that purpose, the following least mean square error function εðzÞ is defined E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 5 ; 3 2 6 ; 6 4 4 where E denotes expectation. Minimizing the error function εðzÞ, using a Wiener deconvolution approach, results in an estimate for sdfðzÞ. In analogy to sdfðzÞ, there exists another impulse response rdfðzÞ of a similar form, which describes the dispersion-induced degradation of the system resolution. For calculating the GVD, it is not necessary to explicitly derive the rdf since only its width is required, which can be estimated from E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 6 ; 3 2 6 ; 5 3 6 where d 0 is the system resolution and d sdf ðzÞ is the width of the sdf at depth z. The result of the convolution of the rdf with the OCT image is a degraded image with resolution width d d ðzÞ given by E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 7 ; 3 2 6 ; 4 4 7 since the convolution of two Gaussians, the psf and the rdf, results also in a Gaussian with a width that is the root mean square (rms) of the widths of the two original functions. Given this width d d , the GVD can be calculated using Eq. (3).

Using the Cross Correlation of Speckle to Estimate the GVD
Another approach to estimate the GVD, without the need for strong distinct reflectors, is the use of the speckle cross correlation to estimate the walk-off shift between OCT images at different wavelengths. The walk-off shift can be estimated from the cross correlation of A-scans from corresponding regions of the two half-spectrum images at different center wavelengths created as described in Sec.

Experimental Methods
For the experimental verification and comparison of the various GVD estimation methods, a swept source OCT system, with a center wavelength of 1300 nm and a resolution of 12 μm in air, was used. Samples of various glasses, a collagen gel (2 g of collagen in 10 ml of water let to solidify and then cut into slices) as well as fresh ex vivo sections of porcine muscle and porcine adipose tissue were used. All the samples were placed over a reflector, which served as a reference for the actual thickness and system resolution measurements, and eight images (5 mm × 4 mm) were acquired from different regions for each type of sample. The glass samples were chosen so that they spanned a wide range of GVD values. Their well-characterized properties were used to verify the validity of the OCT techniques. The glass dimensions were 12.5-mm diameter and 2-mm thickness in all cases. The biological samples were chosen based on their scattering properties, ranging from minimal scattering (collagen) to very scattering (adipose tissue). They were sliced manually so their thickness varied between 1.5 to 2.5 mm. Fortunately, since the actual thickness was measured for each lateral location (see Sec. 4.1), the thickness variation was included in the GVD calculation and did not affect the comparison. In Fig. 3, a list of all the samples with photographs as well as the values of their index of refraction and GVD, obtained from the literature, is provided.
The OCT data were processed in MATLAB TM . Initially, an automated algorithm detected the top surface of the sample as well as the mirror peak location in free-space and below the tissue. The tissue thickness, index of refraction, and mirror reflection Gaussian widths were subsequently calculated assuming a Gaussian envelope shape. From those measurements, the GVD was estimated by combining the values from 250 A-scans in each image. The GVD of each sample was taken as the median of those 250 values. The standard deviation of the GVD values obtained from the eight images of each sample, as well as the mean error between the estimated and measured GVD were used as indicators of the precision and accuracy of each of the different methods described above. Furthermore, the standard deviation of the GVD measurements of single A-scans within an image was used as a measure of each technique's robustness. Given that value, the minimum number of measurements that must be averaged in order to get a GVD estimate with an error E of 10% or less with a confidence level α of 95% was calculated by E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 8 ; 3 2 6 ; 6 3 1 where α ¼ 0.05, E ¼ 0.1, σ is the standard deviation, and Z α∕2 ¼ 1.96 is the critical value of the normal distribution at α∕2. The results are compared in the following sections, to evaluate the precision, accuracy, and robustness of each method. Figure 4 shows typical OCT images of the samples utilized in this study. Using the reflector location (blue lines) as reference, the sample thickness and the group index of refraction, at each lateral location, were calculated using the technique described in the literature, 16 i.e., E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 1 9 ; 3 2 6 ; 4 4 0

Imaging Results
The signal-to-noise-ratio (SNR), defined as the maximum sample intensity divided by the rms noise of the background, was an  Fig. 3 Samples used for GVD estimation with photographs and references values for n and GVD. [24][25][26][27] Journal of Biomedical Optics 046003-4 April 2019 • Vol. 24 (4) average of 70 dB with a 7.4-dB standard deviation for the original images and 68 dB with a 6.5-dB standard deviation for the half-spectrum images. One concern might be that the walk-off shift estimation might be affected by the lower SNR. However, a closer look at the relationship between the error in the GVD estimation and SNR reveals that there is no correlation between the two, at least not in the range of 60 to 80 dB (Fig. 5).   Figure 8 shows an example of the calculation of the GVD using the phase difference method. The real part of the FT of each A-scan was multiplied by a Gaussian envelope located at the location of the reflector behind that tissue thus isolating that single peak [ Fig. 8(b)]. The spectrum of the single peak was obtained from the inverse FT of the single-peak interferogram [ Fig. 8(c)] resulting again in a full spectrum but which now contained only the frequencies corresponding to the single peak. Applying Eq. (13) to this clean spectrum, the GVD was calculated as a function of wavelength [ Fig. 8(d)]. Unwrapping of the phase was performed by changing absolute phase jumps greater than pi to their 2 Ã pi complement. The second derivative was approximated numerically by taking the difference of adjacent values twice. The measurement of the GVD using the phase derivative did not produce accurate results for highly scattering samples. This was due to the presence of strong discontinuities in the phase [Figs. 9(b) and 9(e)], as a result of the scattering discontinuities, producing erroneous GVD estimations [ Fig. 9(f)]. This is consistent with the literature, which predicts minimum and nonminimum phase discontinuities from Mie scatterers. 28    Figure 11 shows an example of the estimation of the GVD from the walk-off shift of half-spectrum images using the cross correlation of corresponding A-scans. Corresponding regions, just above the bottom of the sample [ Fig. 11(a), red Lines], were selected from each half-spectrum image [Figs. 11(b) and 11 (c)]. The cross correlation of corresponding A-scans was calculated and the first peak in the cross correlation was detected. The walk-off shift was estimated from the distance of the peak from  the zero lag location [ Fig. 11(d)] and the GVD was calculated using Eq. (7). The walk-off shift estimation was more robust when there was ample speckle in the images to provide a better cross-correlation approximation. Figure 11(d) shows some typical cross-correlation curves with the arrow pointing to a miscalculation of the walk-off shift due to a weak cross correlation between A-scans. This phenomenon is more common in clear samples such as the collagen gel used in the ex vivo experiments.   Table 1 summarizes the results of the GVD measurements using the techniques described above. The left part of Table 1 lists the results of the GVD calculations over entire samples, estimated from the median of the GVD of 250 A-scans from each image. The accuracy of each technique (i.e., how close the results are to their expected values) for each different type of sample is quantified by the percentage error, the difference from the expected values of Fig. 3. Methods that are inaccurate are highlighted in italics in Table 1. The precision of each technique (i.e., how concentrated the results are around their mean) is described by the standard deviation of the values of the results from complete images (median of 250 A-scans). The right part of Table 1 lists the standard deviations of the GVD values from the  individual A-scans within each image. Larger values imply that the technique is not as robust and that more averages are required to get a good estimate of the GVD of the entire sample. This is also evident from the minimum number of averages required to get an estimate of the GVD with an error of 10%, or less, with 95% confidence, which is listed in the last column of Table 1. Methods that are not robust are listed in bold font in Table 1. Further discussion of the results follows in Sec. 5.

Summary of Results
Another important consideration, when calculating the GVD using OCT, is the size of the fast Fourier transform (FFT) to be used for the reconstruction of the images. In order to get a precise measurement of the degradation of the PSF or walk off, each peak in the OCT A-scan must be adequately sampled to avoid sampling errors. The effect of the FT sampling on the GVD was estimated from the images acquired experimentally, by changing the measurements by AE1 pixel and is shown in Fig. 12(a). These calculations indicate that for a small sampling  Fig. 12(b)].
The most computationally intensive steps in the processing of the data were the required FFTs. To create an image from the raw interferogram, 512 FFTs of length at least 2 16 were required. In the case of the phase difference method, 2 19 FFTs were required in order to accurately reconstruct the spectra from a single peak. The complexity of the speckle degradation method was burdened by an additional 190 two-dimensional (83 × 250) Weiner-type deconvolutions, whereas the final method, speckle cross correlation, required 250 A-scan cross correlations. The average execution time on a PC (i7, quad core, 2 GHz) was ∼1.13, 1.15, 8.79, 1.97, and 1.23 s per image for the psfdegradation, walk-off shift, phase difference, speckle degradation, and speckle cross-correlation methods, respectively.

Discussion
All three mirror-based techniques performed well in lowscattering and fairly uniform samples, such as the glasses or the collagen gel, with the exception of the walk-off shift and the phase difference when interrogating BaF 2 , which is characterized by very low GVD. Based on these observations, only the PSF degradation method can be used to estimate GVD values below 20 fs 2 ∕mm. The measurements from all collagen gel samples were comparable and exhibited little variation (σ ∼ 4% to 6%). In muscle tissue, which is significantly more scattering, only the PSF degradation and walk-off shift methods produced sufficiently consistent results with σ ∼ 12% to 14%, while the phase derivative resulted in σ ∼ 36%. The measurements of the width and walk-off shift were degraded by the presence of speckle noise. Finally, in adipose tissue, which is even more scattering, the PSF degradation and walk-off shift method results were accurate resulting in a σ of ∼11% to 19% while the phase derivative was completely wrong for the reasons described earlier. In general, the methods that did work were sufficiently accurate with an error always <10%. The specklebased techniques were not as precise as the reflector-based methods (σ ∼ 10% to 20%) but they are the only methods applicable to in vivo imaging. With the exception of the speckle walkoff method when applied to low-scattering samples (collagen gel), the remaining estimates were accurate with an error ∼0.3% to 5%. In addition to the methods that result in erroneous estimates (described above and italized in Table 1), there are other methods that are impractical due to the unfeasibly high number of individual GVD values that must be averaged to get a good estimate of the GVD of the sample (bold font in Table 1).

Conclusions
GVD is presented in all tissues and could potentially provide diagnostically useful information. OCT can be used to estimate the GVD and, therefore, enhance the effort for early diagnosis of serious diseases such as cancer. Given the results presented above, for ex vivo GVD estimation, the resolution degradation method is the best choice since it is less sensitive to tissue scattering properties in contrast to the other techniques evaluated. This method performs particularly well for larger GVDs with an error of <0.3%. If the GVD is to be used to provide sensitive diagnostic information from highly scattering human tissues in vivo, it would be preferable to use the speckle degradation as an estimator of GVD. Given that tissue GVD is usually >100 fs 2 ∕mm, this method is expected to perform well with an error of ≤5%. Whichever the case, the use of the GVD as a disease marker is an exciting prospect which should be further investigated.

Disclosure
The authors have no relevant financial interests in this article and no potential conflicts of interest to disclose.