The problem of elasticity imaging (elastography) in optical coherence tomography (OCT) has been attracting a great deal of attention since the end of 1990s (e.g., Refs. 12.3.–4) following a similar trend in ultrasonic medical imaging.5,6 Although in medical ultrasonics the elastography regime has already been implemented in several scanner platforms,6 in OCT still there are no similar commercially available scanners.
For most of the biological tissues (excluding such tissue types as dental enamel, bones, and cartilages) the Poisson’s ratio is very close to the value typical of liquids. In such a case, the material stiffness can be conveniently characterized by the shear modulus determining the velocity of shear-wave propagation or the Young modulus , for which in the limit , and the simple relation is known.7 In medical ultrasonics, one of the approaches to evaluating the tissue stiffness is based on excitation of a shear wave, the velocity of which is measured using ultrasonic methods. In the problems typical of OCT, usually excitation of the bulk shear waves is very problematic so that the surface Rayleigh waves can be used instead. For a Rayleigh wave, its velocity is close to that of bulk shear waves and weakly depends on the Poisson’s ratio of the material. It can be shown,8 in particular, that the Rayleigh-wave velocity . Thus, the Rayleigh-wave velocity can be used for estimation of the tissue stiffness, see e.g., Ref. 9 and references therein.
Other conventionally considered in OCT approaches to elastographic mapping are mainly based on one or another way of the initial reconstruction of particle-displacement field in a strained tissue using special processing of conventional OCT images.1 In particular, correlation processing of OCT images is often discussed for reconstructing displacements produced by compression and/or shearing of the studied tissue. Then, in such a displacement-based (DB) approach, the displacement field can be differentiated to determine the local elastic strains and, correspondingly, to estimate the shear-modulus distribution.1–3 This approach is very similar to the widely known digital image correlation technique10 used in engineering applications. In the DB approach to elastography, the reduction of cross-correlation between fragments of consequently obtained images (“decorrelation noise”) caused by distortions of scatterer patterns in deformed tissues is the main negative factor that reduces the accuracy of the displacement-field reconstruction11 and strongly complicates its practical implementation. Besides cross-correlation methods (usually implying one or another variant of super-resolution in order to detect subpixel displacements), other methods, e.g., ones based on measurement of optical-wave phase, can be used for estimating the displacement and strain fields in the tissue. A common feature of all such DB approaches is that the significant distortion of the scatterer pattern in the tissue leads to a strong increase in the measurement errors and is considered as a negative factor for solving the problem of elasticity mapping.
However, by analogy, with the use of nonlinear distortions in nonlinear-acoustic approach12 to nondestructive testing, just the image distortions causing reduction in cross-correlation between the deformed tissue images can be used for evaluation of the tissue stiffness. Indeed, for regions with higher stiffness that better conserve their shape, the distortions of the scatterer pattern are smaller and, correspondingly, the reduction in the cross-correlation between the consequent OCT images is also smaller. Thus, the comparison of the correlation stability of different parts of OCT images can be directly used for evaluation of spatial distribution of tissue stiffness avoiding the intermediate stages of the displacement-field reconstruction. Besides, the procedure of numerical differentiation of the experimentally determined displacement field is very sensitive to measurement inaccuracy and may introduce strong errors in the reconstructed strain field. According to the basic idea of evaluating the relative stiffness of different tissue regions via the correlation stability of compared OCT images, in what follows, we call such a method of characterizing the relative tissue elasticity the “correlation stability approach” (CS approach).
In our previous works,13,14 we discussed the implementation of elastographic CS mapping based on cross-correlating the compared OCT images at the level of morphological features of the deformed tissue. For implementation of this variant of the CS approach, the speckle-level structure of OCT images should be essentially eliminated. To this end, procedures of spatial filtering and temporal averaging of consequently obtained OCT B-scans can be used as demonstrated in Refs. 13 and 14. However, speckle structure of OCT images is a very important intrinsic feature of OCT images. The speckle-level decorrelation usually occurs for significantly smaller strains than the decorrelation due to purely geometrical distortion of the patterns of morphological features of the tissue. Although the distortion of OCT images at the speckle level is essentially different from purely geometrical distortions, the basic idea of the CS approach that the stiffer regions exhibit better cross-correlation remains valid. Thus, the speckle-level decorrelation can efficiently be used for elasticity mapping even for tissues without pronounced morphological features. In the following sections we discuss both the general features of the CS approach to elasticity mapping and the differences between realization of CS mapping either at the morphological-structure level or at the speckle level. Results of numerical simulations and experimental results for both variants of the CS approach are also presented using both phantoms and in vivo obtained images.
Basic Correlation Procedures to be Used in CS-Based Elasticity Mapping
In the proposed CS approach, the basic procedure of cross-correlation between the reference and deformed images is similar to the conventional calculation of the two-dimensional correlation coefficient using a sliding window.1–3 Usually such a window corresponds to rectangular subsets of elements taken form the compared pixelized images, so that the cross-correlation coefficient in plane is given by1) with similar in size areas from the deformed image; the quantities and are the respective mean values over the correlated areas in the reference and deformed images. If the compared images are identical, coefficient (1) reaches its maximum unity value for . If the compared image is distorted, the coordinates and are found, for which coefficient (1) reaches maximum within a certain search area around the initial position of the window in the reference image. The coordinates of the area corresponding to the maximum cross-correlation correspond to the new position of the same group scatterers from the subset in the reference image. Figure 1(a) (the solid-line rectangle) shows the initial subset taken from the reference image, and the subset from the deformed image [the dashed-line rectangle in Fig. 1(b)] corresponds to the new position of the same group of scatterers. Usually the scatterers experience not only their translational displacement but also their mutual positions are somewhat distorted due to the tissue deformation [compare the solid-line rectangle in panel (a) and the dashed-line rectangle in panel (b) of Fig. 1]. Stiffer regions experience smaller distortions than the surrounding softer regions and therefore demonstrate higher cross-correlation maxima. The maximum values of correlation function (1) characterize the degree of distortion of the scatterer patterns from different subsets and, therefore, the map of the so-found maxima (CS map) characterizes the relative stiffness of differently deformed regions of the tissue.
We emphasize that when performing the CS mapping, we do not try to minimize the deformation-induced worsening of the cross-correlation, but intentionally use the difference in the reduction of cross-correlation between images of different tissue regions to characterize their relative stiffness. Thus, the error-sensitive procedures of the displacement-field reconstruction and its numerical differentiation are avoided. Furthermore, for clear distinction between the areas with different stiffness, the scatterer pattern in the softer areas should already be sufficiently strongly distorted. Correspondingly, the correlation coefficient (1) in those areas should be significantly reduced. This requirement is in contrast with the conventional DB approach, for which the tissue deformation is supposed to be fairly small to ensure that the cross-correlation does not decrease significantly in order to achieve sufficiently high accuracy of the reconstruction of the displacement field. These qualitative arguments indicate that the CS approach intrinsically requires significantly larger strains than the conventional DB approach. The larger operability range in terms of strain should be favorable for free-hand implementation of the CS approach, which is the most interesting for practical applications.
It should also be emphasized that the CS mapping differs from the so-called correlation mapping OCT (cmOCT).15 In cmOCT, correlation is evaluated without moving the subset window over the search area, so that eventual translational motion of the scatterers also contributes to the reduction in the cross-correlation and masks the effect of mutual motion of the scatterers within each subset. Thus cmOCT does not characterize relative rigidity of different regions, but can efficiently detect regions with moving scatterers (i.e., it is intended for mapping microcirculation and flows10).
It is clear from the above-formulated basic idea of the CS approach to mapping the relative elasticity that the key point is the sufficiently strong reduction in the correlation coefficient (1). Understanding of the reasons for such reduction makes it possible to evaluate the operability region of the CS approach in terms of strains. In this context, as emphasized in Sec. 1, two main reasons of the correlation-coefficient reduction should be clearly distinguished.
The first one is related to the purely geometrical distortion of the pattern of structural elements in OCT images. Under the “structural elements” here, we understood the visualized morphological tissue features with the characteristic sizes exceeding the OCT-system resolution (typically 5 to 20 μm), for instance, such structures as hair roots of similar features of skin images. Those morphological structural elements in biological tissues may experience deformations thus contributing to a decrease in the cross-correlation of the compared images.
The other reason of the cross-correlation reduction is the peculiar phenomenon of speckle “blinking” or “boiling,” which is typical of coherence imaging techniques like OCT. The speckle formation is related to the interference of scattered fields from subresolution scatterers. Since deformation of the tissue affects the mutual positions of not only macroscopic (resolved) image features but also subresolution scatterers, the interference of the respective scattered fields may change from constructive to destructive (or vice versa). This results in the peculiar speckle blinking and boiling. However, it is clear that in stiffer regions such distortions are weaker, so that the decorrelation due to speckle blinking/boiling should also be stronger in softer regions.
It is essential that quite different strains may be required for significant decorrelation due to purely geometrical distortions of resolved morphological features and due to speckle boiling related to mutual displacements of subresolution scatterers. However, both effects in principle may be used for implementation of the elastographic CS mapping. It the following sections we consider both variants.
CS Approach Based on Purely Geometrical Distortions of Morphological Features in OCT Images
Qualitative Arguments and Numerical Simulations
First, we discuss decorrelation of the compared OCT images due to purely geometrical distortions of macroscopic (i.e., resolved) morphological features as if they were visualized by incoherent light. The most straightforward way to obtain such speckle-free images is the application of spatial filtering (that averages neighboring speckles) combined with temporal averaging (filtering) of consequently obtained OCT images of gradually deformed tissue by a sufficiently fast OCT scanner. The latter procedure averages speckle boiling in time. Certainly, such processing somewhat worsens the resolution but makes it possible to obtain much stabler images, for which the cross-correlation is determined by geometrical distortions of morphological features rather than by blinking and boiling of individual speckles.
We assume that for such an image, the correlation window in size contains independent structural features. Then for nonoverlapping windows containing completely independent sets of such structural elements, the average absolute value of cross-correlation is (see e.g., Ref. 16). For example, for , the background correlation is , which may correspond to the correlation window size to 40 and the characteristic scale of independent elements of 5 to 7 pixels.
The condition of significant reduction of the cross-correlation due to purely geometrical distortions is intuitively clear. Indeed, the dissimilarity between the cross-correlated images becomes noticeable if the image fragment corresponding to the correlation window is deformed by a value of the order of the size of one element. Assuming for simplicity that the window is rectangular, , this condition is reached for a characteristic strain , where the proportionality factor is on the order of unity. Its exact value depends on the type of straining as well as on the assumed level of the cross-correlation reduction and can be more accurately evaluated by a fairly simple numerical modeling. For example, Fig. 2 shows such results of numerical simulations of the cross-correlation reduction due to geometrical distortions of the scatter pattern. In Fig. 2(a), the horizontal axis corresponds to the number of independent scatterers that can be placed along one side of a rectangular correlation window (such that their total number in the window is ) and the vertical axis shows the characteristic strain . The plot is obtained for the vertical compression of the tissue with free boundaries in the lateral directions, which is close to the case of compressing the tissue by a rigid cylinder piston (like OCT probe) under the assumption that the characteristic depth is not much smaller that the piston radius. The slope of the log–log plot in Fig. 2(a) is close to minus unity and agreescwell with the expected law with the proportionality factor about 1.4 for the chosen reduction in the cross-correlation by a factor of .
Next, Fig. 2(b), for the shear-type staining of the tissue, shows the simulated dependence of correlation coefficient given by Eq. (1) on the created strain using a rectangular window with the size corresponding to . This plot also shows that for the discussed decorrelation of purely geometrical origin, significant reduction in the correlation coefficient requires strains on the order of several tens of percents. Evidently, the stronger the contrast in the stiffness between different parts of the imaged region, the wider is the range of strains, for which the softer parts are already significantly decorrelated, whereas the stiffer ones still remain fairly weakly distorted and keep high correlation. Therefore, for quite reasonable values , the discussed characteristic strains can be on the order of (20% to 60%) and even more, which is favorable for free-hand implementation of the CS approach.
More specifically, the procedure of utilization of geometrical distortions for elastographic CS mapping of a layer with contrasting stiffness can be demonstrated using simulated images with a typical for OCT size of . The initial pattern of “scatterers” in the image is formed by putting a random value in each pixel. Then, Fourier filtering is applied to smooth the image and obtain intensity inhomogeneities with correlation properties similar to those typical of real OCT images subjected to the above described procedures of speckle averaging and filtering. The main requirement is to ensure the same average absolute value of background correlation , as for real OCT images. For the latter, quite a typical value of is about 0.3 to 0.25 for correlation windows 20 to 40 pixels in size (similar to those used in Refs. 2 and 3). In agreement with the above-presented arguments, this value corresponds to characteristic scales of independent structural features about 5 to 7 pixels.
The so-prepared simulated images then have to be properly distorted to simulate tissue deformation. In real experiments, the OCT probe acts as a rigid piston producing in the tissue rather complex strain fields combining axial and lateral distortions. In Ref. 14, the idea of CS mapping was illustrated by simulated images subjected to shear deformation. Here, we present another example that is closer to deformation produced by a rigid piston pressed onto the tissue in the axial direction. For simplicity, we consider a structure containing three layers which are compressed by uniaxial stress in the normal direction. We also assume that the layers can freely slide along their interfaces. This allows one to use a simple analytical solution even for inhomogeneous tissue containing a layer with a different value of the shear modulus , where is the shear modulus of the surrounding tissue and factor is the contrast. If the total sample thickness is , the stiffer layer with the thickness is located between layers with thicknesses and (see Fig. 2) and the applied vertical stress is , then the vertical displacements and lateral displacements inside the sample are given by2), we also took into account that in biological tissues, for which the Poisson’s ratio is close to , the Young modulus . The patterns of the scatterers in the deformed sample can be found using standard procedures of image resampling with the proportions of the compression/stretching defined by solution (2).
Figure 3 shows the simulated initial and distorted according to Eq. (2) images, in which the stiffer layer is not yet visible in the pattern of the scatterers [although in the simulated image Fig. 3(b), it is noticeable due to smaller lateral expansion]. An example of visualization of this stiffer layer using the CS mapping is shown in Fig. 4 for different contrast coefficients . The CS maps obtained by cross-correlating the initial and distorted images show that the stiffer layer is clearly seen due to increased correlation coefficient for the stiffer layer as compared with the surrounding softer regions. The softer regions experience stronger distortions, so that the cross-correlation coefficient for them can be reduced down to the background value to 0.2. The latter condition is reached for strains 30% to 50% and even higher, although within the stiffer layer the cross-correlation still remains fairly high (0.6 to 0.8). As discussed in more detail in Refs. 13 and 14, for such strains, the DB approach, implying cross-correlation procedures for initial reconstruction of the displacement field, does not work because the deformation-produced “decorrelation noise” introduces a huge amount of errors in reconstruction of the displacement field.
In Vivo Demonstrations of Elastographic CS Mapping Based on Geometrical Mechanism of Decorrelation at the Level of Morphological Structural Features
The experimental demonstrations of the above-described approach were made using a recently developed spectral-domain OCT scanner based on principles described in Refs. 17 and 18. It acquires images with a rate of 21 fps, axial resolution 10 μm, and lateral one 20 μm. The fiber-optic OCT probe is 2 mm in diameter and acts as a piston producing deformation of the tissue during freehand operation. The compression of the tissue by the probe is schematically shown in Fig. 5. The probe has been designed first of all for gastro-enterological applications, so that the fiber-optic probe body is flexible, whereas the height of the solid forward-looking head is short enough (15 mm) to ensure compatibility with standard gastrofibroscopes. In view of the strict limitations on the probe-head size, a compact optical scheme has been chosen in which the lateral scanning is ensured via tilting of the OCT beam (see Fig. 5). Consequently, the resultant B-scans have a sector shape similar to that typical of ultrasound scanners. The so-obtained OCT images were subjected to filtering, normalization of the amplitudes and averaging of a few frames for reduction of the speckles and retaining only larger-scale morphological features. Then, the compressed-tissue images are cross-correlated with a reference one.
For experimental demonstrations of the CS approach based on processing of morphological features in B-scans, we decided to use in vivo images in which natural morphological features with contrasting stiffness were present for sure. In particular, for such examples, we chose images of human cheek skin, where hair roots served as stiffer inclusions. Another example is the mucous tissue of a human lip in which a stiffer inclusion has been also revealed (see Fig. 6). In both cases, attentive examination of the series of consecutive B-scans made it possible to discern that the image regions, for which the correlation processing revealed stabler cross-correlation, indeed better maintained their shapes. However, without the preliminary CS mapping, such purely visual comparison of the shape stability for different image regions was practically impossible. The resulting elastographic CS images in Fig. 6 are shown in a color palette overlapped with the intensity OCT images (by analogy with similar overlapping of elastographic and conventional images in ultrasonic scanners). The stiffer features in Fig. 6 are clearly seen as regions of increased cross-correlation that are shown by the red color by analogy with the simulated CS maps in Fig. 4. The average strain of the tissue is about 20% to 50% for the images in Fig. 6. The tissue straining is produced by the OCT probe in the free-hand mode and the characteristic numbers of processed B-scans is several tens.
Correlation Stability Approach Based on the Speckle-Level Decorrelation of OCT Images of Deformed Tissues
Main Features of Speckle-Level Decorrelation Essential for CS Mapping and Simulation Examples
The peculiar speckle structure intrinsic to OCT images is notoriously known19 and is usually considered as an interference factor that should be suppressed for improving the OCT image quality, although in some works just the speckle structure is considered as an informative factor. For example, starting from the earliest publications on elastography in OCT (e.g., Ref. 1), the use of speckle tracking is discussed for estimations of local stains. Measurements of speckle variance can be used as a basis for mapping microcirculations in OCT.20,21 In this section, we consider feasibility of the CS approach to elasticity mapping using the speckle-level cross-correlation stability of compared OCT images. This approach can be especially attractive for examination of tissue regions in which pronounced morphological features are absent, so that after speckle averaging and filtering there would be no sufficiently clear image features that could be cross-correlated. Nevertheless, even if the subresolution scatterers are distributed in the studied region fairly uniformly, the interference of optical fields scattered by them should produce the peculiar speckle structure of the OCT image. In such a context, the speckle-level heterogeneity of the image plays a useful role and can also be used for implementation of the elastographic CS mapping.
First, we recall that the origin of speckles is related to the interference of optical fields scattered by subresolution scatterers located within one resolution volume as is schematically shown in Fig. 7. The difference in the optical path lengths for the light scattered by the scatterers located within the same resolution volume is evidently determined, first of all, by the difference in the axial coordinates of the scatterers. If the tissue is strained, so that the local strain is , the variation in the character of interference of such scattered fields is determined by the variation in the optical-path lengths . Here, the factor 2 accounts for the propagation of the optical wave first towards the scatterers and then back to the photodetector. Evidently, the axial separation of such scatterers in one resolution volume is limited by the axial resolution determined by the coherence length of the light source, i.e., . Significant change in the character of the interference (i.e., the change from constructive interference to destructive one or vice versa) occurs if the strain-induced variations in the optical paths lengths for the subresolution scatterers become close to , where is the optical-wave length. Consequently, the characteristic strain corresponding to the onset of speckle blinking is determined by the condition . For typical OCT scanners with values to 20, the corresponding characteristic strain is to 2.5%. This value is at least an order of magnitude smaller than that discussed in the previous section values of required for significant decorrelation due to purely geometrical distortion of the image. Therefore, if the speckle structure of OCT images is not smoothed, the role of purely geometrical distortions in the decorrelation is usually negligible in comparison with the speckle blinking/boiling. However, just like for decorrelation of purely geometrical origin, the speckle-level decorrelation is also determined by the local strains and thus the stiffer regions should demonstrate stabler correlation, albeit for significantly smaller strains. Consequently, the speckle-level correlation stability can be used to distinguish the difference in the local stiffness even for samples without pronounced morphological structure.
The above-presented arguments can be supported by both numerical simulations and experimental demonstrations. The numerical simulations are interesting by themselves because they allow one to show the effect of mutual displacement of the subresolution scatters, either accompanied by geometrical distortions or demonstrated by the most pure (although artificial) form of the effect. Namely, this assumes retaining only mutual phase variation for the signals from the subresolution scatterers and completely excluding the macroscopic geometrical distortions of intensity spots at the image. The latter effect is expected to influence the decorrelation only for much greater deformations. For the simplest form of simulating the speckle blinking, we assumed the presence of two subresolution scatterers with identical scattering properties in one resolution volume. Correspondingly we initially prepared a pair of identical random matrices such as described in the previous section. Then, the elements of one of the matrices were multiplied by a random factor imitating the phase difference between the contributions of the spatially separated scatterers, where the uniformly distributed random value corresponded to the random initial difference in the axial coordinates of the subresolution scatterers and is the local strain. Besides, the spatial spectrum of so-obtained matrices was subjected to filtering, although in the discussed case, only the minimal filtering was made to ensure the Nyquist–Kotelnikov sampling criterion.
The superposition of the two prepared matrices for imitates interference of the fields scattered by each pair of unperturbed subresolution scatterers and produces the initial image with characteristic speckle-like inhomogeneities very similar to the ones observed in real OCT images. Then increase in the strain results in gradual transformation of the image with peculiar speckle blinking and boiling. This blinking is random in appearance although its reason is the strain-induced regular variation in the separation of the subresolution scatterers. The difference in the local stiffness of different regions of the imaged tissue can also readily be introduced in the simulation by simple substitution of the strain with for the stiffer region, where the factor has the same meaning of the stiffness contrast as in the previous section. Then for the so-obtained images corresponding to gradually increasing strain, cross-correlation coefficient (1) can be calculated.
Figure 8 shows results elucidating some important features of the speckle-level decorrelation. Figure 8(a) demonstrates the difference in the decorrelation that is calculated first for pure geometrical distortion of the speckle spots and then accounting for the joint influence of the speckle blinking and geometrical distortions. The correlation-window size is in both cases. The plot clearly demonstrates that for the same characteristic size of intensity spots, the effect of speckle blinking strongly dominates in the resulting decorrelation. For comparison, in Fig. 8(b) the decorrelation curve is calculated retaining only speckle blinking (i.e., only the effect of mutual dephasing of the contributions of the subresolution scatterers), whereas the influence of geometrical distortions of the speckle spots is excluded in the simulation. The curve in Fig. 8(b) for pure speckle blinking is very close to the curve in Fig. 8(a) for the joint effect of the speckle blinking and geometrical distortions, which confirms that usually . Figure 8 also demonstrates that even for randomly distributed scatterers, the above-made simple estimate of gives quite a correct value of the characteristic strain , for which the speckle-level cross-correlation becomes significantly reduced. It can be also noted that with increasing strain, peculiar partial recurrence of the cross-correlation is observed. This is explained by the fact that the character of interference (constructive or destructive) is repeated if the axial distance between the scatterers located at the opposite sides of the resolution volume is changed by , where is an integer number. However, for increasing strain, in average the correlation is gradually reduced, because the closer located subresolution scatterers give increasing contribution to the decorrelation for monotonically increasing strain. Depending on the density and strength of the subresolution scatterers this partial recurrence can be more or less pronounced in different parts of the simulated images. It is interesting to mention that in some real physical experiments with phantoms, such a partial recurrence of the cross-correlation was sometimes noticeable, which confirms that such a simple model indeed fairly well captures the origin of the speckle blinking/boiling.
Figure 9 demonstrates a simulated example of the speckle-level elastographic CS mapping for a two-layer sample with a twice stiffer lower layer. For this example, in the simulation the local strain in the lower half of the image is made smaller than in the upper part by the stiffness-contrast factor . For a single frame with simulated speckle pattern shown in Fig. 9(a), the difference in the speckle structure for the two layers is not visible at all, whereas in the color CS map shown in Fig. 9(b) and corresponding to the average strain , the lower stiffer layer is perfectly visible due to better correlation stability of the speckle pattern.
Finally, Fig. 10 shows the examples of the reconstructed displacements of the tissue particles found by searching the maximum of the cross-correlation coefficient (1) using a moving correlation window. For detecting subpixel displacements, a super-resolution technique (based on parabolic approximation of the correlation-function maximum) was applied. The curves in Fig. 10(a) (for vertical and horizontal displacement components) correspond to the processing of the simulated images, for which the speckle blinking is excluded and only purely geometrical distortions of the intensity spots are retained. Although deviations of the reconstructed curves from the expected perfect linear dependences are noticeable, these curves fairly well agree with the theoretically expected displacement for the strains used in the simulation. In Fig. 10(b), the displacements for the same point of the sample are found using the same correlation procedures but for the images distorted geometrically with allowance for speckle blinking. It is clear that in Fig. 10(b), the actual displacements are practically completely masked by apparent displacements related to the speckle blinking/boiling.
We emphasize that the reconstructed apparent displacements shown in Fig. 10(b) are fake ones because the actual strain-produced displacement of the speckle spots is almost completely masked by the speckle blinking/boiling induced by the same strain. Furthermore, this statement relates not only to relatively large strains but also to fairly small strains below the characteristic strain of speckle blinking. This conclusion is important in view of the widely accepted opinion that, for deformations within the first 1% to 2%, the tissue displacements should be measurable via speckle tracking using correlation methods. Such speckle tracking should be feasible if translation-like collective motion of scatterers dominates the influence of their relative displacements, e.g., in the case of a fluid flow that moves a group of scatterers without significant change in their relative positions. The extreme situation of completely masked actual displacements shown in Fig. 10(b) corresponds to the most intense speckle blinking in the case of a pair of identical scatterers within each resolution volume. If the speckle blinking is less intense (e.g., for a larger amount of different-strength scatterers within a resolution volume), the speckle tracking for probably may become feasible. However, a priori the density and properties of subresolution scatterers are not known and the examples shown in Fig. 10 demonstrate that even at the very initial stage of speckle-pattern decorrelation, the notion of speckle tracking should be applied very carefully and its feasibility cannot be taken for granted even for small strains.
Furthermore, in real biological tissues containing a significant portion of liquids, local strain-induced microflows and such factors as Brownian motion of subresolution scatterers may very significantly modify speckle patterns even without regular macroscopic deformations and may additionally complicate speckle tracking. The proposed CS approach does not relied on the measurement of displacements and local strains and suggests an attractive alternative way for elasticity mapping. Certainly, the above mentioned factors causing decorrelation of OCT images (like stress-induced microflows and time-dependent decorrelation due to Brownian motion) also reduce the contrast of CS mapping and in such a context play a negative role. Nevertheless, feasibility of CS mapping should be possible by using a sufficiently fast OCT system that is able to acquire B-scans, for which the direct strain-produced decorrelation of speckles occurs earlier than the decorrelation due to the above-mentioned rheological effects unrelated to macroscopic deformations.
Experimental Demonstrations of Elastographic CS Mapping Based on Speckle-Level Decorrelation using Phantoms and an In Vivo Example
In this section, we present some examples of the CS approach realization using phantom samples made of silicone layers with different stiffnesses. In the structure of such optically turbid phantoms, the subresolution scatterers were distributed fairly uniformly. Thus, after speckle smoothing/averaging, the remaining image would not contain sufficiently pronounced macroscopic inhomogeneities that could allow for realization of CS mapping based on purely geometrical distortions of morphological features discussed in Sec. 3. In such a situation, the direct use of speckle-level decorrelation is practically the only feasible variant of the CS approach. A similar situation occurs for certain biological tissues (e.g., some mucous ones) that does not contain pronounced morphological features except for eventual layers that may be visible due to slightly different scattering properties of subresolution scatterers.
To illustrate such a case we prepared layered phantom samples made of two types of silicon. The thickness of the unstressed layers was about 0.7 mm each. A corresponding OCT image of a two-layer silicone phantom is shown in Fig. 11(a).
In this example, the softer silicone layer is located in the upper part of the image. The interface between these layers is well visible due to somewhat stronger scattering from the second stiffer layer. Near the image bottom, a less pronounced but still distinguishable other boundary corresponds to the interface with a third deeper layer of an optically turbid plastic with much higher stiffness. The presence of the interface with the much stiffer material at the bottom made it possible to directly estimate the contrast in the stiffness for the two upper silicone layers. This was done by gradually increasing the pressure produced on the sample by the OCT probe supported by a clamp, the axial displacement of which was controlled by a screw. The variation in the thicknesses of the compressed silicone layers could be directly measured in the OCT images. The so-estimated contrast in the stiffnesses of the two layers was about 3 times. Assuming that in the lateral directions, the material can expand fairly freely, the meaning of the so-determined “stiffness” (which has no rigorous mechanical definition) corresponds to the Young modulus . We recall that for soft rubber-like materials with the Poisson’s ratio close to 0.5 (which is also the case of soft biological tissues), the following relationship is valid: , where is the shear modulus entering Eq. (2). Consequently, in the case of the silicon phantoms, the contrast in the shear moduli of the two layers is equal to the contrast in their Young moduli.
For the used method of the sample compression by a finite-size OCT probe, the assumption of the free material expansion in lateral directions looks reasonable, although rigorously speaking it is not perfectly correct. The main inaccuracy originates from the fact that at the interface with the compressing OCT probe the compressed material does not freely slide. Consequently, the expansion in the lateral direction is not completely free and the effective stiffness of such a layer becomes somewhat greater than the Young modulus . If a layer of a rubber-like material is compressed without lateral sliding between two fairly rigid boundaries (e.g., glass plates) and if the thickness of the layer is much smaller than the lateral sizes of the compressing plates, the deformation of such a thin layer is closer to the case of forbidden displacements in the lateral directions. In such a situation, the effective resistance of the material to the change in the thickness is determined not by the Young modulus, but by a many times greater modulus corresponding to the deformation typical of a compressional elastic wave. Consequently, the effective layer stiffness in such an experimental configuration is mostly determined by the bulk modulus and can be much greater than the actual Young modulus. At a smaller degree, a similar phenomenon should also take place in the discussed configuration of the compression produced by a finite-size OCT probe. Indeed, it was experimentally verified that if the silicon softer layer in the above-described configuration was squeezed down to values significantly smaller than the OCT-probe radius (say, down to 100 to 200 μm), the effective (apparent) stiffness of that layer could significantly increase compared to the initial value for the same layer in a weakly compressed state. For a moderate stiffness contrast between the layers, this effect may become rather important in practice. Indeed, such a thin upper layer made of a material with twice smaller Young modulus may paradoxically look even stiffer than the underlying thicker layer made of a material with a larger Young modulus. The apparent higher stiffness of the actually softer material of the upper layer is absolutely real in such a configuration. In such a case, the tissue particles in the upper layer indeed experience smaller displacements and produce smaller image distortions than in the deeper located material despite the actually smaller stiffness of the upper layer. It should be clearly understood that this remark relates not only to the discussed CS approach but also to any other methods of estimating displacements and strains in compression-type schemes of elastography.
Returning to the examples shown in Fig. 11, it should be noted that the presented records were obtained using the same OCT scanner as for the morphological variant of the CS approach discussed in Sec. 3. However, in contrast to images with smoother and larger-scale morphological features, the speckle-level correlation processing requires more exact reproducibility of the scanning light-beam motion. Therefore, only the B-scans corresponding to the same direction of the scanning beam motion were used for speckle-level cross-correlation (i.e., B-scans of the same parity, either only odd or even ones, were compared with each other). In view of this, the repetition rate of the compared B-scans was twice reduced, down to . Examples of the so-obtained cross-correlation maps are shown in Fig. 11 [plots (c) through (f)] for several consequent B-scans recorded during gradual deformation of the sample. We point out that the sequence of real CS maps demonstrates very similar features to those observed for simulated maps that are shown in Fig. 9 and is based on the proposed simple model of speckle blinking caused by the material straining. Finally, Fig. 11(b) shows the elastographic CS map obtained by averaging 10 consequent CS maps calculated for five B-scans before and five after the chosen reference frame. The averaged CS map allows one to clearly see that the layers with different stiffness in the CS map exactly correspond to the two layers of the silicone phantom that can be distinguished in the raw intensity OCT image shown in Fig. 11(a).
The next demonstration of CS mapping is based on in vivo obtained images processed using the same speckle-level cross-correlation approach. As an example, the skin of a human finger pad was chosen, for which the outer corneous layer was expected to serve as a stiffer region overlaying the deeper softer tissue. Figure 12(a) shows a raw OCT image (B-scan with number 65 chosen as a reference one) and plots (c), (d), and (e) in Fig. 12 demonstrate three CS maps obtained for the immediate neighboring B-scans with numbers 64 and 66, and the more distant frame with number 70. The sequence of the images was recorded in the free-hand mode, so that the compression could not be well-controlled and was strongly affected by the natural tremor of muscles, heart beating, and other mechanical perturbations. Besides, the decorrelation was additionally enhanced by rheological processes in the tissue, the influence of which determined stricter requirements to the acquisition rate of B-scans. However, Fig. 12 demonstrates that even the available modest rate of 10.5 fps already allowed for free-hand implementation of the elastograpic CS mapping, although the correlation rapidly decreased with the distance of the compared frame from the reference one. The near-surface stiffer corneous layer of the finger-pad skin is quite clearly seen in the averaged CS map in panel (b), as well as in panels (c) and (d) in Fig. 12 corresponding to the frames immediately preceding and following the reference one (i.e., separated by about 100 ms). For more distant frames separated from the reference one by 300 to 400 ms, the cross-correlation coefficient is already reduced almost down to the background level. A several times higher acquisition rate would make the free-hand operation simpler by ensuring a larger number of frames demonstrating gradual decrease in the cross-correlation. However, even for the available modest acquisition rate, the CS approach realized at the speckle-level already allowed us to demonstrate a free-hand realization of the approach.
The presented results of numerical simulations and experimental examples obtained for both phantom samples and in vivo acquired images of biological tissues demonstrate that the proposed CS approach to elastographic mapping suggests an attractive simple alternative to other discussed variants of elastorgraphic mapping in OCT. This approach may be considered as a variant of compressional elastography in the sense that it is based on the observation of gradual decorrelation between compared images under the action of an external stress. The CS approach has several important advantages compared with earlier discussed methods:
1. Its very attractive feature is elimination of the conventionally discussed intermediate procedures of the displacement-field reconstruction.
2. Even more important is elimination of such a very error-sensitive procedure as numerical differentiation of the experimentally measured displacements. In fact, the issue of errors introduced by numerical differentiation of reconstructed displacements relates not only to cross-correlation techniques1,2,3 supplemented with one or another super-resolution technique4 but also to alternative phase/Doppler techniques.22,23 The problem is not only in the necessity to increase the sensitivity of measuring local displacement of the tissue particles [say, the vertical component ] but also in the appropriate increase in the accuracy of determining the initial location of the displaced scatterers, which is determined by the system resolution. Indeed, even for a rather high displacement sensitivity (say, up to nanometer scale), much higher uncertainty in the initial positions of the scatterers may introduce a rather strong uncertainty in determining the genuine local strains . Simple avoidance of differentiation by substituting the genuine local strain by the average ratio (like in Refs. 17 and 18) does not yet resolve the problem of clear mapping the local stiffness.
3. The operability range in terms of strains for the proposed CS approach is significantly wider than for methods based on direct determination of displacements. The evident reason for this difference is that for conventionally discussed methods, the similarity between the compared OCT signals from the reference and deformed samples should necessarily be sufficiently high to ensure the high-accuracy measurements. In contrast, the CS approach intrinsically implies that the image distortions (does not matter at the speckle level or at the level of morphological features) should be sufficiently strong to significantly reduce the cross-correlation between a significant portion of the compared images. Consequently, realization of the CS mapping intrinsically assumes creation of significantly larger strains.
4. The wider operability range makes the CS mapping very attractive for implementation in the genuine free-hand mode, which is the most interesting variant for practical applications. Besides, the possibility of using the OCT probe itself for the tissue straining also simplifies practical realization of this approach in comparison to some other interesting approaches that require the usage of additional technical means, e.g., for excitation of vibrations24 or surface waves9 in the studied tissue.
5. It is worthwhile to especially emphasize that the speckle-level correlation processing made it possible to successfully apply the CS approach even to the samples without pronounced morphological inhomogeneities, for which speckle suppression would remove practically all image features to which the correlation processing could be applied.
The main drawback of the elastograpic CS mapping described in the present paper is the fact that in such a simplest form, the obtained maps give only noncalibrated relative representation of the contrast in the stiffness of different tissue regions. In this context, it can be said that such qualitative elastographic CS maps may be compared with palpation of macroscopic tissue areas. For such a palpation, the very fact of detection of a stiffer region without quantification of its stiffness is already a very useful result. Then, other methods (e.g., biopsy) are used to determine whether the so-detected inclusion is malignant or not. Thus, implementation of CS mapping even in the simplest qualitative form should be very beneficial for extending possibilities of OCT scanners.
Besides, the above-described CS mapping can be combined with elements of DB approaches (e.g., for evaluation of strains in the already revealed areas with different stiffness). As is seen from the above-presented examples, often such regions are contoured by fairly well visible boundaries, for which the displacements and strain can be relatively easily evaluated (in contrast to “blind” estimation of the local strain distribution over the entire image). Such procedures can help to quantify the relative difference in the stiffness between the visualized regions with different properties. Finally, there remains a possibility to combine the OCT probe with a built-in force sensor, so that even absolute values of the tissue stiffness can be additionally evaluated. Such prospects are the subject of further development of the proposed CS-based elastography. Another interesting way of development of the CS approach is its extension for evaluating rheological properties of biological tissues. In the context of the present paper the decorrelation of rheological origin was considered as an interference factor complicating the elasticity mapping. However, intentional use of decorrelation due to rheological properties of biological tissues may open additional interesting diagnostic possibilities for clinicians.
The study was supported by RFBR grants Nos 13-02-00627 and 13-02-97131, and grant No 14.B25.31.0015 of the Russian Federation Government. Matveev L.A. acknowledges support of grant No MK-4826.2013.2 of the President of the Russian Federation for young scientists.