## Abstract

Exploration of nanoscale tissue structures is crucial in understanding biological processes. Although novel optical microscopy methods have been developed to probe cellular features beyond the diffraction limit, nanometer-scale quantification remains still inaccessible for *in situ* tissue. Here we demonstrate that, without actually resolving specific geometrical feature, OCT can be sensitive to tissue structural properties at the nanometer length scale. The statistical mass-density distribution in tissue is quantified by its autocorrelation function modeled by the Whittle-Mateŕn functional family. By measuring the wavelength-dependent backscattering coefficient *μ _{b}(λ)* and the scattering coefficient

*μ*, we introduce a technique called inverse spectroscopic OCT (ISOCT) to quantify the mass-density correlation function. We find that the length scale of sensitivity of ISOCT ranges from ~30 to ~450 nm. Although these sub-diffractional length scales are below the spatial resolution of OCT and therefore

_{s}*not resolvable*, they are nonetheless

*detectable*. The sub-diffractional sensitivity is validated by 1) numerical simulations; 2) tissue phantom studies; and 3)

*ex vivo*colon tissue measurements cross-validated by scanning electron microscopy (SEM). Finally, the 3D imaging capability of ISOCT is demonstrated with

*ex vivo*rat buccal and human colon samples.

©2013 Optical Society of America

## 1. Introduction

Studying nanoscale structural features within live intact biological tissue, while extremely challenging, has profound significance in understanding biological systems. Sub-micron structural changes are known to play a fundamental role in physiology and pathology. For example, the reorganization of the higher-order chromatin structure in epithelial cells has been associated with gene expression and carcinogenesis [1, 2]. Structural change in the extracellular matrix (ECM), which supports and interacts with cells, can also be used as a substantial marker for various types of epithelial cancers [3–5].

Despite the great importance of these nano-scale changes, there is as yet no imaging method that can provide a quantitative measure of sub-diffractional tissue structures in a label-free and non-invasive manner. For centuries, researchers observing biological samples via microscopy were unable to resolve structures smaller than approximately half the wavelength of the illumination, due to a physical limitation called the diffraction limit. During the last three decades, numerous advanced optical techniques have been demonstrated to break the diffraction limit and visualize sub-diffractional structures [6–13]. However, such optics or instrumentation are only applicable to individual cells or thinly sliced tissue rather than intact, *in situ* tissue. Most of these techniques require specific fluorophores that could potentially interfere with endogenous biological processes and only enable imaging of molecular species for which fluorescence dyes exist. Other *in vivo* optical techniques based on the light reflectance [14], such as polarized enhanced backscattering (EBS) [15] and polarization gated spectroscopy (PGS) [16], could be sensitive to nanoscale perturbation, while they are not able to provide a geometrical image of tissue. No existing technique allows for three-dimensional (3D) imaging of tissue with sub-diffractional sensitivity via wholly elastic interactions of light with the tissue.

Optical coherence tomography, on the other hand, is an ideal 3D imaging modality applicable to *in vivo* tissue [17]. The penetration depth is on the order of millimeters while the resolution is on the level of several microns, similar to that of conventional histology. A typical OCT system adopts laser scanning illumination so that the size of the focal spot determines the transverse resolution. The depth resolution is determined by the temporal coherence length of the incidence source, which is a function of center wavelength and the bandwidth [18]. Efforts at improving OCT resolution have pushed the axial resolution down to the sub-micron level using ultra broadband illumination [19], and also improved the transverse resolution using modified illumination methods [20, 21]. However, the resolution limit still persists, and sub-micron structural features remain inaccessible to OCT. Yet the physical reality is that structures smaller than the diffraction limit still affect scattering events. This manifests itself as a phenomenon called “speckle” - appearing as a granular artifact in OCT images. The underlying principle of speckle is that random interference within a resolution voxel such that a dark pixel in the image does not necessarily indicate a weaker-scattering structure [22, 23]. Thus, there may be information revealing the sub-diffractional structural properties encrypted in the speckle signal. The question then becomes: *can one extract such information at nanometer scale from OCT?*

At the fundamental level, biological tissues are composed of complex macromolecules organized in a complicated manner. According to Gladstone-dale’s equation, the mass density of these macromolecules is linearly proportional to the refractive index (R.I.) such that a denser structure gives higher R.I [24]. The complex interconnected structure of tissue, even at nanometer scale, leads one to consider tissue as a continuously varying R.I. medium [25–28]. Given the random and continuous nature of the R.I. distribution, the most comprehensive way to quantify such a medium is by its autocorrelation functions (R.I. correlation function). Recently, we have demonstrated a technique called inverse spectroscopic optical coherence tomography (ISOCT) to quantify the tissue R.I. correlation function [29]. The forward model under the first-order Born approximation was established, and the inverse problem was solved by measuring optical scattering properties.

Here we show that the length scale of sensitivity of ISOCT, particularly the power-law exponent of the backscattering spectrum, is from ~30 to ~450nm. We first briefly introduce the principle of recovering the R.I. correlation function by ISOCT. Then numerical and phantom experiments are presented to verify sensitivity. We report measurements taken from an *ex vivo* colon biopsy samples, and deduce the mass-density correlation functions from the epithelial and stoma compartments. At the same time, SEM images were taken, and the nano-scale image correlation functions were calculated for comparison. We found that ISOCT can indeed quantify the sub-micron structural variations from these two distinct tissue types.

## 2. Inverse spectroscopic optical coherence tomography (ISOCT)

