1.

Introduction

Proper action potential conduction through the heart are required for normal heart function. Fiber orientation, the alignment of myocytes within the heart wall, greatly influences wavefront propagation. Therefore, abnormal fiber orientation increases the likelihood of abnormal cardiac rhythms, or arrhythmias.¹

Current methods for measuring fiber orientation range from histology^{2, 3} to diffusion tensor magnetic resonance imaging.⁴ It is known that fiber angle varies in a nearly linear fashion from the epicardium to the endocardium.^{2, 3} Karlon and colleagues developed an automated image analysis method for measuring fiber orientation from histology slices.⁵ The method uses intensity based gradients within an image to calculate fiber orientation and angular standard deviation.

Optical coherence tomography (OCT) is a nondestructive, noncontact imaging modality that can generate 3-D images rapidly with high spatial resolution.⁶ OCT detects minimally scattered light and utilizes intrinsic scattering contrast within samples to generate an image. OCT is capable of imaging to depths of 1 to $2\mathrm{mm}$ in cardiac tissue.^{7, 8, 9} In this letter, we present and validate an automated method for quantifying fiber orientation from structurally intact excised cardiac tissue preparations using intensity-based gradients applied to OCT en face images.

2.

Methods

Three-dimensional OCT image sets of the endocardial surface of an isolated right ventricular free wall (RVFW) preparation from a New Zealand white rabbit were used in this study. The protocol was approved by the Institutional Animal Care and Use Committee at Washington University. The experimental procedures have been previously described.¹⁰ The RVFW was dissected, stretched, and pinned epicardial side down onto a silicon disk. The sample was placed in 3.7% formaldehyde for one day and 20% sucrose solution for an additional two days. This dehydration step improves the visibility of the fibers under OCT imaging.

A microscope-based OCT system¹¹ was used to image volumes of the sample. The axial and lateral resolution of the system was approximately $10\mathrm{\mu }\mathrm{m}$ (in air). The three data sets presented in this work vary in structural complexity [Fig. 1a ]. Each volume was $4.5\mathrm{mm}\times 4.5\mathrm{mm}\times 3.16\mathrm{mm}$ , corresponding to a pixel size of $11.25\mathrm{\mu }\mathrm{m}\times 11.25\mathrm{\mu }\mathrm{m}\times 4.86\mathrm{\mu }\mathrm{m}$ . To calculate fiber orientation from 2-D OCT images, a modification of the intensity-based gradient algorithm described by Karlon ⁵ was used.

Fig. 1

3-D OCT data sets of the right ventricular free wall (RVFW). (a) 3-D OCT images of the RVFW, data sets (a1) RVFW1, (a2) RVFW2, and (a3) RVFW3. The arrow in A1 points to cross-sectional (en face) slices. (b) View of three en face slices in depth within the 3-D OCT data sets. The fiber structure is visible within the en face slices. The uneven sample surface (b1) and trabeculations [(b2) and (b3)] induce low-frequency changes and produce shadows within subsequent en face slices, inducing a gradient that may be greater than the gradient produced by the fibers.

Within en face OCT images [Fig. 1b], uneven sample surface topology (e.g., RVFW1) and shadows cast by trabeculations (e.g., RVFW2 and RVFW3) cause low spatial frequency changes in the background intensity. These artifacts introduce unwanted intensity gradients within the image. A 2-D second-order Butterworth high-pass filter was used under the assumption that the surface topology and shadowing artifacts have lower spatial frequency components compared to visible fiber structures. The high-pass filter was convolved with the en face OCT image to suppress variations in background intensity. The high-pass filtered en face OCT image was convolved with a Wiener filter for noise reduction.

Two-dimensional $3\times 3$ Sobel filters were used to estimate local gradients in the image. $G_x$ and $G_y$ are defined as the convolution of the horizontal and vertical Sobel filters with the 2-D en face filtered OCT images. For each pixel, the magnitude of the gradient, $G(i,j)=(G_x^2+G_y^2{)}^{1∕2}$ and the gradient direction, ${\theta }^{\prime }=\mathrm{atan}(G_y∕G_x)$ , was calculated.

Within a small local window of the image, $W$ , the dominant local direction of the gradient was computed by taking the maximum of the angular distribution function, $A_{\theta }^W$ , a function of $G(i,j)$ and ${\theta }^{\prime }(i,j)$ , as described by Karlon ⁵ The angular distribution function is a fit of a radial normal distribution to the distribution of angles within the local window.

The directions of the cardiac fibers were assigned as perpendicular to the direction of the dominant local gradients. Two criteria were used to reject invalid fiber orientation assignments. First, the algorithm identified windows with no tissue present by using a threshold on the average pixel intensity within the window. Second, to measure the confidence in fiber orientation measurements produced by the automated algorithm, a D’Agostino-Pearson ${\kappa }^2$ (normality) test was conducted on the angular distribution of calculated orientations within each window. High values of ${\kappa }^2$ indicate that the angular distribution function is not a normal distribution, and therefore, the confidence in the fiber orientation assignment in that window is low. Threshold values for ${\kappa }^2$ (0.02) and average intensity values within a window (80, 1.5 times the noise floor) were selected. The automated algorithm was implemented using the software package MATLAB 7.3.0.267 R2006b ( $\copyright$ 1984–2006, The Mathworks, Inc.).

