Subvoxel Processing: A New Method for Alleviating Partial Volume Blurring in MR Images of Trabecular Bone (2000)

Introduction MRI provides a noninvasive method for quantifying trabecular bone architecture for the purpose of assessing fracture risk in subjects with osteopenia (1, 2). However, as trabecular thickness is on the order of the achievable voxel size in vivo (100-150μm), more sophisticated processing and analysis techniques are required for accurate extraction of structural parameters. Linear interpolation is commonly applied to increase the apparent resolution of digital images. In one dimension, for example, the bone volume fraction (BVF) at a spatial location between the centers of two adjacent voxels would be computed as the average of the two voxels. Therefore, additional values calculated in this manner will never increase beyond the original values and, thus, contradict the notion that smaller voxels are more likely to contain larger fractions of bone. The subvoxel processing technique presented here reduces partial volume blurring effects, building on a previously developed method for estimating voxel BVF’s from images acquired in the limited spatial resolution regime (3). Theory Subvoxel processing rests on two assumptions: (1) smaller voxels are more likely to have high BVF and (2) bone is generally in close proximity to more bone. Each voxel is partitioned into eight subvoxels. The algorithm enforces strict conservation of bone mass, i.e. the total BVF in the original voxel is divided among the subvoxels. The precise amount allotted to a subvoxel is determined by the amount and location of bone outside the voxel but adjacent to the subvoxel. Thus, bone tends to be sequestered in the area of the voxel which is closest to other bone. The principle is illustrated in 2D with a 3×3 array of voxels in which the central voxel (BVF=0.2) has been partitioned into four subvoxels (Fig. 1). The weight for subvoxel a is the sum of the three adjacent voxel BVF’s (i.e., 0.3+0.4+0.6=1.3). Similarly, the weights for subvoxels b and c are 0+0.5+0.6=1.1 and zero, respectively. Subvoxels are assigned a nonzero weight only if bone is “attracted” to both outer sides by voxels containing bone. For example, subvoxel b has two voxels with BVF=0.6 and BVF=0.5 attracting bone to the right. The voxel with BVF=0.5 also attracts bone towards the bottom. In contrast, subvoxel d has no voxels attracting bone from the left and is thus assigned zero weight. Therefore, 1.3/(1.3+1.1)=0.54 and 1.1/(1.3+1.1)=0.46 of the bone is allotted to subvoxels a and b, respectively, yielding 0.54•0.2•4=0.43 and 0.46•0.2•4=0.37 (Fig. 1B). The factor of four is necessary since the volume of the subvoxel is four times smaller than that of the voxel. If voxel BVF is high, a subvoxel may be assigned an unrealistic BVF>1. In such a case, the subvoxel is assigned a BVF of 1 and the excess bone is assigned to the other subvoxels. The algorithm is easily applicable to all three spatial dimensions. Materials and Methods Subvoxel processing of an in vivo 3D MR image of the distal radius of a 60-year-old woman (Fig. 2) was compared to trilinear interpolation. The image was acquired at 1.5T using FLASE (4) at a voxel size of 137×137×350μm. A voxel BVF map was generated, after applying iterative deconvolution to estimate the noiseless histogram of voxel intensities (3). A cylindrical VOI (31×31×22 voxels) was arbitrarily chosen (Fig. 1B). Subvoxel processing and trilinear interpolation were applied to increase the apparent resolution in all three dimensions.