Spatial coordinate corrected motion tracking for optical coherence elastography

: We investigate a spatial coordinate correction (SCC) method to track motion with high accuracy for optical coherence elastography (OCE). Through SCC, we refer the displacement ﬁeld tracked by optical coherence tomography (OCT) in the loaded sample to individual material points deﬁned in a ﬁxed coordinate system. SCC allows OCE to perform spatially and temporally unambiguous tracking of displacement and enables accurate mechanical characterization of biological tissue for cancer diagnosis and tumor margin assessment. In this study, we validated the eﬀectiveness of motion tracking based on SCC using experimental OCE data obtained from ex vivo human breast tissues.


Introduction
Knowledge on tissue mechanical properties is highly significant in different fields of biomedicine [1]. For breast cancer management, breast lesion in macro-scale or meso-scale is known to have larger stiffness compared to normal breast tissue. The stiffness of breast tissue has been assessed through manual palpation and in elastography techniques including ultrasound elastography [2,3], magnetic resonance (MR) elastography [4] and optical coherence elastography (OCE) [5][6][7][8]. OCE is a functional extension of optical coherence tomography (OCT) that is a 3D microscopic imaging modality based on low coherence light microscopy [9]. Compared to other medical tomographic imaging technologies, OCT and OCE can perform imaging and sensing through a miniature probe to access deep breast tissue with minimal invasiveness [10][11][12]. In recent years, OCE technology has advanced significantly for breast cancer applications [13][14][15][16]. Particularly, our quantitative OCE (qOCE) technology features a fiber-optic probe (2 mm in diameter) and an inline force sensor, and has the unique advantage for in situ tissue characterization. OCE evaluates tissue mechanical properties by tracking spatially resolved displacement and quantifying local deformation. Therefore, spatially and temporally accurate displacement tracking is critical for effective breast tissue characterization using OCE.
In OCE measurement, Doppler analysis or speckle analysis is performed on OCT signals to quantify the deformation of the loaded sample [17][18][19][20]. Various motion tracking methods have been developed for OCE. Most of these methods assumed that the same OCT pixel observes the same material particle throughout the entire mechanical loading process. In other words, OCE has been conventionally tracking motion for individual OCT pixels rather than for individual material points. It has been noted that an OCT pixel in general sees different material points during the process of sample loading in OCE measurement [21,22]. However, this issue has not been discussed in detail before and has not been systematically addressed.
Under mechanical loading, material points that have different mechanical properties are anticipated to have different instantaneous velocity and strain rate. Moreover, the accuracy of motion tracking results depends on the magnitude of OCT signal. When a pixel sees highly scattering material points and generates OCT signal with large magnitude, the motion tracking based on phase analysis is accurate and unbiased [19,23]. When the same pixel sees a different group of material points and generates OCT signal with small magnitude, the motion tracking results can be highly biased. Without referring the results of motion tracking to individual material points, the observed temporal variation of mechanical behavior is inevitably convolved with the spatial heterogeneity of the sample. Therefore, motion tracking referred to OCT pixels can lead to reduced accuracy and compromised signal quality, particularly for mechanically and/or optically heterogeneous sample. On the other hand, tissue boundaries with discontinuity in optical and mechanical properties are of great significance for OCE's biomedical application. OCE directly measures depth resolved displacement and estimates local strain by taking differential operation on the displacement. Hence accurate assessment of displacement is critical for subsequent analysis of OCE data and there is a strong need for spatially and temporally accurate displacement tracking, to allow unambiguous OCE measurement for effective tissue characterization.
In this manuscript, we describe a spatial coordinate corrected (SCC) motion tracking method that refers the results of displacement tracking to individual material points defined in a fixed coordinate system. SCC is critical when large displacement is introduced to deform the sample. Introducing large displacement slowly in loaded sample eliminates the need for tissue preloading, ensures quasistatic characterization of elastic behavior of the sample, and generates sufficiently large displacement and force. To be specific, for material points at tissue boundary, SCC is needed if the displacement is larger than the depth interval of an OCT pixel (5µm for our OCT system). Considering a 0.5 mm imaging depth, SCC is relevant for a small strain of 1%. SCC motion tracking generates unambiguous motion field to enable accurate OCE assessment of mechanical heterogeneity for cancer detection and tumor margin assessment. Moreover, the motion tracking method described in this manuscript can be utilized to perform quantitative imaging of dynamic scenes using OCT in other applications such as blood flow visualization and quantification.
The manuscript is organized as the follows. First, we describe the principle of the SCC method for OCE motion tracking. Afterwards, details on OCE hardware used in this study and measurement protocol are presented. We then present experimental results, and finally summarize this study. The effectiveness of the SCC motion tracking method was demonstrated using experimental data obtained from compression OCE measurement of ex vivo human breast tissue samples.