The principle of ISOCT is illustrated in Fig. 1. Conventional Fourier domain OCT (FDOCT) performs an inverse Fourier transform (IFT) on the entire bandwidth of the interference spectrum to an A-scan depth profile. Due to the limited bandwidth and resolution, the A-line profile misinterprets the actual R.I. fluctuation (ΔR.I.), and thus the A-line correlation function contains significant error vs. the actual one. By way of comparison, ISOCT measures the optical properties, including the backscattering coefficient spectrum *µ _{b}(λ)* and the reflection ratio

*α*(defined as the ratio of backscattering and scattering coefficients,

*α = μ*), to much more accurately deduce the R.I. correlation function [29]. The advantage of using OCT is that the spectrum can be extracted from a localized volume sized around a 10 µm cube, where the scattering problem can be largely reduced to an integration of single scatterings from all the scatterers inside the cube. Thus, the inverse problem dramatically simplified under the first-order Born approximation. Another advantage is that the 3D capability of OCT provides guidance to investigate a particular region of interest.

_{b}/ μ_{s}#### 2.1 Whittle-Mateŕn correlation function

We modeled the tissue R.I. correlation function by the Whittle-Mateŕn (W-M) functional family with three physical parameters: normalization factor *N _{c}*, length scale

*l*, and a unitless “

_{c}*functional factor*”

*D*defining the function type [28]:

_{v}(.) denotes the modified Bessel function of the second kind, and

*ρ*is the spatial displacement between any two points inside tissue. An advantage of using the Whittle-Mateŕn family is that, through a combination of the three parameters, it covers most of the realistic types of correlation functions, from inverse power-law functions to exponentials to Gaussian. Examples of the different functional forms encompassed by the model are shown in Fig. 2(a). When 0 <

*D*< 3, the tissue is organized as a mass fractal, and the correlation function has the form of an inverse power law for

*ρ*<

*l*. In this case,

_{c}*D*equals the mass fractal dimension ranging from 0 to 3 in Euclidean space. When

*D*= 3, the corresponding spectral density shares the identical form with the often-used Henyey Greenstein phase function. When 3 <

*D*< 4, the correlation function forms a stretched exponential. When

*D*= 4, the correlation function evolves into an exponential. It approaches a Gaussian function when

*D*approaches infinity. Qualitatively speaking, a smaller

*D*value indicates a more drastic change of R.I. at a smaller length scale, or in other words, a smaller length scale of tissue heterogeneity.

The parameter *l _{c}* determines the physical length scale of the correlation function, beyond which the function drops off quickly, as shown in Fig. 2(a). This parameter can be conceptualized as the cut-off length scale of the correlation function. The precise physical meaning of

*l*depends on

_{c}*D*. For example, for

*D*< 3,

*l*defines the upper range for which

_{c}*C*(

_{n}*ρ*) is a power law with the fractal dimension

*D*: it defines the upper extent of the mass fractal range. For

*D*= 4,

*C*(

_{n}*ρ*) is an exponential function and

*l*is the correlation length, where

_{c}*C*(

_{n}*l*) =

_{c}*C*(0)

_{n}*e*

^{−1}. An extreme case is when

*l*approaches 0, and the correlation function is compressed into a delta function acting as a dipole scatterer for which the scattering is fairly isotropic. On the other hand, biological tissue usually satisfies

_{c}*kl*>> 1 [30], so that the scattering is highly anisotropic and forward directed.

_{c}*k*=

*2π/λ*is the wave number. The shape of the Whittle-Mateŕn function can be fully described by

*D*and

*l*is a scaling factor of the correlation function quantifying the variance of the R.I. fluctuation. A definition of an R.I. correlation function requires a normalization such that at the origin (

_{c}. N_{c}*ρ*= 0), the function should be equal to the variance of the R.I. fluctuation

*σ*. The Whittle-Mateŕn correlation function approaches infinity at the origin when

_{n}^{2}*D*≤ 3. When

*D*> 3, this function has finite value and

*N*can be analytically associated with

_{c}*σ*[28]. To establish a relationship between the normalization factor

_{n}^{2}*N*and

_{c}*σ*, one option is to perform a conditional normalization:

_{n}^{2}*ρ*is the lower length scale of sensitivity. Since a system is no longer sensitive to structures smaller than

_{l}*ρ*, one may assume that the value of the correlation function remains flat when at

_{l}*ρ*.

_{≤}ρ_{l}We note that quantifying the R.I. correlation function enables us to study the statistical tissue mass-density distribution. Barer *et al*. established a linear relationship between the R.I. and local molecular mass-density [24]:

*ε*is the local concentration of the solid material (e.g., macromolecules) and

*χ*is the R.I. increment.

*n*is the R.I. of water. It has been experimentally verified that this relationship holds in the physiological range, and the value of

_{0}*χ*is ~0.18cm

^{3}/g with very little variation among different substances and up to at least 50% concentrations [31]. This linear equation allows one to investigate the mass-density properties by quantifying the R.I. correlation function.

#### 2.2 W-M correlation function and Henyey Greenstein phase function

Since its publication in 1941, the Henyey Greenstein (H-G) phase function [32] has been the primary theoretical tool to quantify the anisotropic factor *g*. Although it was empirically introduced, its accuracy has been verified in the biophotonics community over decades, and it is still being widely used for the inverse problem in diffuse reflectance, integrating sphere, goniometer measurements, *etc*. As described above, when *D* equals 3 (corresponding to a limiting case of a mass fractal structure that fills the entire space and, at the same time, a transition from a mass fractal regime to a stretched exponential correlation function), the corresponding spectral power density of the W-M correlation function shares an identical form with the H-G phase function, with $k{l}_{c}=\sqrt{g}/(1-g)$*(see the detailed derivation in Appendix 1)*. This provides a physical explanation for the H-G phase function.

#### 2.3 Forward model for ISOCT

