Edge Aberration in MRI

Abstract

Direct application of conventional models for sub-voxel edge detection to modalities with intricate image formation like MRI results in systematic edge dislocations on a sub-voxel scale (edge aberration). By quantitative experimental analysis of this effect, a simple correction term can be calibrated, which is demonstrated to improve edge localization precision by a factor of 2.5 to surpass voxel size by 2 orders of magnitude.

Introduction

In raster image formation, a continuous so-called underlying function ρ is discretized into an image function f which maps voxels to (for example) gray scale values. Following a widely-used formulation, f is determined by f (x) = [A ∗ ρ] (x) at x being a voxel center (1) where, via the convolution ∗, A(x), called voxel aperture (VA), plays the role of a weighting and averaging function for ρ, thus substantially determining (small scale) image semantics-besides the nature of ρ and imaging geometry.

Due to the convolution with A, an underlying function ρ exhibiting an edge separating two regions within one voxel generally results in some intermediate gray value for this voxel. This so-called partial volume effect can be exploited to locate the edge in ρ more accurately than the grid resolution of f , which is referred to as sub-pixel or sub-voxel edge detection (Fig. 1).

However, such techniques assume an ideal case where A is a box function and its carrier is the voxel itself. Their application to differently formed images (nonbox A) such as MRI datasets therefore is in principle erroneous, but common practice. In this paper, we shall assess the systematic sub-voxel edge dislocations thereby introduced (edge aberration) and present a simple heuristic correction. 1.1

State of the Art

Initially, sub-pixel edge detection has been developed regardless of a specific application [ 1 ]. It proofed to be very useful and accurate in optic, notably aerial, imaging [ 2 ], where image formation is close to the ideal case.

In MRI, sub-voxel edge detection techniques are typically employed for segmentation of small structures like knee cartilage [ 3, 4 ] and the brain [ 5 ]. Also, edge enhancing filters may implicitly perform sub-voxel edge detection [ 6 ].

In these contributions, an MRI-specific VA is taken into account only in [ 5 ], albeit in the form of simple filtering of a super-sampled putative segmentation in a discrete domain. The authors of [ 6 ] explicitly assume a box-shaped VA.

Implicitly or explicitly, sub-voxel edge detectors first estimate a continuous approximation of ρ [ 1 ] or ∇ρ [ 3, 4 ], often consisting of a polynomial [ 1, 4 ]. The approximation is then analyzed and (up to) one edge point for each voxel is computed [ 1, 2, 4 ], or contours are optimized based on the approximation [ 3, 5 ].

All cited edge-point and contour techniques for volume data use a slice-byslice approach, in which edges are located in 2D slices and assembly of surfaces in the volume constitutes a post-processing step. 2

Materials and Methods

Applying our previously-reported edge detector [ 4 ] to MRI slices of an edge phantom, we evaluate the difference in measured edge locations when the phantom is moved by fractions of a voxel along one of the image axes. This leads us to a correction term for (1D) locations of single edge points, extrapolatable to 2D.

In what follows, we shall w.l.o.g. assume x = 0 at the center of the voxel under consideration and unit voxel width, denoted 1 vx. 2.1

Edge Detection Algorithm

The employed 2D edge detector [ 4 ] first computes the gradient image g(x) = |∇f (x)| via the Sobel operator. Next, g is locally approximated in every voxel’s 3 × 3 neighborhood by a second-order bivariate polynomial gˆ(x). Let gˆ(t) denote gˆ evaluated along the straight line given by the local image gradient, parameterized in t. Then if gˆ(t) has a local maximum inside the voxel under consideration, the location of the maximum is accepted as an edge point. Thus up to one edge point per voxel is obtained (Fig. 1).

Quanti cation of Systematic Errors and Error Correction

Practical considerations led us to an experimental setup in which the phantom is fixed in the scanner and the scanner’s field of view (FOV) is moved by small offsets ∆x. While this facilitates relative positioning of the phantom with respect to the FOV, the true position of the edge, x, remains unknown.

However, we can extract edge positions xm1 and xm2 from the image before and after the translation of the FOV. We define ∆xm = xm2 − xm1 and xm = (xm1 + xm2)/2. As the error in ∆xm is small and has a component systematic in xm, we can view the ratio β = ∆x/∆xm as a “function” of xm. We observe β(xm) = ∆x dx ∆xm ≈ dxm =⇒ x(xm) ≈

β(xm)dxm ∫ xm 0:5 vx which gives us a correction function for measured edge positions xm. We call α(xm) = x(xm) − xm the additive correction term to xm.

Using the normalized local image gradient (nx, ny) = ∇f (x, y)/ ∥∇f (x, y)∥, a measured 2D edge point (xm, ym) can be corrected via (x, y) = (xm + nxα((nx)xm), ym + nyα((ny)ym)) (2) (3) 2.3

