The biomechanical properties of tissues can provide crucial information for assessing tissue health and integrity. For example, tumors¹ and atherosclerotic plaques² exhibit distinct biomechanical properties from surrounding healthy tissues. Elastographic techniques such as ultrasound elastography (USE)³ and magnetic resonance elastography (MRE)⁴ are clinically valuable tools for detecting diseases, such as breast cancer⁵ and hepatic fibrosis.⁶ However, their spatial resolutions are limited, and these techniques require relatively large displacements to produce a detectable signal. Thus, USE and MRE are not generally used to characterize the biomechanical properties of small and thin samples.

Optical coherence elastography (OCE) is a rapidly emerging technique that can noninvasively measure localized mechanical contrast in tissue with spatial resolution at the micrometer scale.^7,8 Previous studies have demonstrated the use of OCE to characterize the biomechanical properties of various tissues, such as hyaline cartilage,⁹ cornea,¹⁰ crystalline lens,¹¹ breast malignancies,¹² prostate tissue suspect of cancer,¹³ and other soft tissue tumors.¹⁴

However, accurately quantifying biomechanical properties from OCE measurements requires the use of an appropriate mechanical model. We have previously developed a modified Rayleigh-Lamb frequency equation (mRLFE) to quantify the viscoelasticity of the cornea based on OCE measurements of an air-pulse induced elastic wave.¹⁵ This model assumes that the sample is a thin plate with a fluid–solid interface at the posterior surface. However, the geometry of many tissues such as the skin means that this model is not applicable.

In this study, we demonstrate that the feasibility of using the Rayleigh wave model (RWM) with OCE measurements to quantitatively assess the viscoelasticity of soft tissue. By utilizing phase velocities calculated from OCE-measured displacement profiles of an air-pulse induced elastic wave, the elastic wave dispersion curves were fitted to the RWM to quantify the Young’s modulus and shear viscosity.

The air-pulse OCE system was comprised of an air-pulse delivery device¹⁶ and a phase-stabilized swept source optical coherence tomography (PhS-SSOCT) system.¹⁷ The air-pulse system delivered a short-duration focused air-pulse ($\le 1\mathrm{ms}$) to the excitation point, which was at the central region of the sample upper surface, to induce a small-amplitude deformation (order of micrometer) that propagated transversely as an elastic wave. The elastic wave was imaged by the OCT system by synchronizing multiple air-pulse excitations with the OCT system frame trigger.¹⁰ However, the penetration depth of the OCT probe beam was only $\sim 1\mathrm{mm}$ in the tissue, with a total sample thickness of over 10 mm. The phase data were corrected for the surface motion and refractive index mismatch.¹⁸ The phase velocities, which are intrinsic mechanical characteristics of the sample, were calculated by fitting the phase shifts $\mathrm{\Delta }\theta$ obtained by fast Fourier transform (FFT) of the axial temporal displacement profiles to the elastic wave propagation distances ($\sim 1$ to $\sim 2\mathrm{mm}$), $\mathrm{\Delta }r$, by $c_p=\omega \mathrm{\Delta }r/\mathrm{\Delta }\theta$, where $\omega$ was the angular frequency.¹⁹ The phase velocity calculation was then repeated for each in-depth layer of the sample, and then averaged depth-wise (over $\sim 0.2$ to $\sim 0.5\mathrm{mm}$ depending on sample optical properties) to obtain a dispersion curve for a given sample.

The viscoelastic properties were obtained by fitting the analytical solution of the air-pulse induced elastic wave to the OCE-measured phase velocities. Given the limited depth penetration of the elastic wave, the wave can be modeled as a Rayleigh wave²⁰

This model requires that the sample upper surface is free (no stress), and that the sample is of infinite thickness. In contrast, our previously developed mRLFE model for the cornea required finite thickness, a free upper surface, and a fluid–solid interface at the posterior surface.¹⁵ In Eq. (1), $k=\omega /c_p$ is the wave number at angular frequency $\omega =2\pi f$, with $f$ being the frequency, and