Under the Born approximation, the spectral power density of a medium can be obtained via Fourier transformation, and the optical scattering properties can be derived analytically [28, 29]. The backscattering coefficient *μ _{b}* is equal to 4π times the backscattering cross section per unit volume [29], and can be written:

*kl*>> 1,

_{c}*μ*spectrum exhibits a power law with its exponent (often called the scattering power, SP) equal to 4

_{b}*-D*,

*kl*<< 1, the medium shows Rayleigh scattering since

_{c}*l*is much smaller than the illumination wavelength, and

_{c}*μ*has the wavelength dependence of

_{b}*k*

^{4}(

*λ*

^{−4}). Figure 2(b) shows the spectrum of

*μ*at various

_{b}*D*values. Since the scattering at a particular wavelength depends on the length scale of the scatterers, we normalized the wavelength by

*l*so that

_{c}*μ*is a function of

_{b}*λ*/

*l*, i.e. (

_{c}*kl*)

_{c}^{−1}. For showing the power law dependence, Fig. 2(b) plotted in log-log scale with the x-axis labeling

*kl*in an ascending order. The y-axis labeled

_{c}*μ*normalized by

_{b}*k*and

*N*. We see that in the Rayleigh regime (

_{c}*kl*<< 1),

_{c}*μ*spectrum is independent of

_{b}*D*since they all behave as dipole scatterers; while when

*kl*>>1,

_{c}*μ*spectrum shows power law dependence (linear curve in log-log plot) with different

_{b}*D*values.

Also the expression of α can be given as:

*kl*>>1 and

_{c}*D*> 2,

*D*values. When

*kl*<< 1 in the Rayleigh regime, the scattering has the dipole scattering pattern whose backscattering ratio α is 1.5. When

_{c}*kl*>>1, the value of α decreases as

_{c}*D*becomes larger. The physical interpretation of a larger

*D*is that more large structures are present in the medium, resulting in more forward-directed scattering.

#### 2.4 Inverse methods for ISOCT

OCT detects backscattering signals. Assuming a homogeneous medium, the A-line OCT signal in depth can be expressed based on Beer’s law:

where*I*is the illumination intensity;

_{0}*r*is the reflectance of the reference arm;

*n*is the refractive index of tissue (~1.38); and

*z*is the penetration depth;

*L*is the coherence length of the light source. The maximum intensity is proportional to

*μ*and the decay rate along depth is proportional to

_{b}*μ*. By taking the natural log, the above equation becomes a linear function:

_{s}*μ*can be estimated by least square (LS) fitting the A-line signal. The microscopic heterogeneity of

_{s}*μ*(

_{b}*z*) can be recovered with the fitted

*μ*:and thus the reflection ratio α(z). The exponential term compensates the illumination attenuation along the penetration depth.

_{s}The next step is to take the Short-time Fourier transform (STFT) and extract *μ _{b}* spectrum. At each wave band selected by a scanning Gaussian window, the process in step 1 is iterated to obtain the

*μ*(z) spectrum. We empirically chose the window width to be

_{b}*k*= 0.18μm

_{w}^{−1}(~15nm at 710nm) granting ten independent points on the spectrum. By another LS fit on the spectrum, the exponent of the spectrum and thus

*D*(

*z*) is obtained according to Eq. (5).

The final step is to calculate *l _{c}*(

*z*) based on

*D*(

*z*) and α(

*z*) based on Eq. (6) or (7). By combining the B-scan transverse scan, the 2D tomographic

*D*and

*l*map can be obtained.

_{c}#### 2.5 Experimental setting for ISOCT

The experimental setup was described in a previous publication [29]. Briefly, a Fourier-domain OCT was adopted and an illumination band from 650 to 800nm was generated by a supercontinuum source. We used polystyrene-microsphere phantoms to calibrate the measurements. As reported in our previous work, within a 90% confidence interval, the uncertainties were measured as ± 0.05mm^{−1} and ± 0.22mm^{−1} for *µ _{b}* and

*μ*, respectively, and ± 0.2 for

_{s}*D*[29]. For an accurate

*D*measurement on biological tissue, a 1:5 diluted 0.08μm polystyrene-microsphere solution (Thermo Scientific, 10% solids by weight) was used to normalize the OCT spectrum due to its Rayleigh scattering spectrum of

*k*

^{4}.

## 3. Length scale of sensitivity of ISOCT

At different length scales, the R.I. correlation function may change, e.g., the R.I. distribution within mitochondria may not apply to the entire cell. It is then important to provide a range of length scales to which ISOCT is sensitive. We first conducted a numerical simulation to systematically perturb the length scale of the correlation function, which leads to a detectable change in the measurement of *D*. Then liquid phantoms composed of microspheres were measured by ISOCT to test the length scale of sensitivity. Lastly, *ex vivo* colon biopsies were measured with SEM verification to test whether submicron discrepancies in different tissue compartments could be detected by ISOCT.

#### 3.1 Numerical perturbation on W-M correlation function

In order to study the lower length scale of sensitivity *r _{min}*, we used a Gaussian smoothing window

*G*to filter out the fine structural details in the R.I. fluctuations, similar to an imaging process with a finite spatial resolution. The changes in D measurements were recorded and the minimum detectable change in

*D*defined the length scale of sensitivity in ISOCT.

Applying a smoothing window in the spatial domain was equivalent to convolving the R.I. correlation function with the autocorrelation function of that window *(see the detailed derivation in Appendix 2)*:

*ACF*(

*f*) stands for the autocorrelation function of an arbitrary function

*f*. We configured the original

*C*with

_{n}*D*= 2.8 and

*l*= 1μm, which are realistic values according to our previous measurements [29]. The FWHM of the Gaussian window

