1574 IEEE TRANSACTIONS ON MEDICAL IMAGING, VOL. 24, NO. 12, DECEMBER 2005 A Novel Local Thresholding Algorithm for Trabecular Bone Volume Fraction Mapping in the Limited Spatial Resolution Regime of In Vivo MRI Branimir Vasilic* and Felix W. Wehrli Abstract—Recent advances in micro-magnetic resonance imaging have shown the possibility of in vivo assessment of trabecular bone architecture. However, the small feature size and relatively low signal-to-noise ratio (SNR) achievable in vivo cause the intensity histogram to be unimodal. The critical first step in the processing of these images is the extraction of bone volume fraction for each voxel. Here, we propose a local threshold algorithm (LTA) that determines the marrow intensity value in the neighborhood of each voxel based on nearest-neighbor statis- tics. Using the local marrow intensities we threshold the image and scale the intensities of voxels partially occupied by bone to produce a marrow volume fraction map of the trabecular bone region. We show that structural parameters derived with the LTA are highly correlated with those obtained with the previously published histogram deconvolution algorithm (HDA) and that the LTA is robust to image noise corruption. The LTA is found to correctly identify trabeculae with a significantly higher reliability than HDA. Finally, we demonstrate that the LTA is superior in preserving connectivity by showing for 75 in vivo images that the genus of the trabecular bone surface is always higher than when processed with the HDA. Index Terms—Image segmentation, local, threshold, trabecular bone. I. INTRODUCTION L ONG, load bearing bones adapt a structure that is optimal [1] to carry the maximum load at a given mass of bone. Long bones have an outer shell of cortical bone that encloses the marrow and the trabecular bone, an intricate network of struts and plates. Trabecular bone is concentrated mostly at the ends of long bones where it serves to distribute the applied loads with a minimum of bone material. Similarly, trabecular bone is predominant in the vertebrae. While a significant part of bone’s strength % comes from bone density alone, the microstructure of trabecular bone also plays a very important role. In a healthy person bone is dynamically resorbed and re- built [2], [3] through the process of remodeling for the purpose of repairing microdamage and adapting its structure to changing Manuscript received April 1, 2005; revised August 23, 2005. This work was supported in part by the National Institutes of Health (NIH) under Grant R01 AR49553, Grant R01 AR41443, and Grant T32 EB00814. The Associate Ed- itor responsible for coordinating the review of this paper and recommending its publication was W. Niessen. Asterisk indicates corresponding author. *B. Vasilic is with the Department of Radiology, the University of Pennsyl- vania Medical Center, 3400 Spruce Street, 1st Floor Founders, MRI Learning Center, Philadelphia, PA 19104 USA (e-mail: brana7@comcast.net). F. W. Wehrli is with the Department of Radiology, the University of Pennsyl- vania Medical Center, Philadelphia, PA 19104 USA. Digital Object Identifier 10.1109/TMI.2005.859192 loads. The balance of bone remodeling can be affected by var- ious diseases, the most common being osteoporosis [4], leading to a deterioration in bone competence. An ideal way to deter- mine bone quality would be an in vivo measurement capable of detection and quantification of changes in bone microstructure. Advances in magnetic resonance microscopy ( MRI) [5], [6] have made this noninvasive imaging modality a promissing tool that could potentially be used for longitudinal monitoring of tra- becular bone changes in patients. Since MRI pushes the limits of magnetic resonance imaging (MRI) in both resolution and signal-to-noise ratios (SNRs) [7], adequate image processing is crucial in extracting relevant in- formation from acquired images. Fundamental to the analysis of MR images of trabecular bone is to correctly classify voxels as bone or marrow. All subsequent analysis of microstructure depends on the proper voxel classification and it is, therefore, critical to preserve all available information on bone occupancy of voxels. This task is complicated by the following factors: a) Noise leads to misclassification of voxels and has to be dealt with, b) Partial volume effects are always present in MR im- ages of bone since the resolution is often too low to cover trabec- ulae with more than a single voxel, and finally c) Coil shading, especially important with surface coils, modulates signal inten- sities across the image which has to be accounted for by the algorithm used. Coil shading can be corrected before bone oc- cupancy analysis or can be incorporated into the algorithm used for bone volume fraction (BVF) mapping. Bone appears dark in MR images while marrow is bright due to the intense fat signal. However, noise and partial volume ef- fects make the image histograms unimodal, thus hiding the dis- tinction in intensities between bone and marrow. The absence of bimodality in the histogram causes simple thresholding methods to give very poor results so other, more sophisticated methods have to be employed. To this end we have used [8] a histogram deconvolution algorithm (HDA) that identifies voxels partially occupied by bone and iteratively removes the noise in those voxels producing a noiseless BVF map. However, we have ob- served that in low SNR images the HDA often fails to correctly identify thin or weakly connected trabeculae so we explored al- ternative solutions to the problem. In this work we introduce a statistical method which calcu- lates a local intensity threshold value that distinguishes voxels containing pure marrow from those partially occupied by bone. In the local threshold algorithm (LTA) we determine a condi- tional probability that relates the intensity of a voxel with the probability that it has a certain value of the discrete Laplacian. 0278-0062/$20.00 © 2005 IEEE