Principle
We use the laboratory frame of reference to describe the spatial location of material points. A material point is defined as its original spatial location within the undeformed sample. The material point is anticipated to arrive at different spatial coordinates upon compression. To achieve spatially resolved motion tracking, we calculate the Doppler phase shift between complex OCT signals acquired with a time interval, and convert the phase shift to displacement and then velocity [24]. In the subsequent discussion, OCT signal indicates Doppler OCT signal unless otherwise stated. The need for SCC in OCE motion tracking is illustrated in Fig. 1. Consider an OCT pixel P that takes signal from the volume in the vicinity of coordinate r 0 =(x 0 ,y 0 ,z 0 ) in the laboratory frame of reference. Before the loading process starts ( Fig. 1(a)), signal of pixel P comes from the green cube. After the sample has been compressed, the point that is originally at r travels to r 0 . As a result, signal of pixel P comes from the blue cube instead ( Fig. 1(b)). As shown in Fig. 1, the same OCT pixel P observes different material points (within green and blue cubes) during the sample loading process, causing spatial-temporal ambiguities in motion tracking. The objective of SCC motion tracking is to establish a motion field for individual material points using the history of OCT measurement. For 3D motion field established during OCE measurement, the instantaneous velocity (v) at an arbitrary spatial coordinate (r=(x,y,z)) is a 3D vector with x, y and z component in Cartesian coordinate: v(r)=(v x (r),v y (r),v z (r))). In this study, we use a cylindrical indenter to perform uniaxial compression on the sample and acquired Ascan OCT signals over time from the center of the indented volume (r=(0,0,z)) ( Fig. 1(c) where the red line in the right inset represents the displacement magnitude increasing monotonically as depth). The measurement geometry has cylindrical symmetry and the non-zero velocity component at coordinate r=(0,0,z) is along z dimension (v(r)=(0,0,v z (r)))) [25]. Therefore, a scalar function of axial velocity v(z) is subsequently used to represent the motion field obtained from OCT measurement.
Assume the compression starts at time t = 0. OCT data (Ascan Doppler signal) is temporally sampled at an interval of δt and spatially sampled in z direction with an interval of δz. The instantaneous motion velocity is calculated using Doppler phase shift (Eq (1) where λ 0 indicates the central wavelength of the light source and δφ(z,t)=φ(z,t)-φ(z,t-δt). (1) A well-known issue for motion tracking using Eq. (1) is the artifact due to phase wrapping. However, in this study, we used inter-Ascan phase shift to calculate the axial velocity. Our OCT system operated at 1.3µm had an Ascan acquisition rate of 92kHz and could track axial speed without phase wrapping up to approximately 30 mm/s. This speed is much larger than the axial speed that we induced by compressing the sample. Hence phase wrapping will not be discussed in great detail. The principle of SCC method is illustrated in Fig. 2. At time t=δt ( Fig. 2(a)), axial velocities (v OCT (mδz,δt) m = 0, 1, 2, . . . ) are obtained directly from OCT signal at spatial locations δz, 2δz, 3δz, . . . for material points represented by the green circles. At t=δt ( Fig. 2(b), the volume sampled by the m th pixel in the Ascan comes from a material particle located at a different spatial location (l m ) in the undeformed sample (red circles). l m can be calculated using Eq (2), where ε m (δt) indicates local strain and ε m (δt) indicates instantaneous strain rate δt given ε m (0)=0). Therefore, v OCT (mδz,δt) directly obtained from the m th pixel of OCT Ascan in fact corresponds to the motion of material particle located at l m in the undeformed sample (Eq (3)). l m is generally non-uniform ( Fig. 2(b) and (c)). To generate motion field in a uniform sampling grid ( Fig. 2(d)), interpolation (operator H in Eq. (4)) is performed using the non-uniformly sampled measurements (v(l m (δt),δt)). Here v OCT (z,t) indicates motion field (spatially resolved instantaneous axial velocity) defined in the coordinate of OCT image, and v(z,t) indicates motion field defined for specific material points that are non-uniformly sampled. Ascans acquired subsequently will be processed in a similar manner. Assume the compression starts at t = 0 and ends at t = Nδt where δt is the Ascan sampling interval and N is an integer. We obtain Doppler OCT signal v OCT (mδz,Nδt) which originates from material points in the undeformed sample at coordinate shown in Eq. (5). Interpolation is performed to generate an instantaneous motion field in a uniform sampling grid (Eqs. (6) and (7)). Here the subscript 0 denotes the sample is considered as undeformed when t = 0. This procedure repeats until all the Ascans are converted to the actually undeformed coordinate at t = 0. For the material particle located at mδz at the beginning of the compression process, its overall displacement can be calculated using Eq. (8). Figure 3 summarizes the above described procedure for SCC motion tracking.
The above described SCC method relies largely on accurate assessment of local strain. Strain estimation is a differential measure, based on quantifying the change of displacement over space. Displacement obtained from Doppler OCT signal is affected by randomness in phase measurement. Therefore, local strain calculated using the finite difference of the displacement is extremely noisy. In this study, we choose to assess local strain using a weighted least square method that was proven to provide more accurate local strain quantification [19]. In the future, we may implement vector based method for more robust strain estimation [26 27].

OCE system, experiments, and signal processing methods
We used the quantitative optical coherence elastography (qOCE) system developed in our lab to acquire experimental data to validate the effectiveness of SCC motion tracking (Fig. 4). Details about the qOCE system can be found in our previous publications [28,29]. The system uses a spectral domain OCT engine at 1.3µm and uses a graphic processing unit (GPU) for real-time signal processing (wavelength to wavenumber calibration, fast Fourier transform, and Doppler analysis). The OCT system has a 7.5µm axial resolution and a 0.1 nm displacement sensitivity. For tissue mechanical characterization, the qOCE system features a miniature probe (cylindrical with 1.8 mm diameter) that has an integrated force sensor. The qOCE probe compresses the sample, and acquires optical signal from the sample that interferes with reference light to generate interferometric signal for depth resolved displacement tracking. Moreover, the probe deformation is interrogated by a common path OCT signal for force sensing. A spatial-division-multiplexing approach is implemented to simultaneously acquire signals for tissue displacement and force sensing.
We conducted OCE experiments on breast tissue specimens acquired from human subjects. The study was approved by the Rutgers University IRB. All patients gave written informed consent before enrollment in the study. We collected tissue specimens from patients who underwent breast cancer biopsy. Biopsy recommendations had been made based on the standard of care. Biopsy cores were obtained with a 14 gauge core biopsy needle by Dr. Hubbi, and were characterized using qOCE immediately following the biopsy procedure. Afterwards, the tissue was submitted in formalin for standard of care histopathologic analysis.
To acquire experimental OCE data, we placed the tissue specimen on the rigid surface of a motorized linear stage that moved along the axial direction, and brought the qOCE probe in contact with the surface of the tissue specimen. To compress the sample, we commanded the motor to translate for a specific depth at a constant speed. During the loading process, we obtained instantaneous motion velocity using Doppler OCT Ascans. Structural and Doppler OCT signals were displayed in real-time and saved for subsequent analysis.

Results
SCC is needed for OCE characterization of heterogeneous sample, because the same OCT pixel sees different material points during the loading process. We placed adipose human breast tissue on a rigid substrate, compressed the tissue and acquired OCT signals. Figure 5 shows OCT Ascans (structural) obtained before (black signal) and after (red signal) the compression. Clearly, the same OCT pixel sees different material points because of sample deformation. For example, the peak corresponding to the tissue-substrate interface (blue arrows in Fig. 5) shifted and was sampled by different OCT pixels. Results in Fig. 5 clearly show that OCE measurement not referred to individual material points is susceptible to misinterpretation and artifact. We further demonstrated that SCC accurately tracked the displacement of specific material points over time. We fabricated an elastic material phantom that contained TiO 2 for light scattering. We compressed the phantom, acquired OCE data, and extracted instantaneous motion speed for individual OCT pixels throughout the sample compression process. For the material particle sampled by the m th OCT pixel within an Ascan at t = 0, we determined its spatial coordinate l m (Nδt) at time t = Nδt in the laboratory frame of reference using Eq. (5). For different material points, their displacements were calculated: D 0 (mδz,Nδt)=l m (Nδt)-mδz, and the results are shown as solid curves in Fig. 6. For the same material particle, the displacement increased over time. For different material points, larger displacement was experienced at a larger depth. On the other hand, the sample was made to be optically and mechanically homogeneous. Hence different material points observed by the same OCT pixel during the compression process had the same pattern of motion. Because of the homogeneity, the displacement observed by the m th pixel can be calculated directly using OCT signal: D OCT (mδz,Nδt)= N n=1 [v OCT (mδz, nδt)δt], plotted as dashed curves in Fig. 6. Consistency between D 0 (mδz,Nδt) and D OCT (mδz,Nδt) in Fig. 6 suggests the effectiveness of SCC in accurate determination of spatial coordinate for material points. An ideal sample to validate the effectiveness of SCC is an engineered heterogeneous phantom with known mechanical and optical properties. However, it is very challenging to fabricate a phantom that resembles the mechanical and optical heterogeneity of breast tissue, and to the best of our knowledge such a phantom has not been fabricated before. To demonstrate the capability of SCC for accurate motion tracking, we performed OCE measurement on human breast tissues. Data was obtained from different tissue specimens (adipose breast tissue and nearby tissue of breast lesion) of the same patient. The patient was diagnosed as having benign breast tissue with stromal fibrosis.
Structural images of the adipose tissue specimen is shown in Fig. 7(a) (scale bar represents 0.5 mm). Figure 7(a) has a porous appearance that is typical for adipose breast tissue [30]. The cytoplasmic rims of the adipocytes produce strong OCT signal and vacuoles occupied by the lipid droplets generate dark areas. For OCE measurement, we commanded the motor to move 0.1 mm at a speed of 0.01 mm/s. The OCE probe tracked sample motion from the same lateral position. We tracked the instantaneous axial velocity with and without SCC throughout the compression process. Results in Fig. 7(b) and (c) suggest that SCC accurately registered the motion tracking results to specific material points. The vertical axis of Fig. 7(b) represents the spatial location (depth) of material points in the undeformed sample. Results in Fig. 7(b) indicate that the motion of the same material particle remained consistent over time. On the other hand, results in Fig. 7(c) were obtained without SCC. The vertical axis of Fig. 7(c) represents the spatial location (depth) of OCT pixels from which displacement/velocity was extracted. In Fig. 7(c), the motion observed at the same OCT pixel varied significantly over time, because different material points are seen by the same OCT pixel. With depth resolved instantaneous velocity (v 0 ) extracted for individual material points (Fig. 7(b)), we computed the temporally accumulated depth resolved displacement (D(z m , and show the result as the black curve in Fig. 7(d). Displacement obtained with SCC has multiple peaks, because the magnitude of structural OCT signal fluctuates above and below noise level and generated unbiased (peaks) and biased (valleys) velocity estimation. For comparison, we computed the temporally accumulated displacement using velocity directly seen by individual OCT pixels Fig. 7(d)). Structural characteristics for adipose breast tissue are no longer discernable, because motion from different parts of the sample got smeared in temporal averaging. Figure 7(e) shows a structural Ascan signal obtained from the OCE probe. The peaks in Fig. 7(e) are well aligned with peaks in the black curve of Fig. 7(d).
The capability of SCC in achieving spatially and temporally unambiguous motion tracking is further illustrated in Fig. 7(f). Notably, without SCC, it is challenging to scale the displacement field using the bulk compression of the sample, because the strain field is in general not uniform, as demonstrated in our previous study [31]. Figure 7(f) shows how the displacement accumulates over time during the compression process. Two depths were chosen (blue at a depth of 0.24 mm and red at a depth of 0.78 mm). The magnitude of displacements obtained with SCC (solid signals) in general increases overtime. Figure 7(f) also shows displacements obtained without SCC seen by OCT pixels at 0.24 mm (blue) and 0.78 mm (red) over the compression process (dashed curves). For the red dashed signal, it is clear that the slope of displacement against time (velocity) changes abruptly at approximately 4s. As shown in Fig. 7(c), the pixel at the depth of 0.78 mm saw material points generating large OCT signal at the beginning part of compression process. Therefore, the displacement tracking derived from this pixel was accurate at the beginning part of the compression process. When the sample was further deformed after t = 4s (indicated by the arrow in Fig. 7(f)), the pixel at the depth of 0.78 mm saw material points generating small OCT signal and provided incorrect estimation of motion velocity (approximately 0) and displacement (approximately constant). Without SCC, motion tracking result shown in the red dashed curve was incorrectly interpreted as time varying characteristics of the sample's mechanical behavior, while SCC motion tracking eliminated such ambiguity.
We also performed OCE characterization on the tissue of the breast lesion (stromal fibrosis) and show results obtained in Fig. 8. Figure 8(a) (scale bar represents 0.5 mm) shows the Bscan structural image with a homogeneous speckle pattern that is characteristic for uniformly distributed scatterers within diseased breast tissues. Figure 8(b) and (c) show depth resolved instantaneous velocity obtained with and without SCC. In Fig. 8(d), the accumulated displacement extracted with SCC is shown as the black curve while the accumulated displacement directly seen by individual OCT pixels is shown as the red curve. Figure 8(e) shows a structural Ascan signal obtained from the OCE probe. Motion tracking results in Fig. 8(d) show less pronounced difference as compared to results in Fig. 7(d), because of the spatial homogeneity of the optical properties for the breast lesion. Nevertheless, the displacement obtained with SCC (black curve) shows a significant peak (indicated by the green arrow in Fig. 8(d)) at the interface between tissue specimen and rigid surface on which the sample was placed (referred as interface I), and the displacement of this interface was accurately quantified (0.1 mm). On the other hand, results obtained without SCC (red curve) does not have this peak, because interface I was seen by different OCT pixels during the compression process (Fig. 8(c)), and temporal averaging diminished spatial resolution. In addition, Fig. 8(f) shows how the displacement accumulates over time during the compression process, at 0.24 mm (blue) and 0.78 mm (red) over the compression process, with (solid curves) and without (dashed curves) SCC. Solid and dashed signals in Fig. 8(f) are highly consistent, suggesting higher relevance of SCC in spatially heterogeneous sample (adipose tissue in Fig. 7(a)) compared to homogeneous sample (lesion tissue in Fig. 8(a)). Using the displacement obtained through SCC tracking of material points, we obtained the strain ( ) of the tissue specimens (adipose tissue shown in Fig. 7 and breast lesion tissue shown in Fig. 8) thought weighed least square estimation. The reaction forces (F) were also obtained from the force sensing channel of the OCE probe. The apparent Young's moduli (E) of the tissue specimen could be calculated: E = F/(πr 2 )/ and compared in Table 1. Although the apparent Young's moduli were calculated without considering measurement geometry, the results provides a relative measure to compare stiffness of different types of tissues. As shown in Table 1, stromal fibrosis breast tissue had larger stiffness compared to adipose tissue, which is consistent with results reported in literature [32,33]. However, the absolute values extracted from our measurement are different from some published values [34]. It is worth mentioning that OCE and many other elastography technologies provide an apparent assessment of tissue stiffness rather than Young's modulus in a rigorous sense. Young's modulus is the coefficient that linearly correlates the strain and stress. OCE and most other elastography techniques measure the total reaction force (F) to derive the stress (σ=F/A). However, the stress field established within the tissue is essentially unknown and largely depends on the boundary condition of the measurement. Therefore, the apparent Young's modulus (E apparent =(F/A)/ ) is different from the actual Young's modulus, and can vary significantly when obtained from different measurement platforms. Hence, a robust inverse algorithm considering the geometry of the measurement system is needed to achieve absolute measurement of tissue mechanical properties [35]. (c) instantaneous velocity extracted without SCC during the compression process; (d) depth resolved displacement accumulated during the compression process with (black) and without (red) SCC; (e) Ascan (structural) obtained from the site where OCE measurement was performed; (f) cumulative displacement of material points originally located at 0.24 mm (blue) and 0.78 mm (red) over the compression process (solid curves), and cumulative displacement seen by OCT pixels at 0.24 mm (blue) and 0.78 mm (red) over the compression process (dashed curves).

Conclusion and discussion
In this manuscript, we described and validated motion tracking based on spatial coordinate correction. Our results suggest that SCC allows spatially and temporally accurate tracking of material points within the sample under compression OCE measurement. To map the motion field established in a spatially heterogeneous sample, SCC is essential for motion tracking. Moreover, as shown in Fig. 7(b) and (f), Fig. 8(b) and (f), SCC accurately tracks temporal evolution of the displacement field. This is not only relevant for elastic mechanical characterization, but is also critical to investigate the viscoelastic behavior of the tissue specimen under sudden mechanical loading. We validated the SCC motion tracking method using one-dimensional data (Ascans). However, SCC can be easily generalized for 3D motion tracking. With history of instantaneous speed extracted for individual spatial locations in 3D, a series of 3D interpolations can be performed to refer the measurement to individual material points within the sample, given the 3D space is sampled sufficiently dense. The principle described in this manuscript can be used to correct spatial coordinate and track material points (rather than OCT pixels), as long as the history of motion at OCT pixel is obtained through Doppler analysis or speckle decorrelation analysis.
In SCC, we defined material points using their spatial coordinates when the sample is undeformed. For data obtained at an arbitrary OCT pixel during the compression process, we retrieved the material point from which the signal originated. Referring material points to their spatial location before compression is an arbitrary choice. Alternatively, material points can be defined using their spatial coordinates when the compression is done, or at any time point during the compression process. To suppress the noise due to decorrelation in motion tracking when the same OCT pixel sees different scatterers before and after sample compression, V. Y. Zaitsev et al developed a hybrid method of strain estimation [22]. In comparison, the objective of SCC is to determine the coordinate for each material particle during the compression process, using the history of observation made at a set of static observation coordinates.
The effectiveness of SCC motion tracking relies on accurate estimation of the trajectory of a material point over time, using the history of OCT data (Eq (2) and (5)). We determined the incremental displacement using estimated strain, and the strain was obtained using weight least squares involving multiple pixels within an Ascan assuming uniform strain distribution [19]. Instead of strain and strain rate estimated from a group of pixels, we have used the instantaneous velocity (v OCT (mδz,δt)) directly extracted from OCT data to assess the displacement of a material point within time interval δt. However, suboptimal results were obtained because v OCT (mδz,δt) was extremely noisy. The displacement extracted from an area with low scattering takes a value close to 0 and does not provide information on tissue properties. The accuracy of SCC for tissue underneath the low scattering layer may be affected, because of discontinuity in tracking the spatial coordinate for material points.