Here, $\rho$ is the material density, and $\lambda$ and $\mu$ are the Lamé constants. For an elastic material, $\lambda =E\nu /[(1+\nu )(1-2\nu )]$ and $\mu =E/[2(1+\nu )]$, where $E$ is the Young’s modulus and $\nu$ is the Poisson’s ratio. Assuming the sample is a Kelvin–Voigt viscoelastic material, the dynamic shear modulus is complex: ${\mu }_D=\mu +\mathrm{i}\eta \omega$, where $\eta$ is the shear viscosity, and $i$ is the imaginary unit. It should be noted that the widely used surface wave equation, $c_R=[(0.87+1.12v)/(1+v)]{(\mu /\rho )}^{0.5}$, is the approximate solution for Eq. (1) when viscosity is neglected.²⁰

We first performed preliminary OCE measurements on 10% gelatin phantoms (Type A gelatin, 250 Bloom/8 Mesh, PB Gelatins/PB Leiner, Iowa) with the addition of 0%, 10%, 20%, and 30% castor oil (Walgreens Castor Oil, Walgreens, Illinois) ($\mathrm{w}/\mathrm{w}$, $n=3$ for each type) to artificially control viscosity.²¹ Milk was added to increase scattering. Each phantom sample had identical dimensions: diameter $D=50\mathrm{mm}$ and height $H=11\mathrm{mm}$. This height was greater than the longest wavelength of recorded elastic waves ($\sim 9\mathrm{mm}$), relatively satisfying the infinite depth assumption for a Rayleigh wave.^18,22

Measurements were taken at the central region of the upper surface. OCT images of a gelatin phantom without oil and with 20% oil are presented in Fig. 1. The addition of oil increased back-scattering due to the presence of more refractive index mismatches. The viscoelasticity of phantoms was obtained by an iterative error minimization procedure, where the difference between the OCE measurements and analytical solution of the RWM was minimized by altering the Young’s modulus and shear viscosity. To assess the accuracy of the elasticity assessment with OCE and the RWM, uniaxial mechanical compression testing (Model 5943, Instron Corp., Norwood, Massachusetts) was performed on the phantoms ($n>9$ for each concentration) after the OCE measurements, where the effect of sample geometry was considered.²³

Fig. 1

OCT images of (a) a gelatin phantom and (b) a 20% oil-in-gelatin phantom sample. The elastic wave analysis region is highlighted in yellow, and the air-pulse excitation point is marked by the red arrow.

Figure 2 shows the gelatin phantom viscoelasticity assessments by OCE and the RWM. Poisson’s ratio was 0.4999 to account for the near incompressibility of the phantoms, and the density was $1000\mathrm{kg}/{\mathrm{m}}^3$. Figure 2(a) plots the phase velocities as measured by OCE and fitted analytical solution of the RWM for one typical sample with the addition of 0%, 10%, 20%, and 30% ($\mathrm{w}/\mathrm{w}$) oil. The Young’s moduli of the phantoms with different oil concentrations were estimated as 10.32, 11.25, 12.15, and 12.70 kPa, and the shear viscosities were 2.27, 2.41, 3.04, and 4.03 Pa s for 0%, 10%, 20%, and 30% oil, respectively. The elastic wave wavelength was estimated as $\sim 3$ to $\sim 9\mathrm{mm}$. The approximate lateral resolution in these phantoms was $\sim 4.5\mathrm{mm}$ (a half-wavelength of the longest wavelength).²⁴ Other factors such as displacement signal-to-noise ratio also influence the spatial resolution.²⁵

Fig. 2

(a) Typical examples of the air-pulse induced elastic wave dispersion calculated from OCE measurements for 10% gelatin phantoms with additional 0%, 10%, 20%, and 30% oil and the fitted analytical solution of the RWM. (b) Young’s modulus of the gelatin phantoms with oil concentrations of 0%, 10%, 20%, and 30% as assessed by OCE and the RWM ($n=3$ for each concentration) as compared to uniaxial mechanical compression testing (MT) ($n>9$ for each concentration). (c) Shear viscosity of the phantoms as assessed by OCE and the RWM ($n=3$ for each concentration).