_{c}*G*was set to be

*r*[14]

_{min}*For the upper length scale of sensitivity*

_{.}*r*, we retained fine structures by subtracting the above filtered R.I. fluctuation with the original. The equivalent correlation function is:where

_{max}*δ*is a delta function in the spatial domain. The spectrum of

*µ*from 650 to 800nm was then numerically generated from

_{b}*C*

_{n}’,*D*was measured according to Eq. (4).

The perturbed correlation function *C _{n}’* at small length scales leveled off below

*r*as in Fig. 4(a) when the finer structural detail was eliminated. In Fig. 4(b), the larger length scale correlation dropped to 0 much faster than the original. Examples of perturbed R.I. fluctuation was shown in [14]. Given the perturbed correlation function, the relative change of

_{min}*D*values in terms of the original is shown in Figs. 4(c) and 4(d) with systematically increasing

*r*and

_{min}*r*. Within the given wavelength range, we found that a change in

_{max}*D*could only be detected at length scales larger than

*r*= 35nm and smaller than

_{min}*r*= 450nm, given a 5% sensitivity threshold. This shows that ISOCT is primarily sensitive to the shape of the R.I. (and therefore the mass-density correlation function) for length scales between

_{max}*r*= 35nm and

_{min}*r*= 450nm. The measured

_{max}*D*value increased with increasing

*r*since removing finer structures equivalently manifested as a larger length-scale medium. On the other hand, with the

_{min}*r*perturbation,

_{max}*D*measurements decreased since the medium appeared to possess more small structures.

*One should note that, although the length scale of sensitivity was from ~35 to ~450nm, the calculated l*

_{c}could extend to larger length scales to match the measured optical properties as if the correlation functional form maintains the same form when ρ>450nm.#### 3.2 Phantom study

We next conducted a phantom study to test the length-scale sensitivity of measuring *D*. Aqueous phantoms were composed of a suspension of polystyrene microspheres in a series of diameters ranging from 30nm to 4.3μm as shown in Fig. 5(a). The volume fraction of the microspheres followed an inverse power law with respect to their microsphere diameters (power = −1.2), so that the correlation function was approximated to that of a W-M function with *D* = 2.8 assuming infinite numbers of sphere diameters [33]. The scattering coefficient of the phantoms was maintained at the physiologically relevant value of 9 mm^{−1}.

Two groups of experiments were performed. The first group was intended to determine the sensitivity of ISOCT to subdiffraction-size microspheres and the second group was intended to determine the sensitivity of ISOCT to microspheres of size the order of one wavelength and larger. For each group of experiments, ISOCT yielded the D value. Each *D* value was averaged over 100 A-scans and taken from the surface of each phantom. For comparison, we also used Mie theory to calculate the *µ _{b}* spectrum and further

*D*values for each phantom.

In order to visually compare the difference within phantoms, pseudo-color tomographic *D* maps were created in HSV (Hue-Saturation-Value) color space. *D* encoded Hue and the image intensity encoded Saturation and Value. Ten frames of B-scan images were averaged to overcome the Brownian motion. Moreover, each wavelength-dependent B-scan image was smoothed by a moving window 32µm x 30µm in z (depth) and x (transverse), and then tomographic *D* maps were calculated.

Consider first the methodology and results of the first group of experiments. Here, a series of experiments wherein individual phantoms were constructed by successively excluding the smallest microspheres included in the previous phantom. Specifically, initial phantom (No. 1-1) was comprised of a suspension with microsphere diameters from 30nm to 1μm. The next phantom (No. 1-2) excluded the 30nm microspheres, and was comprised of a suspension with microsphere diameters from 40nm to 1μm. The third phantom (No. 1-3) excluded both the 30 and 40nm microspheres, and was comprised of a suspension with microsphere diameters from 60nm to 1μm. Additional phantoms were constructed with successively larger microspheres excluded.

Phantoms No. 1-1 and 1-2 exhibited the measured *D* value of 1.57 ± 0.04 and 1.51 ± 0.05 respectively, i.e. no observable change within the experimental uncertainty of 0.2. However, phantom No. 1-3 exhibited a measured *D*-value of 1.78 ± 0.06, which was a detectable perturbation relative to phantoms No. 1-1 and 1-2 above the experimental uncertainty. Successive phantoms continued the trend of increasing *D*-value perturbations. Figure 5(b) plots the measured *D* changes for difference phantoms. The x-axis labels the largest particle size excluded in each phantom. As noted above, the change in D first became appreciable, i.e., above the experimental uncertainty, for phantom No.1-3. The difference became more significant with larger particles excluded. Without the small particles, the R.I. correlation function was flattened only at small length scales, verifying our above numerical simulation.

Consider next the methodology and results of the second group of experiments. Here, a series of experiments wherein individual phantoms were constructed by successively adding larger microspheres to the previous phantom. Specifically, the initial phantom (No. 2-1) was comprised of a suspension with microsphere diameters from 0.08 to 0.36μm. The next phantom (No. 2-2) was comprised of a suspension with microsphere diameters from 0.08 to 0.65μm. The third phantom (No. 2-3) was comprised of a suspension with microsphere diameters from 0.08 to 1μm, and so on. Additional phantoms were constructed with successively larger microspheres included.