In order to validate the method quantitatively, an investigator blinded to the results of the automated algorithm manually measured fiber orientation angles on the en face OCT images analyzed by the automated algorithm. En face images were analyzed in increments of $25\mathrm{\mu }\mathrm{m}$ in depth for all three data sets. Results from the automated algorithm were compared to manual measurements by analyzing the mean and standard deviation of fiber orientation assignments for each depth and orientation as a function of depth. This comparison was made using several window sizes, but the mean of the absolute difference was the lowest for a window size of $563\mathrm{\mu }\mathrm{m}\times 563\mathrm{\mu }\mathrm{m}$ (almost four myocytes in length, assuming that an adult cardiac myocyte has an average length of $150\mathrm{\mu }\mathrm{m}$ ). Results using a $563\mathrm{\mu }\mathrm{m}\times 563\mathrm{\mu }\mathrm{m}$ window were resampled to obtain 256 vectors per image. Under these settings, the 2-D fiber orientation algorithm runs in less than $7\mathrm{s}$ per image using a Workstation with an Intel processor running at $2.66\mathrm{GHz}$ and 2-GB SDRAM memory, $533\mathrm{MHz}$ with Windows.

3.

Results

Accurate fiber orientation measurements were obtained from all three data sets up to $500\mathrm{\mu }\mathrm{m}$ below the sample surface. Figure 2 shows example vector plots of fiber orientations overlaid on raw en face OCT images. The low-frequency background intensity change and shadows created by endocardial trabeculations are apparent within the raw OCT image. Preprocessing of the data effectively reduced the gradient contribution of these features, producing accurate measurements of fiber orientation for samples of varying structural complexity.

Fig. 2

Fiber orientation in the rabbit right ventricle. An automated algorithm enables quantification of the fiber orientation in the plane parallel to the wall surface for samples of varying structural complexities: (a) RVFW1, (b) RVFW2, and (c) RVFW3. Although there are gradients within the images due to uneven surfaces (a) and shadows created by trabeculations [(b) and (c)] the two-step filtering process effectively reduce their contribution to the gradient magnitude calculations.

It has been well established that normal fiber orientation varies nearly monotonically from the epicardium to the endocardium. Figure 3 shows fiber orientation measurements as a function of depth for representative $281\mathrm{\mu }\mathrm{m}\times 281\mathrm{\mu }\mathrm{m}\times 500\mathrm{\mu }\mathrm{m}$ volumes from each of the three data sets. In all three cases, there is a nearly linear change in fiber orientation with depth. Comparing the slopes for the RVFW1 volume shows that the automated measurements correspond very well to the manual measurements [Fig. 3a]. A quantitative comparison of fiber orientation assignments was conducted on an $1125\times 1125\times 500\mathrm{\mu }\mathrm{m}$ volume randomly selected within the RVFW1 data set [Fig. 3a]. Fiber orientation assignments made by the automated algorithm correlated very well to manual measurements, with a 1.002 slope and a 0.823 correlation coefficient.

Fig. 3

Fiber orientation as a function of depth. (a) Fiber orientation as a function of depth for RVFW1 comparing automated (solid squares) and manual (open squares) measurements for a $281\times 281\times 500\mathrm{\mu }\mathrm{m}$ volume. (a1) Comparison of fiber orientation measurements for an $1125\times 1125\times 500\mathrm{\mu }\mathrm{m}$ volume within RVFW1. Open squares—RVFW1 manual; solid squares—RVFW1 automated. The automated fiber orientation measurements agree very well with the manual fiber orientation measurements. (b) Fiber orientation as a function of depth for RVFW2 (solid diamonds) and RVFW3 (solid triangles) for a $281\times 281\times 500\mathrm{\mu }\mathrm{m}$ volume. The fiber orientation changes monotonically with depth.

4.

Discussion and Conclusion

In summary, an automated algorithm was developed and validated for quantifying cardiac fiber orientation within OCT image sets. During this study, the sample surface was made relatively flat through the use of pins, and dehydration agents enhanced the visibility of fibers in OCT image sets. These are procedures that cannot translate directly to in vivo imaging. Therefore, developing an OCT scanner with high axial and transverse resolution may alleviate the need for dehydration and fixation protocols. Maintaining a relatively constant transverse resolution with depth will increase the visibility of fibers in depth. In addition, there is a need to extend this algorithm to quantify fiber orientation in three dimensions to accommodate samples with irregular fiber orientation patterns and samples with curved surfaces. Quantifying fiber orientation in structurally intact excised preparations using OCT can potentially be used to correlate fiber orientation with electrical conduction properties (e.g., conduction velocity) and mechanical properties (e.g., strain analyses) of the myocardium in a variety of heart disease and arrhythmia models.