Measurements

A test tube filled with Gd-based contrast agent served as edge phantom, producing a disk in MR slices. 6 volume datasets of 8 slices each were acquired with a T2-weighted fat-saturation spoiled spin-echo 2D sequence in a Siemens “TrioTim” 3 T whole-body scanner. In successive acquisitions, the FOV was translated by about 0.1 mm in-plane along x (row direction, phase encoded), while the phantom remained fixed. Voxel size was 0.625 × 0.625 × 3 mm3 at 3 mm slice thickness, TE = 11 ms, TR = 3.5 s, and SNR = 141. Ground truth data for ∆x was obtained from the center-of-mass shift of the noise-corrected slices.

First, 3 slices of each volume were analyzed. From each slice, 3 voxel rows perpendicular to the border of the tube were selected, where x-coordinates of the edge points from the left and right border of the tube were taken for each of the 6 FOV positions (Fig. 2). For each of the 3 voxel rows in each of the 3 slices, 5 edge offsets for ∆x ≈ 0.1 mm, 4 offsets for ∆x ≈ 0.2 mm, and so on down to 1 offset for ∆x ≈ 0.5 mm were thus measured, for left and right edges respectively.

Next, the measured offsets for ∆x ≈ 0.1 mm were used to approximate a correction function according to ((2)). Applying this, the remaining 5 slices were analyzed in a similar fashion to evaluate the correction technique. 3

Results

Plotting the ratio β against the measured intra-voxel position xm clearly reveals a systematic component that can be fitted by a polynomial. For this study, we find an additive correction term α(xm) = 21.3xm7 + 1.46xm6 − 10.9xm5 − 0.301xm4 + 1.13xm3 − 0.0830xm2 + 0.0677xm + 0.0168, with α and xm in vx. (Fig. 2.3.)

Edges to which the correction has been applied are shown in Fig. 2. The smoothing effect visible in the Fig. indicates that precision has been increased.

In numbers, we find an RMS absolute error in edge offset measurement for arbitrary displacements of 0.039 vx before and 0.016 vx after correction. The overall average signed offset error is 0.0022 vx and −0.0024 vx before and after correction, respectively. Fig. 4 shows the error distribution for different offsets. 4

Discussion

Although [ 5 ] addresses the MRI-specific VA in simulation of MR images, both MR physicists and the computer vision community seem to be unaware of the effective VA in MRI and the assumed VA in edge detection not coinciding-clearly a semantic mismatch of image formation and image analysis.

We have shown the resulting errors to be systematic and in principle amenable to correction by means of a simple additive term to the edge detector, calibrated by measurement of an edge phantom. Our focus was not on building a robust edge detector, and neither did we devise a generally-valid formula. Also, noise as well as the heuristic extrapolation from 1D to 2D demand further attention.

The proposed correction technique does not significantly change overally accuracy of about 1/500 vx because it is unaffected by edge aberration due to quite balanced over- and underestimation (Fig. 2.3). Promisingly however, precision of offset measurement is improved by a factor of 2.5, resulting in an estimated single edge localization precision of 0.016 vx/√2 = 0.011 vx =b 6.9 µm which is 2 orders of magnitude better than the image grid and much closer to ground-truth precision than before correction (Fig. 4). This is most relevant for applications relying on single edge points, e. g. for detection of subtle changes in local thickness.

References

1. Matalas I , Benjamin R , Kitney R. An edge detection technique using the facet model and parameterized relaxation labeling . IEEE Trans Pattern Anal Mach Intell . 1997 ; 19 ( 4 ): 193 - 205 .

2. Steger C. An unbiased detector of curvilinear structures . IEEE Trans Pattern Anal Mach Intell . 1998 ; 20 ( 2 ): 113 - 25 .

3. Cohen ZA , McCarthy DM , Kwak SD , et al. Knee cartilage topography, thickness, and contact areas from MRI: in-vitro calibration and in-vivo measurements . Osteoarthritis Cartilage . 1999 ; 7 ( 1 ): 95 - 109 .

4. Ko¨nig L , Groher M , Keil A , et al. Semi-automatic segmentation of the patellar cartilage in MRI . Proc BVM . 2007 ; p. 404 - 8 .

5. Bosc M , Heitz F , Armspach JP . Statistical atlas-based sub-voxel segmentation of 3D brain MRI . Proc ICIP . 2003 ;II: 1077 - 80 .

6. Salvado O , Hillenbrand CM , Wilson DL . Partial volume reduction by interpolation with reverse diffusion . Int J Biomed Imag . 2006 ; 2006 (92092): 1 - 13 .