Figure 2(b) compares the elasticity values of the gelatin phantoms as estimated by OCE and the RWM and as measured by uniaxial mechanical compression testing (MT). The results demonstrate that the Young’s modulus as estimated by OCE and the RWM generally agreed with the stiffness as measured by MT. Figure 2(c) shows that the shear viscosity of the 10% gelatin phantoms with 0%, 10%, 20%, and 30% oil was $2.18\pm 0.12$, $2.72\pm 0.28$, $2.99\pm 0.05$, and $3.96\pm 0.22\mathrm{Pa}\mathrm{s}$ as assessed by OCE and the RWM, respectively. The viscosity of the phantoms increased with oil concentration, which corroborates with the literature.²¹ Therefore, the combination of OCE and the RWM can quantify the viscoelasticity of thick samples.

Next, OCE experiments were performed on in vitro chicken liver, where the dimensions of the smallest sample were $D\times H=30\mathrm{mm}\times 10\mathrm{mm}$. For each of the three samples, five OCE measurements were made across the central region of the sample, and all three samples were optically similar. During the calculation, the material density of the liver was assumed as $\rho =1100\mathrm{kg}/{\mathrm{m}}^3$ (close to water density) while the Poisson’s ratio was $\nu =0.4999$ to account for the near incompressibility of biological tissue.²⁶

Figure 3 plots the viscoelasticity assessment by the RWM of three chicken liver samples. Figure 3(a) shows the OCE-measured phase velocities of two selected samples and viscoelasticity quantification by the RWM. The Young’s modulus of one measurement point in the sample 1 was $E=3.10\mathrm{kPa}$ and shear viscosity $\eta =1.40\mathrm{Pa}\mathrm{s}$ (with $R^2=0.9956$). In sample 3, the RWM estimated $E=1.52\mathrm{kPa}$ and $\eta =1.33\mathrm{Pa}\mathrm{s}$ ($R^2=0.9905$). Similar to the phantoms, the lateral resolution in these samples was approximated from half of the longest wavelength as $\sim 3\mathrm{mm}$. Furthermore, the infinite thickness assumption was relatively satisfied as the thinnest sample was 10 mm, greater than the longest wavelength of $\sim 6\mathrm{mm}$.

Fig. 3

(a) The estimated ($E$) Young’s modulus and ($\eta$) shear viscosity of two chicken liver samples as assessed by the OCE and the RWM. (b) Comparisons of the Young’s modulus and shear viscosity of each of the three liver samples ($n=5$ OCE measurements for each sample).

Figure 3(b) shows the average viscoelasticity for each of the three samples as assessed by the OCE and the RWM. The Young’s moduli of the three samples were $3.19\pm 0.39$, $1.38\pm 0.17$, and $1.56\pm 0.08\mathrm{kPa}$. The large degree of intersample variability was most likely because the samples were from different animals.²⁷ The shear viscosity [striped bars in Fig. 3(b)], was estimated as $1.27\pm 0.02$, $1.09\pm 0.01$, and $1.23\pm 0.01\mathrm{Pa}\mathrm{s}$. Correlation between the analytical solution to the RWM and OCE-measured phase velocities for the three samples was $R^2=0.94\pm 0.06$, $0.96\pm 0.03$, and $0.98\pm 0.01$. The averaged Young’s modulus was estimated by OCE, and the RWM for all three samples was $E=2.04\pm 0.88\mathrm{kPa}$ and shear viscosity was $\eta =1.20\pm 0.13\mathrm{Pa}\mathrm{s}$ with an average $R^2=0.96\pm 0.04$. The error is intrasample standard deviation. The Young’s modulus of human liver as measured by transient elastography is between 2 and 10 kPa (Ref. 27) and between 2 and 5 kPa as measured by MRE,²⁸ which is in the range of our measurements.

In summary, we have utilized the RWM in conjunction with air-pulse OCE measurements to quantify the viscoelasticity of chicken liver tissue by utilizing OCE-measured phase velocities and the RWM. The accuracy and feasibility of this model was verified with gelatin phantoms without and with oil. The RWM was able to accurately quantify the elasticity as compared to uniaxial mechanical compression testing, which corroborates with our previous work.²⁶ Furthermore, the RWM was able to assess the increased viscosity due to the addition of oil in the gelatin phantom. To test the feasibility of the presented method on biological tissue, the viscoelasticity of chicken liver was quantified. The combination of the RWM and OCE may be a potentially useful method for quantifying the viscoelasticity of tissues for early disease detection.