Phantom No. 2-1 and 2-2 measured the initial *D* value of 2.78 ± 0.08 and 2.42 ± 0.06, which was detectable given the experimental uncertainty of 0.2. When spheres larger than about 2.1µm were added, the *D* value perturbation reached a plateau which indicates that those large spheres did not further alter the SP due to the small volume fraction, as in Fig. 5(c). We should note that the second group of experiments differs from the numerical simulation for estimating the upper length scale of sensitivity [Fig. 4(d)], in which we only eliminated larger length scale while retaining the fine structures. For a sphere, the R.I. correlation function is close to a triangular function that covers length scale from 0 up to the diameter of the sphere [32]. When large particles were added to the solution, the entire correlation function including the small length scales were altered. Further, in our phantom there were only finite numbers of sphere sizes. The large sphere added to the phantom introduced strong oscillatory spectrum which resulted in the decrease in *D* value in Fig. 5(c) instead of an increase as in Fig. 4(d). Nonetheless, the data in Fig. 5(c) still illustrates that *D* can only be sensitive up to a certain upper length scale, above which the value of the correlation function is small and the scattering contribution is negligible.

In tomographic OCT images [Fig. 6(a)] for the first series of phantom study, no apparent difference can be observed across phantoms; while the increment of *D* can be visualized from phantom No. 1-3 and become more dominant in phantom No. 1-4. Although we observed variation within the D maps, the overall difference was distinguishable [Fig. 6(b)]. Similarly for the second series of phantom study, the decrease of *D* values could be observed by the tomographic *D* maps while gray scale OCT images failed to discriminate different phantoms [Figs. 6(c) and 6(d)].

#### 3.3 Experimental verification on biological tissue ex vivo

To demonstrate that ISOCT can be sensitive to tissue structural differences at a sub-diffractional regime, we performed ISOCT measurement on an *ex vivo* human colon sample and compared with the high-resolution SEM images. The structure of colon mucosa consists of a folded epithelial cell layer (crypts) surrounded by the lamina propria (extracellular matrix surrounding crypts), as shown in Fig. 7(a). An SEM image shows distinct bundle structures of collagen in lamina propria (LP), while in the epithelial (Epi) layer, the structures appear more random and fine-grained. We calculated the image autocorrelation function on two tissue components, and Epi has a steeper decline in the correlation function in the sub-micron region than LP [Fig. 7(b)]. Meanwhile, thanks to the capability of localizing the OCT spectrum, *D* and *l _{c}* can be quantified along the penetration depth [Fig. 7(c)]. The correlation functional forms, which were measured by ISOCT, show a similar result as SEM image analysis, i.e., the Epi has a sharper correlation function than LP [Fig. 7(d)]. Both

*D*values in Epi and LP are mostly between 2 and 3 indicating that they are in a mass fractal regime, and it can be seen in Fig. 7(d) that the R.I. correlation function exhibits a power law at sub-micron length scales [30].

The high-resolution SEM image analysis qualitatively verified the structural difference between Epi and LP measured by ISOCT. However, the calculated image correlation functions should not match those of ISOCT. In SEM, the structure was imaged into a two-dimensional (2D) space instead of 3D as in ISOCT, and the intensity of the pixel was not necessarily linearly proportional to the R.I. Furthermore, the sample preparation altered the R.I. distribution from the *ex vivo* intact sample, e.g., fixation changed the R.I. of the tissue composition. Yet, its validity for providing a qualitative structural comparison between the two tissue components remains.

## 4. Three dimensional (3D) ISOCT imaging of biological tissue

As an illustration of 3D ISOCT imaging capability, Fig. 8 shows an example of mapping *D* on an *ex vivo* rat buccal sample. There are three dominant layers: keratinized epithelium (KE), stratified epithelium (SE) and sub-mucosa (SM) [Fig. 8(a)]. The intensity discontinuity at layer conjunctions allows us to segment them at their boundaries so that the *D* from each layer can be measured separately. Specifically, we smoothed the A-line signal by 15µm moving window and calculated the derivative. The first positive peak of the derivative marked the surface of the sample, the second dominant positive peak marked the surface of SM. We then smoothed the surface curves to remove noises. We numerically flattened the sample surface by shifting the A-lines, and averaged the image across transverse direction. By taking a derivative, the averaged thickness of KE was obtained by the distance between the first positive and the first negative peaks (boundary of KE and SE). We then marked the surface of SE by shifting the sample surface by the averaged thickness of KE. The wavelength-dependent attenuation from KE and SE was compensated before the *D* calculation in SM as in Eq. (10) [34]. The resulting *D* map was plotted in HSV space, where *D* encoded Hue and the image intensity encoded Saturation and Value [Fig. 8(b)]. The *D* value histograms at each layer are shown in [Figs. 8(c)-8(e)]. The averaged *D* ± standard deviation from KE, SE and SM were calculated as 2.12 ± 0.59, 2.73 ± 0.49 and 4.82 ± 1.6, respectively. A similar increase of *D* in SM can be observed in sub-mucosa which also contains a large portion of extracellular matrix such as the lamina propria in the above colon sample.

With 3D imaging capability, we were able to detect structural alterations in a region of interest (ROI) demonstrated in Fig. 8(f) with *ex vivo* human colon biopsies. The 3D pseudocolored *D* map and the gray scale morphology were volume-rendered and overlaid concurrently. The cross sectional B-scan *D* map is shown in the movie (Media 1).

## 5. Discussion and conclusions

We have demonstrated results indicating that our spectroscopic analysis by ISOCT is sensitive to R.I., and thus the mass-density correlation function, at length scales from ~35 to ~450nm. The forward model is based on the W-M correlation function, which provides excellent flexibility covering a wide range of functional types. We found that, overall, the value of *D* is ~3 for biological tissues [29]. This verified the Henyey-Greenstein phase function, which shares the identical form of the power spectral density with the W-M function at *D* = 3.

The physical principle of nanoscale sensitivity in ISOCT is straightforward. According to the Wiener-Khinchin theorem, the power spectral density is the Fourier transform of the R.I. correlation function over the entire spatial regime including sub-diffractional length scales. ISOCT essentially collects and analyzes a finite range of power spectral density (backward direction, finite bandwidth. Although the resolving power is limited by the finite bandwidth, the sub-diffractional information is encrypted in the spectrum, e.g., by the *k*^{4} spectral dependence. It is well known that when the size of the scatterers decreases and reaches the Rayleigh regime, the scattering cross-section is proportional to *k*^{4} [35]. When the tissue structure contains more Rayleigh scatterers, the SP is closer to 4, which yields a smaller *D*, if quantified by ISOCT. When the size of the particles increases, the spectrum becomes flatter, leading to a larger *D*. Thus, the way that ISOCT analyzes the spectrum, i.e., by quantifying the SP, is conceptually quantifying a relative scattering weight of the sub-diffractional structures. Higher *D* values indicate a relatively low power ratio of the sub-diffractional scatterers. The introduction of the W-M correlation function further elucidates the physical reality and grants us a quantification method.

Our model treats tissue as a continuous R.I. fluctuation medium having fluctuation. Another useful model of tissue is the discrete particle model - one of the particle shapes is the homogeneous sphere whose scattering field can be analytically solved by Mie theory [35]. Light scattering by some organelles such as mitochondria and nuclei have been successfully modeled by Mie theory [11, 36, 37]. However, there are several drawbacks to the discrete particle model. The cellular structures are mostly interconnected rather than isolated or discrete. Also there are often inner structures or building blocks within a particular organelle which exhibit high heterogeneity, e.g., within a cell nuclei, there are nucleolus, euchromatin and heterochromatin. On the other hand, a continuous R.I. fluctuations model better represents these highly heterogeneous and interconnected structures. It does not require structures to be discrete, allows R.I. variances at any length scales and more comprehensively quantifies the variation by its autocorrelation function. Although there is a fundamental difference between the discrete particle model and the continuous R.I. fluctuation model, we can still bridge their understanding by means of the R.I. correlation function. Under the Born approximation, the tissue R.I. correlation function can also be decomposed into an integration of the correlation function of spheres in infinite numbers of sizes, so that the two models should converge. As mentioned above, the correlation function from a sphere is a triangular-like function which contains entire length scales from 0 up to its diameter [38]. Alterations of a particular size of particles (e.g. increase volume fraction, increase of R.I. *etc.*) in the discrete model lead to a change of the R.I. correlation function from the sphere diameter and below. We should note that the particle size can be bigger than 450nm, but it is still the sub-micron length scales that prominently determine *D* measurement [Fig. 4(c) and 4(d)].

We demonstrated the three-dimensional capability of measuring *D* from biological tissue. The accuracy for *in vivo* applications relies on two central assumptions of a continuous R.I. fluctuation model: the Born approximation and a continuous medium. In ISOCT, the spectrum is extracted from a coherence volume (~10^{3} μm^{3}), within which the Born approximation is largely satisfied. However, in a highly scattering medium where multiple scattering dominates, the Born approximation no longer holds and measurements can be inaccurate. Previously, we have found that in the superficial layer, the backscattering coefficient has much higher tolerance to the stronger scattering medium (*μ _{b}* up to 10mm

^{−1}) while the scattering coefficient

*μ*started to break the Born approximation above ~12mm

_{s}^{−1}[29]. Thus, superficial compartments of tissue such as epithelial layers have a larger range of applicability for ISOCT.

The other key assumption is that tissue is modeled as a continuous medium. At microscopic level, this assumption is largely valid [26, 27]. Tissues also often possess stratified components such as epithelium, sub-mucosa *etc*. as in the buccal tissue presented in Fig. 8(a). There are distinct inter-layer boundaries where Fresnel reflections also contribute to the backscattered signal. We can treat each layer as a homogeneous medium within which the R.I. fluctuations are continuous and apply the algorithm. However, the deeper tissue is illuminated with wavelength-dependent attenuated incidence, which can be compensated for given the measurement of *µ _{s}(λ)*. Once the penetration depth is over several times the mean free path length (MFP)

*l*(

_{s}*l*= 1/ µ

_{s}_{s}), the multiple scattering has a significant effect [39], so that the compensation may not be sufficient and the first assumption may be defied. As a result, the

*D*measurement may be inaccurate.

To conclude, in this paper we demonstrate that OCT can be sensitive to nanoscale structural alterations by ISOCT analysis. Numerical studies showed that the length scales of sensitivity of our approach range from around 35nm to 0.45μm, given a 5% sensitivity threshold. Phantom studies and *ex vivo* human colon biopsies provided experimental verification of this nanoscale sensitivity.

## Appendix:

1. Relationship between Henyey-Greenstein (H-G) phase function and the spectral power density of W-M function when *D*=3.

Under the Born approximation, the spectral power density is the Fourier transform Φ of the W-M correlation function:

*k*is the momentum transfer whose value is equal to $\left|{k}_{s}\right|=2k\mathrm{sin}(\theta /2).$

_{s}*k*is the wave number and θ is the scattering azimuthal angle.

*D*,

*l*and

_{c}*N*follow the same definition in W-M correlation function. When

_{c}*D*= 3 and after plugging in the absolute value of

*k*, the above equation becomes

_{s}*p*(

*θ*) has the same shape as formulated above with the normalization so that the integration of the phase function is unity:

When$k{l}_{c}=\sqrt{g}/(1-g)$, the above equation can be further expressed as

which has the identical form as H-G phase function.2. Derivation of Eqs. (11) and (12)

For the lower length scale analysis, the original tissue R.I. fluctuation ${n}_{\Delta}$ (r) is convoluted with a blurring Gaussian window denoted here as *G*(*r*)

*l*denotes the lower length scale. According to Wiener–Khinchin theorem, the Fourier transform of the correlation function is equal to power density spectrum.

*FT*denotes Fourier transform and

*ACF*denotes autocorrelation function. Thus,

For retaining the small length scale and eliminating only the large length scale, the original fluctuation ${n}_{\Delta}$(r) is subtracted by the above filtered distribution:

where the subscript*h*stands for the upper length scale. A similar formulism in Eq. (18) is applied for

_{${n}_{\Delta}^{h}$}:

## Acknowledgments

The authors would like to acknowledge our financial support from NIH grant R01CA128641, R01EB003682 and NSF grant CBET-0937987.

## References and links

**1. **M. Cremer, K. Küpper, B. Wagler, L. Wizelman, J. von Hase, Y. Weiland, L. Kreja, J. Diebold, M. R. Speicher, and T. Cremer, “Inheritance of gene density-related higher order chromatin arrangements in normal and tumor cell nuclei,” J. Cell Biol. **162**(5), 809–820 (2003). [CrossRef] [PubMed]

**2. **C. Lanctôt, T. Cheutin, M. Cremer, G. Cavalli, and T. Cremer, “Dynamic genome architecture in the nuclear space: regulation of gene expression in three dimensions,” Nat. Rev. Genet. **8**(2), 104–115 (2007). [CrossRef] [PubMed]

**3. **O. Nadiarnykh, R. B. LaComb, M. A. Brewer, and P. J. Campagnola, “Alterations of the extracellular matrix in ovarian cancer studied by second harmonic generation imaging microscopy,” BMC Cancer **10**(1), 94 (2010). [CrossRef] [PubMed]

**4. **S. M. Pupa, S. Ménard, S. Forti, and E. Tagliabue, “New insights into the role of extracellular matrix during tumor onset and progression,” J. Cell. Physiol. **192**(3), 259–267 (2002). [CrossRef] [PubMed]

**5. **N. Théret, O. Musso, B. Turlin, D. Lotrian, P. Bioulac-Sage, J. P. Campion, K. Boudjéma, and B. Clément, “Increased extracellular matrix remodeling is associated with tumor progression in human hepatocellular carcinomas,” Hepatology **34**(1), 82–88 (2001). [CrossRef] [PubMed]

**6. **A. Lewis, M. Isaacson, A. Harootunian, and A. Muray, “Development of a 500 Å spatial resolution light microscope: I. light is efficiently transmitted through λ/16 diameter apertures,” Ultramicroscopy **13**(3), 227–231 (1984). [CrossRef]

**7. **E. Betzig, G. H. Patterson, R. Sougrat, O. W. Lindwasser, S. Olenych, J. S. Bonifacino, M. W. Davidson, J. Lippincott-Schwartz, and H. F. Hess, “Imaging intracellular fluorescent proteins at nanometer resolution,” Science **313**(5793), 1642–1645 (2006). [CrossRef] [PubMed]

**8. **S. W. Hell and J. Wichmann, “Breaking the diffraction resolution limit by stimulated emission: stimulated-emission-depletion fluorescence microscopy,” Opt. Lett. **19**(11), 780–782 (1994). [CrossRef] [PubMed]

**9. **M. J. Rust, M. Bates, and X. Zhuang, “Sub-diffraction-limit imaging by stochastic optical reconstruction microscopy (STORM),” Nat. Methods **3**(10), 793–796 (2006). [CrossRef] [PubMed]

**10. **M. G. L. Gustafsson, “Surpassing the lateral resolution limit by a factor of two using structured illumination microscopy,” J. Microsc. **198**(2), 82–87 (2000). [CrossRef] [PubMed]

**11. **I. Itzkan, L. Qiu, H. Fang, M. M. Zaman, E. Vitkin, I. C. Ghiran, S. Salahuddin, M. Modell, C. Andersson, L. M. Kimerer, P. B. Cipolloni, K. H. Lim, S. D. Freedman, I. Bigio, B. P. Sachs, E. B. Hanlon, and L. T. Perelman, “Confocal light absorption and scattering spectroscopic microscopy monitors organelles in live cells with no exogenous labels,” Proc. Natl. Acad. Sci. U.S.A. **104**(44), 17255–17260 (2007). [CrossRef] [PubMed]

**12. **Z. Wang, L. Millet, M. Mir, H. Ding, S. Unarunotai, J. Rogers, M. U. Gillette, and G. Popescu, “Spatial light interference microscopy (SLIM),” Opt. Express **19**(2), 1016–1026 (2011). [CrossRef] [PubMed]

**13. **C. Joo, T. Akkin, B. Cense, B. H. Park, and J. F. de Boer, “Spectral-domain optical coherence phase microscopy for quantitative phase-contrast imaging,” Opt. Lett. **30**(16), 2131–2133 (2005). [CrossRef] [PubMed]

**14. **A. J. Radosevich, J. Yi, J. D. Rogers, and V. Backman, “Structural length-scale sensitivities of diffuse reflectance in continuous random media under the Born approximation,” Opt. Lett. **37**, 5220–5222 (2012). [CrossRef] [PubMed]

**15. **A. J. Radosevich, N. N. Mutyal, V. Turzhitsky, J. D. Rogers, J. Yi, A. Taflove, and V. Backman, “Measurement of the spatial backscattering impulse-response at short length scales with polarized enhanced backscattering,” Opt. Lett. **36**(24), 4737–4739 (2011). [CrossRef] [PubMed]

**16. **A. J. Gomes, S. Ruderman, M. DelaCruz, R. K. Wali, H. K. Roy, and V. Backman, “In vivo measurement of the shape of the tissue-refractive-index correlation function and its application to detection of colorectal field carcinogenesis,” J. Biomed. Opt. **17**(4), 047005 (2012). [CrossRef] [PubMed]

**17. **D. Huang, E. A. Swanson, C. P. Lin, J. S. Schuman, W. G. Stinson, W. Chang, M. R. Hee, T. Flotte, K. Gregory, C. A. Puliafito, and et, “Optical coherence tomography,” Science **254**(5035), 1178–1181 (1991). [CrossRef] [PubMed]

**18. **A. F. Fercher, W. Drexler, C. K. Hitzenberger, and T. Lasser, “Optical coherence tomography-principles and applications,” Rep. Prog. Phys. **66**(2), 239–303 (2003). [CrossRef]

**19. **B. Povazay, K. Bizheva, A. Unterhuber, B. Hermann, H. Sattmann, A. F. Fercher, W. Drexler, A. Apolonski, W. J. Wadsworth, J. C. Knight, P. S. Russell, M. Vetterlein, and E. Scherzer, “Submicrometer axial resolution optical coherence tomography,” Opt. Lett. **27**(20), 1800–1802 (2002). [CrossRef] [PubMed]

**20. **L. Liu, J. A. Gardecki, S. K. Nadkarni, J. D. Toussaint, Y. Yagi, B. E. Bouma, and G. J. Tearney, “Imaging the subcellular structure of human coronary atherosclerosis using micro-optical coherence tomography,” Nat. Med. **17**(8), 1010–1014 (2011). [CrossRef] [PubMed]

**21. **J. Yi, Q. Wei, H. F. Zhang, and V. Backman, “Structured interference optical coherence tomography,” Opt. Lett. **37**(15), 3048–3050 (2012). [CrossRef] [PubMed]

**22. **A. Fercher, “Inverse Scattering, Dispersion, and Speckle in Optical Coherence Tomography,” Optical Coherence Tomography 119–146 (2008).

**23. **J. M. Schmitt, S. H. Xiang, and K. M. Yung, “Speckle in optical coherence tomography,” J. Biomed. Opt. **4**(1), 95–105 (1999). [CrossRef] [PubMed]

**24. **R. Barer and S. Tkaczyk, “Refractive index of concentrated protein solutions,” Nature **173**(4409), 821–822 (1954). [CrossRef] [PubMed]

**25. **J. M. Schmitt and G. Kumar, “Turbulent nature of refractive-index variations in biological tissue,” Opt. Lett. **21**(16), 1310–1312 (1996). [CrossRef] [PubMed]

**26. **M. Moscoso, J. B. Keller, and G. Papanicolaou, “Depolarization and blurring of optical images by biological tissue,” J. Opt. Soc. Am. A **18**(4), 948–960 (2001). [CrossRef] [PubMed]

**27. **C. J. R. Sheppard, “Fractal model of light scattering in biological tissue and cells,” Opt. Lett. **32**(2), 142–144 (2007). [CrossRef] [PubMed]

**28. **J. D. Rogers, İ. R. Capoğlu, and V. Backman, “Nonscalar elastic light scattering from continuous random media in the Born approximation,” Opt. Lett. **34**(12), 1891–1893 (2009). [CrossRef] [PubMed]

**29. **J. Yi and V. Backman, “Imaging a full set of optical scattering properties of biological tissue by inverse spectroscopic optical coherence tomography,” Opt. Lett. **37**(21), 4443–4445 (2012). [CrossRef] [PubMed]

**30. **İ. R. Capoğlu, J. D. Rogers, A. Taflove, and V. Backman, “Accuracy of the Born approximation in calculating the scattering coefficient of biological continuous random media,” Opt. Lett. **34**(17), 2679–2681 (2009). [CrossRef] [PubMed]

**31. **H. G. Davies, M. H. F. Wilkins, J. Chayen, and L. F. La Cour, “The use of the interference microscope to determine dry mass in living cells and as a quantitative cytochemical method,” Quarterly J.Microscopical Sciences (New York) **3**, 271–304 (1954).

**32. **L. G. Henyey and J. L. Greenstein, “Diffuse radiation in the galaxy,” Astrophys. J. **93**, 70–83 (1941). [CrossRef]

**33. **C. J. R. Sheppard, “Fractal model of light scattering in biological tissue and cells,” Opt. Lett. **32**(2), 142–144 (2007). [CrossRef] [PubMed]

**34. **C. Xu, J. M. Schmitt, S. G. Carlier, and R. Virmani, “Characterization of atherosclerosis plaques by measuring both backscattering and attenuation coefficients in optical coherence tomography,” J. Biomed. Opt. **13**(3), 034003 (2008). [CrossRef] [PubMed]

**35. **H. C. van de Hulst, *Light Scattering by Small Particles* (Dover publications, 1981).

**36. **J. Yi, J. Gong, and X. Li, “Analyzing absorption and scattering spectra of micro-scale structures with spectroscopic optical coherence tomography,” Opt. Express **17**(15), 13157–13167 (2009). [CrossRef] [PubMed]

**37. **F. E. Robles and A. Wax, “Measuring morphological features using light-scattering spectroscopy and Fourier-domain low-coherence interferometry,” Opt. Lett. **35**(3), 360–362 (2010). [CrossRef] [PubMed]

**38. **C. Matzler, “Autocorrelation functions of granular media with free arrangement of spheres, spherical shells or ellipsoids,” J. Appl. Phys. **81**(3), 1509–1517 (1997). [CrossRef]

**39. **R. K. Wang, “Signal degradation by multiple scattering in optical coherence tomography of dense tissue: a Monte Carlo study towards optical clearing of biotissues,” Phys. Med. Biol. **47**(13), 2281–2299 (2002). [CrossRef] [PubMed]