A Comparison of Peripheral Imaging Technologies for Bone and Muscle Quantification: A Review of Segmentation Techniques

Accepted Manuscript Title: A Comparison of Peripheral Imaging Technologies for Bone and Muscle Quantification: a Review of Segmentation Techniques Author: Andy Kin On Wong, Sarah Lynn Manske PII: DOI: Reference:

S1094-6950(18)30033-7 https://doi.org/10.1016/j.jocd.2018.04.001 JOCD 1036

To appear in:

Journal of Clinical Densitometry

Received date: Accepted date:

13-3-2018 11-4-2018

Please cite this article as: Andy Kin On Wong, Sarah Lynn Manske, A Comparison of Peripheral Imaging Technologies for Bone and Muscle Quantification: a Review of Segmentation Techniques, Journal of Clinical Densitometry (2018), https://doi.org/10.1016/j.jocd.2018.04.001. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

A Comparison of Peripheral Imaging Technologies for Bone and Muscle Quantification: a Review of Segmentation Techniques Andy Kin On Wong1,2,*, Sarah Lynn Manske3 1. Joint Department of Medical Imaging, Toronto General Research Institute, University Health Network, Toronto, ON, Canada 2. McMaster University, Department of Medicine, Faculty of Health Sciences, Hamilton, ON, Canada 3. Department of Radiology, McCaig Institute for Bone and Joint Health, Cumming School of Medicine, University of Calgary, Calgary, AB, Canada

Corresponding Author: * Andy Kin On Wong ([email protected]) 200 Elizabeth St. Toronto General Hospital, 7EN-238, Toronto, ON M5G 2C4 Tel: 905-399-0329 Fax: 905-521-1297

Running Title: Technical Review of Musculoskeletal Image Segmentation Keywords: High-resolution peripheral quantitative computed tomography, peripheral magnetic resonance imaging, bone architecture, intermuscular and intramuscular fat, image foresting, region-growing Conflicts of interests The authors have no conflicts of interests to declare.

1

Page 1 of 44

Abstract Musculoskeletal science has developed many overlapping branches, necessitating specialists from one area of focus to often require the expertise in others. In terms of imaging, this means obtaining a comprehensive illustration of bone, muscle, and fat tissues. There is currently a lack of a reliable resource for end users to learn about these tissues’ imaging and quantification techniques together. An improved understanding of these tissues has been an important progression towards better predicting disease outcomes and better elucidating their interaction with frailty, aging, and metabolic disorders. Over the last decade, there has been major advances into the image acquisition and segmentation of bone, muscle and fat features using computed tomography (CT), magnetic resonance imaging (MRI), and peripheral modules of these systems. Dedicated peripheral quantitative musculoskeletal imaging systems have paved way for mobile research units, lower cost clinical research facilities, and improved resolution per unit cost paid. The purpose of this review is to detail the segmentation techniques available for each of these peripheral CT and MRI modalities, and to describe advances in segmentation methods as applied to study longitudinal changes and treatment-related dynamics. While the peripheral CT units described herein have established feasible standardized protocols that users have adopted globally, there remain challenges in standardizing MRI protocols for bone and muscle imaging.

2

Page 2 of 44

1. Introduction A number of medical specialties and allied health professionals have converged on the need to quantify bone and muscle quality associated with aging, and also with disease states such as osteoporosis, metabolic bone diseases, sarcopenia, as well as changes secondary to diabetes, chronic kidney disease, cystic fibrosis, malignancy, and human immunodeficiency virus. The current standard of care for detecting bone loss is dual-energy X-ray absorptiometry (DXA), which generates two-dimensional images from which areal bone mineral “density” (aBMD) is derived. The same modality has also been used to compute full-body composition measures such as fat mass and lean mass, the appendicular portion of the latter is thought to represent mostly skeletal muscle. The need to better quantify the risk for fractures, frailty, and falls has led to motivation to better quantify bone and muscle composition and architecture. The current frontrunners for musculoskeletal imaging have been primarily peripheral quantitative computed tomography (pQCT), high-resolution (HR)-pQCT, and magnetic resonance imaging (MRI). The first two modalities have now been distributed widely across the globe, but have mainly been used for clinical research. Peripheral MRI has also been successful at serving the niche market of musculoskeletal investigations, and yields high signal-to-noise-ratio images even at lower magnetic field strength. The highly standardized imaging protocols for HR-pQCT have facilitated the generation of evidence for bone quality’s associations with fractures. However, it lags behind in terms of imaging muscle compared to pQCT, which can yield both muscle and bone quality at more proximal sites but operates at a lower voxel size (best 200 µm in plane x 2.0 mm slice thickness) compared to HR-pQCT (82 µm isotropic). MRI has been coined the leader in providing high-contrast images of skeletal muscle, but the precise quantification of detailed

3

Page 3 of 44

tissue features remains weak. Not to mention, the plurality of MR sequences, magnet strengths, and coil configurations limit standardization and comparability across studies.

With a wide gantry capable of accommodating the lower leg, both pQCT and pMRI have been catered towards enabling the quantification of appendicular muscle. With MRI, fat and muscle are more distinct, owing to the large differences in proton spin relaxation times. However, the small difference in k-edge between skeletal muscle and adipose tissue precludes the generation of sufficient contrast to observe inter- and intramuscular fat streaks on CT images. Even at low X-ray energies, the contrast between muscle and fat remains poor, and would require more Xrays or possibly a dual-energy approach, and consequently a higher effective dose on the patient due to the large degree of photoelectric absorption by tissues. Across these modalities, the segmentation of bone, muscle, and fat have led to varying degrees of confidence in measuring differences between disease states, changes over time, and association with clinical endpoints.

The present review will provide a comprehensive overview of the current methods and challenges associated with bone, muscle, and fat segmentation for images derived from these modalities. The reader will gain a better understanding of (i) the quality and tissue contrasts within images derived from each imaging modality, (ii) the current segmentation methods applied to quantify features, (iii) the challenges associated with segmentation and quantification, and (iv) new techniques in development for segmenting additional features within these images. The general audience for this review is researchers in the field of musculoskeletal imaging, but can also be extended to general imaging scientists and physicists who may see opportunities to

4

Page 4 of 44

apply the described image processing and segmentation methods to achieve superior precision and clinical value from other medical images. The structure of the review follows with an overview of segmentation methods applied to musculoskeletal imaging in Section 2. In Sections 3, 4 and 5, bone and muscle segmentation in HR-pQCT, pQCT, and MR images will be discussed, respectively. Section 6 will follow with a highlight of more advanced segmentation methods including ways to measure dynamics over time. Finally, Section 7 will end with the current state-of-the-art in segmentation techniques, and future directions in developing algorithms. 2. Overview of image segmentation algorithms Although bone contrasts well against muscle and fat in CT-based images, it appears as void signals in MR images, the geometry of which could be affected by magnetic susceptibility differences. Fat on the other hand is less distinct within CT images, and while clearly shown on MR images, they need to be separated from water signals. As a result, multiple segmentation methods need to be applied to musculoskeletal images in order to separate bone, muscle and fat features. Threshold-based segmentation involves defining one or more cut points for distinguishing an object from background or other tissues, and is useful when the image has uniform brightness and particularly when it is calibrated against a phantom. Thresholds can be selected based on reference or biological standards, or by identifying local minima on histograms (1). Unfortunately, due to the combination of signal inhomogeneities and the lack of a calibration standard across different MR sequences, it is difficult to generate phantom-calibrated MR images that can rely solely on threshold-based segmentation. Alternatively, edge detection marks contrasts in image pixel intensities and utilizes operators that rely on derivative gradient functions. Estimated edges are chained together and weaker edges that contain low gradient

5

Page 5 of 44

amplitudes are removed by thresholding. Performance is dependent on the clarity and gradient of the border between tissues, which could be affected by motion and resolution (1). By virtue of high contrast edges among bone, muscle, fat and cartilage, it is possible to apply edge detection procedures on both MR and CT images. However, when edges are ill-defined due to poorer resolution, motion, artefactual distortion, or simply the inherent characteristics of the tissue, edge detection often fails. Instead, cluster-based segmentation methods overcome the challenges of partial volume artifact and minimize the influence of poorer signal-to-noise ratio (SNR). Generally, this method alternates between segmenting an image and characterizing the properties of labelled classes to assign probabilities of class membership without the need for training data. One such method is the fuzzy c-means clustering algorithm, which hinges on the idea that pixels bear a given probability of belongingness to a certain cluster, with pixels on the edge sharing the lowest membership probability (2). Similarly, a cost-minimizing function is more robust against partial volume voxels. It expands seed points outwards, and at each subsequent pixel defines three properties: the position of its predecessor, the cost function of the path drawn and its seed root (3). The image is segmented into regions by virtue of how closely connected in costminimized function the pixels are. These above algorithms have been applied to pQCT and MR images of bone and muscle, and in some instances a combination of techniques is necessary. The power of these tools is tempered by the trouble encountered with image artifacts including motion, which could result in a false representation of bone, muscle, and fat structures. Below is a comprehensive review of specific segmentation methods and challenges associated with quantifying musculoskeletal features within each modality’s images.

3. HR-pQCT segmentation

6

Page 6 of 44

3.1 Periosteal and endosteal border identification & cortical and trabecular compartment separation Standard image processing for bone micro and macro-structural as well as densitometric outcomes is achieved using manufacturer software (Image Processing Language, Scanco Medical AG, Version 5.08b for XCTI and Version 5.42 for XCTII). Two widely utilized approaches have been applied by the manufacturer (Scanco Medical AG) to identify the periosteal border within HR-pQCT images. The original approach, implemented by default in the first-generation scanner (XCTI) requires user intervention to identify the bone of interest. Lines drawn by the user “snap” to the periosteal edge by applying a snake algorithm that minimizes spline energy based on large image gradients (4). The user manually corrects any errors, typically at locations of poorly mineralized cortical bone or cortical porosity. The second approach requires only identifying the bone of interest, and automatically crops the outer soft tissue from the entire bone region based on an auto-contour algorithm previously described by Buie et al (5). This segmentation method originally proposed by Buie et al (5), and refined by Burghardt et al. (6), relies on a series of morphological operations (dilation and erosion) to separate the periosteal and endosteal surfaces (Figure 1). The periosteal boundary of cortical bone is identified within the manually cropped region by applying a Laplace-Hamming filter to smooth (XCTI and XCTII) and enhance edges (XCTI only), followed by a global threshold to distinguish bone from background. A dilation followed by an erosion step is applied to remove Volkmann’s canals that penetrate the cortical bone. Marrow space is removed using a twodimensional connectivity criterion (6). The endosteal boundary is found by overlaying the periosteal mask on an inverted version of the original binary bone structure image. Trabeculae smaller than 2 voxels in thickness are removed using a 2D connectivity criterion and 3D 7

Page 7 of 44

thickness map, followed by dilation and morphological erosion to return boundaries to the endosteal surface (6). The total bone area and trabecular area masks constructed in 3D are combined to create cortical and trabecular masks. The cortical mask is then refined for analysis of cortical porosity (Ct.Po) to include only intracortical pores consistent with Haversian canals, and excludes pores resulting from transcortical foramen, erosions, pores connected to the marrow, and pores smaller than 5 voxels in size likely attributable to noise (Figure 1). This algorithm is implemented in the IPL platform (Scanco Medical), and has been provided on XCTI as an option, and is the default segmentation approach for XCTII.

An alternative approach developed by Zebaze et al (7) identifies the cortical compartment using StrAx1.0 software. HR-pQCT images are segmented into compact-appearing, transitional, and trabecular zones by identifying local bone edges using linear attenuation profiles. By consolidating information across multiple profiles defined at 0.1o intervals perpendicular to the periosteal surface from the centroid, final regions of interests (ROIs) are drawn (Figure 2). The outer edge of the periosteal surface is identified in each ROI as the first peak in the attenuation curve. The transition zone is marked by the first inflexion point after the first peak, and ends at the second inflexion point, which marks the beginning of the trabecular region (7). The Scanco manufacturer software defines a meta and inner trabecular region within the trabecular mask of which the former is defined by a 40% concentric peel from the endosteal boundary. This meta region could encompass some aspect of the transition zone described by Zebaze et al, though not comparable as the same 40% iterative peel is applied equally across individual slices for all individuals.

8

Page 8 of 44

3.2 Segmentation of bone from marrow and pore spaces For XCTI, bone is segmented from marrow and pores by applying a combined Laplacian Hamming filter to sharpen image features (cut-off frequency 0.4, epsilon 0.5). Greyscale voxels within the image are normalized followed by thresholding. In an effort to smooth the image and remove noise, freestanding clusters smaller than 25 voxels in area are removed (8) (Figure 3). For XCTII, bone is segmented from marrow by first applying a low-pass Gaussian filter (sigma 0.8, support 1.0) to smooth the images and to remove noise, and then a fixed global threshold (320 mg HA/cm3 for trabecular bone) is applied (9). From the cortical mask, information about cortical volumetric (v) BMD, cortical thickness (Ct.Th), and cortical area (Ct.Ar) can be obtained. One approach to calculate cortical porosity (Ct.Po) is by volume fraction: Ct.Po = Pore Volume/Total Volume x 100% (6). Alternatively, Ct.Po analysis can be conducted without segmentation of pores from the cortical bone; instead this approach relies on the assumption that fully mineralized adult human cortical bone is 1200 mgHA/cm3 [Zebaze 2013]. Cortical porosity is then calculated as 1 – Ct.BMD/1200 mgHA/cm3. These two approaches were compared against synchrotron radiation microCT (9 µm voxel size) on human cadaveric bone specimens (10). The StrAx density profile method although yielded a high correlation with synchrotron standards (R2=0.939), tended to overestimate porosity by 6.17% to 20.99% (95% confidence interval); whereas the threshold-based method showed a higher correlation (R2=0.977), underestimated the amount of cortical porosity by 2.60% to 10.76%, but performed more accurately for pore sizes larger than 140 µm in diameter (R2=0.983; accuracy within -4.88% to +2.45%) (10). The inaccuracy of the StrAx density profile method may be explained by the assumption of more homogeneous tissue mineral density, which proved

9

Page 9 of 44

not to be the case as demonstrated within compact bone from synchrotron radiation microCT. One further difference in the approaches to assess cortical porosity is that using the Scanco software, Ct.Po is typically calculated across the whole volume of interest whereas StrAx calculates Ct.Po in only the distal 40 slices. However, this latter method has only been applied to XCTI so far.

Trabecular bone analysis is performed differently on XCTI versus XCTII due to the differing spatial resolution. On both generations, trabecular vBMD and trabecular number (Tb.N) are measured directly from within the trabecular masks. On XCTI, bone volume fraction (BV/TV), trabecular separation (Tb.Sp), and trabecular thickness (Tb.Th) are then derived from trabecular vBMD and Tb.N using equations described by Laib et al (11). In contrast, there is sufficient spatial resolution to accurately assess BV/TV, Tb.Th and Tb.Sp in addition to other structural and densitometric measures directly using XCTII without making any derivations (9). Aside from these basic microstructural parameters, finite element analysis could be applied to the same masks to yield bone mechanical properties such as Young’s modulus, ultimate stress, stiffness, and failure load. 3.2 Muscle and fat – region growing and threshold-based segmentation More recent attempts to quantify soft tissue in distal tibia scans led to a soft tissue analysis (STA v1.0) package that quantifies muscle and myotendinous tissues (MT) (12). After removing skin and bones from downscaled images (82 to 164 µm) (manually in v1.0 and using thresholding in v.1.1), the algorithm applies a threshold-guided region growing technique whereby seeds are planted within MT and fat as defined by empirically-determined thresholds (MT: 34.22 to 194.32

10

Page 10 of 44

mg HA/cm3, Fat: -238.60 to -84.9 mg HA/cm3). Islands smaller than 80 voxels are removed (slice-wise in v1.0 and through a 3D thickness map in v1.1) prior to applying 20 iterations in which seed volumes are expanded by 1 voxel at a time. Overlapping MT and fat regions are set as undefined and ultimately distinguished by applying a global threshold of -25.28 mg HA/cm3, above and below which are classified as MT and fat, respectively. From the volumetric stack, MT volume (MTV, mm3), MT density (MTD, mg HA/cm3), and MT mean cross-sectional area (MTCSA, mm2) are computed. These measures have demonstrated reliability with test-retest measurement and have correlated with the same parameters obtained by pQCT but at the more proximal, 66% site location (12). The most recent update of this software has also improved upon muscle area segmentation as judged only preliminarily by visual inspection (preliminary data). No data have been published to validate this improvement to date. While the XCTII scanner is now capable of imaging more proximal regions containing true muscle, this soft tissue analysis technique has yet to be evaluated for accuracy and reliability at these locations. 4. pQCT segmentation 4.1 Trabecular density - iterative edge detection and concentric peeling Image processing for bone macro-structural and densitometric outcomes is achieved using manufacturer software (Stratec, Version 6.20c) by cropping a region surrounding the cortical shell. Several algorithms are applied to differentially examine trabecular and cortical features. To examine purely trabecular features, Contour mode 2 (iterative contour detection with automatic thresholding), and Peel mode 1 (concentric peel with fixed percentage) at 45%, are applied to segment bone from soft tissue, and to isolate a purely trabecular region, respectively (13). This peeled trabecular region is analogous to the “inner” region defined by 40% iterative peeling obtained on HR-pQCT (end of section 3.1). For iterative periosteal contour detection, the

11

Page 11 of 44

software automatically scans the image for a voxel that best represents cortical bone using a bone density histogram. This seed voxel, when identified with a second neighbouring voxel satisfying the same condition, then follows an iterative contour algorithm that proceeds in a clockwise fashion to identify nearest neighbouring voxels using an eight-voxel pattern. The iteration is repeated until the start voxel is met, closing a path of minimum curvature defining the periosteal contour. A 3x3 median filter is applied to this contour to eliminate noise along the cortical bone boundary (14). Out of the total bone area confined within this boundary, a concentric peel of 45% of pixels from the outer edge are discarded to give a region that consists purely of trabecular bone and marrow. These discarded pixels represent cortical bone, an intermediate region containing potentially trabecularized cortical bone, and a small medullary area containing trabecular bone. This assumption was shown to function optimally when the trabecular region was more regularly shaped (13). Since some of the truly trabecular region is also removed in this procedure, only trabecular bone density is computed from this mask. For a full trabecular envelope segmentation, Contour mode 1 and Peel mode 2 is used in conjunction with an outer threshold of 280 mg/cm3 to separate muscle from the cortex, and an inner threshold of 400 mg/cm3 to separate trabecular bone from cortical bone. While this method can yield a trabecular mask similar to HR-pQCT, it is sensitive to motion artifact and is affected by trabeculae in the transition zone, which could yield a smaller trabecular envelope and a lower trabecular vBMD. 4.2 Cortical density – thresholding-based segmentation To compute cortical bone density, a global threshold of 710 mg/cm3 is applied to the image using Cortical Separation Mode 1 (simple threshold algorithm), removing all pixels below this value both outside the cortex and within the trabecular envelope. This value was recommended by the manufacturer as it distinguishes between pixels that are completely filled by cortical bone from

12

Page 12 of 44

those partially filled by trabecular bone or soft tissue, two cases that would otherwise underestimate true cortical density (15). However, the simple threshold algorithm used by Cortical Separation Mode 1 underestimates cortical thickness, particularly when the cortex is lower density or when it is thin. Alternatively, Cortical Separation mode 31, an iterative contour search plus threshold algorithm can handle thinner cortices and includes more bone in the cortical segmentation, but its validity has not been evaluated. It was recommended previously that areal and density measurements cannot be obtained accurately using the same thresholds (16). While these cortical segmentation methods are feasible for density analysis, cortical thickness can be better computed using a dual threshold segmentation and contour search method that is available on OsteoQ software suite by Inglis Software Solutions (17). 4.3 Bone microstructure - region-growing and threshold-based segmentation The pQCT OsteoQ Software (Inglis Software Solutions, ON, Canada) has so far been the only tool that can yield apparent bone microstructural measures from pQCT images at 0.200 mm inplane pixel size images obtained at 10 mm/s for a higher SNR (See technical review for more information on acquisition (18)). The software uses a region-growing algorithm beginning with a seed point equal to the highest linear attenuation value within cortical bone, which is then expanded to one of eight nearest neighbouring voxels (Figure 4). The lower bounds of the grown region are guided by a density value equal to two SDs above the mean soft-tissue signal computed from soft tissue ROIs from a minimum of 15 images within the image set to be examined (for wrist images: typically 187.95-189.62 mg/cc; for ankle images: typically 206.07209.22 mg/cm3) (19). Voxels lying beyond the cortical bone are trimmed, and islands of unconnected bone within the slice are filled within the region. Bone regions are dichotomized from non-bone by applying a voxel-intensity histogram equalization procedure to increase the

13

Page 13 of 44

contrast in bone signal dynamic range, followed by a high-pass Laplacian (3x3 kernel) filter and a global threshold, derived using the Otsu method applied to the image’s bimodal histogram (19). For connectivity measurements of the trabecular bone, binary images of bone are further skeletonized using the parallel thinning algorithm reported by Zhang and Suen (20). The same microstructural parameters can be generated from these analyses as for HR-pQCT except for Ct.Po and finite element analysis-based mechanical properties. Instead, two-dimensional-based mechanical estimates are generated, including: polar strength-strain index (SSIp), and crosssectional moment of inertia (CSMI).

4.4 Muscle and fat – edge-detection and threshold-based segmentation pQCT scans of the mid-leg (typically 66%, 50% or 38%) and arm (typically 33%) can be further segmented to yield muscle properties. Subcutaneous fat is separated to yield muscle and bone areas together by selecting an ROI around the leg and applying threshold-guided (40 mg/cm3) iterative contour search (Contour mode 3 Peel mode 2) with a smoothing filter (F03F05F05). To subsequently isolate bone areas (tibia-fibula, or radius-ulna), a threshold-based contour search with an inner threshold of 400 mg/cm3 and outer threshold of 280 mg/cm3 is used (Contour mode 1 Peel mode 2). The total areas of these bone regions (including the medullary area) are subtracted from the bone and muscle regions from the first segmentation to provide muscle mass (mg), muscle cross-sectional area (MCSA, mm2), and the quotient of these measures, muscle density (mg/cm3) (21). Due to the lower density of fat versus muscle, muscle density can be a surrogate measure for the amount of fat within (intra-) and between (inter-) muscle groups. Unlike MRI, CT-based images fail to depict fascial boundaries; for example, the fascia found

14

Page 14 of 44

just beneath subcutaneous fat and around each bundle of muscle. Consequently, the segmented muscle regions at the mid-leg and mid-arm may exclude pockets of inter-muscular fat, particularly in individuals with greater fatty infiltration into muscle such as those with diabetes or who have been immobilized for long durations. When fat infiltrates between muscle groups and the entire muscle cross-section is no longer circular, edge detection from Stratec software often fails. Instead, a manual boundary must be drawn around the actual muscle region in order to circumvent this challenge (Figure 5). In the rare occasion, this method also fails and an even lower inner threshold must be used (< 0 mg/cm3), but this method is essentially a manual segmentation so requires the ROI boundaries to be traced accurately along the musclesubcutaneous fat perimeter. Alternatively, muscle can be segmented using a watershed algorithm (SliceOmatic, Tomovision, QC) to overcome challenges with failed segmentations due to fat infiltration. The differences between these two techniques have been described by Wong et al who showed that the watershed algorithm was able to rescue precision particularly in those who have a larger amount of fatty infiltration into muscle (21) (Figure 5).

15

Page 15 of 44

5. MRI segmentation 5.1 Bone microstructure – threshold- and map-based segmentation MRI bone imaging has been enabled by the fact that solid bone is devoid of proton signals. Several bone imaging sequences are available including: 3D gradient echo, spoiled gradient recalled echo (SPGR), balanced steady state free precision (bSSFP), fast imaging with steady state precession (FISP), and fast large-angle spin-echo (FLASE) imaging. A full description of acquisition techniques is available in the technical review article of this bone and muscle imaging series (18). The contrast created by differences in proton distributions between bone and surrounding soft tissue can be capitalized by segmentation algorithms. Since MRI signals are not proportional to density, they cannot be calibrated to density equivalent units and therefore densitometric measures are not typically reported from MR image analyses of bone. Indirect assessment of trabecular bone has been reported by measuring susceptibility differences between bone and marrow using T2* (22), phase distribution (R2=0.46 against pQCT vBMDtr) (23), or the ratio of signals from different echo times (R2*, R2 = 0.74 against QCT vBMDtr) (24). However, the associations with bone density are weaker with poorer resolution and vary with the orientation of trabeculae. Despite the lower resolution achievable with MRI compared to CT, direct measurement of trabecular structure remains possible without relying on densitometric measures by capitalizing on subvoxel processing and by interpolating bone (void signal) distributions. A number of research groups previously segmented bone from MR images. Majumdar (25), Newitt and colleagues (26) have mainly followed a threshold-based method whereby histograms of bone (IB) and marrow (IM) peak intensities are used to select a bone fraction-weighted

16

Page 16 of 44

threshold (IO) to segment trabeculae. The mean intercept length method following Parfitt’s model of parallel plates is then applied to compute structural parameters. Wehrli and colleague’s Virtual Bone Biopsy program (27) uses bone volume fraction (BVF) mapping to dichotomize images based on a BVF cut-off of 0.25. Subsequently, digital topological analysis (28) skeletonizes images and assigns rod and plate-like geometries, which are used to compute structural parameters using run-length analysis. Slight variations of each of the above groups’ methods have been observed across publications. Further modifications to these methods include zerofilling k-space, sinc interpolation (27), and signal inversion (29). However, many algorithms are met with the challenge of needing field-inhomogeneity correction (30). The largest limitation of these methods is the fact that individual scans need to be thresholded independently due to scan-specific variations in radiofrequency transmitter power, and receiver amplified gain, resulting in different reference grey level values. Inter-individual differences and changes over time may therefore be less comparable. Ouyang was successful in internally calibrating images using marrow, subcutaneous fat, tendon, bone and air to yield acceptable short- (CV: 2.00±1.30% for BV/TV) and long-term (CV: 4.92±2.56% for BV/TV) precision errors (31). However, in cases where subcutaneous fat and tendon are minimal to absent, this method may be unreliable. 5.2 Bone microstructure – cost-minimizing and cluster-based segmentation Other non-threshold-based methods such as fuzzy c-means clustering (See section 2) have later been applied by Han (32) and Folkesson (2), the latter of whom was also successful in measuring changes in trabecular structure with alendronate therapy using this method (33) but had poorer luck with threshold-based techniques. Wong et al applied a semi-automated image foresting

17

Page 17 of 44

transform algorithm described previously by Falcao et al (3,34) using OsteoQ software (Inglis Software Solutions, ON) (35). This cost-minimizing function relies on signal intensity, gradient and pixel relative position for separating trabecular and marrow regions from the cortical bone and surrounding tissues (Figure 6). The algorithm combines information on the number of bonetissue boundaries, pixels that form each contour, and a Euclidean distance map of square distance from a seed voxel planted in the trabecular region. Similar run-length procedures as applied by Wehrli are used to compute apparent structure, connectivity, and hole geometry.

5.3 Muscle and fat – cluster-based segmentation For macroscopic muscle and fat imaging, basic T1 and T2-weighted fast spin echo (FSE) sequences are already sufficient to yield high contrast between muscle and inter-/intra-muscular fat. Although spin echo techniques have been criticized for their of fat-water separation, specific sequences have been designed to either water suppress or fat-water separate (2- and 3-point Dixon, iterative decomposition of water and fat signals with echo asymmetry and least-squares estimation (IDEAL)). A review of muscle MR imaging techniques, please see the previous technical review article in this series (18). From fat-water separation methods, simple thresholding of fat fraction images (fat signal/fat + water signal) at 50% has been successful at segmenting fat from surrounding tissues. This method has previously been shown to be effective for visceral, subcutaneous and inter- and intra-muscular fat in children. In fact, full automation of this technique is possible (36). In another study of fat-water separated images in older adults, the automated Otsu method defined thresholds in subcutaneous, visceral and inter/intra-muscular fat separation from surrounding tissues (37). Here, the optimal threshold was

18

Page 18 of 44

determined by minimizing the weighted within-class variance and maximizing the betweenclass variance. Similar thresholding techniques could be applied to the Iterative Decomposition of water and fat with Echo Asymmetry and Least-square estimation (IDEAL) algorithm, which could separate fat and water from fast spin echo (FSE) or spoiled gradient recalled echo (SPGR) images while accounting for field inhomogeneity and chemical shift artifacts at inter-tissue boundaries. Wang et al were able to separate fat from abdominal tissues with IDEAL fat images using K-means (k=2) clustering (38). Csapo et al also applied fuzzy clustering but included 3 clusters (muscle, intermediate, fat), with a threshold selected between the intermediate and fat clusters (39).

The most common technique to handle this type of indistinct boundary employs fuzzy logic cmeans (FCM) algorithms that classify voxels (yi) into predefined tissue clusters (C) based on computed membership probabilities (uij), where the level of fuzziness between classes is defined by m, often set as 2. The algorithm operates as a function of vectors representing voxel clusters (V) that define its multi-dimensional image features, and generates a membership matrix (U). Tissue classification is achieved by minimizing the cost function:

Equation 1. Fuzzy c-means cluster cost function where N is the total number of pixels. The algorithm outputs fuzzy maps representing each class of tissues (background and bone; fat and marrow; muscle and skin). Subcutaneous, inter-/intramuscular fat, and marrow fat are isolated within these masks, but no elaborate processing is done

19

Page 19 of 44

to refine muscle and inter-/intra-muscular fat boundaries, or to improve muscle-fat edge detection. 5.4 Muscle and fat – Region-growing and watershed segmentation Two other groups applied manual segmentation of visceral and subcutaneous fat using Sliceomatic software (Tomovision, Magog, QC), which includes both watershed and regiongrowing algorithms (40,41). Unlike fat-water separated images, T1-weighted FSE images have shown challenges in segmentation due to partial volume artifacts, rendering separation of fat and muscle more challenging. Threshold-guided region growing segmentation has been successful in delineating inter- and intra-muscular fat and muscle boundaries. However, the analysis-reanalysis reproducibility suffers due to the inconsistent selection of thresholds to separate muscle and fat (42). Even small shifts in the threshold selected could impose a large penalty in the amount of fat segmented (Figure 7). The method remains reliant on visual inspection of what grey-level intensity is considered true fat and could greatly affect measurement of inter-individual differences, despite applying gamma correction to account for variations in screen displays. Others have used subcutaneous fat as an internal calibration to standardize thresholding (31), but this method is not robust against field inhomogeneity, which causes over or under segmentation in regions of the image that exhibit shifts in signal intensity distributions. For peripheral small-bore MR magnets, this challenge is not often observed.

20

Page 20 of 44

5.5 Muscle fascial boundary – Snake algorithm The definition of muscle boundaries remains to be a subjective challenge when the fascial border is ill-defined, even on MRI. Tracing the fascia manually could introduce a large amount of precision error within and between raters, particularly as lines within subcutaneous fat often appear like the fascial boundaries. Positano (43) and Orgiu (44) both previously used the gradient vector flow snake algorithm to trace around the fascia lata. This algorithm begins with a seed contour derived from a smoothed muscle mask that is evolved iteratively towards the muscle fascia by virtue of minimizing the snake function:

Equation 2. Snake function where internal energy (Eint) of the snake is governed by tension (α) and stiffness (β) coefficients, x is the relative position of voxels, c’ and c’’ are the first and second order derivatives, respectively, and k is a weighting coefficient to dictate how strongly the external forces act on the snake;

Equation 3. Internal energy of snake and where external energy acting on the snake (Eext) is guided by the original T1-weighted image with fat and muscle contrasts minimized (I(c(x)) by assigning all muscle voxels to equal the mean intensity of adipose tissue determined from fuzzy fat masks.

21

Page 21 of 44

Equation 4. External energy of snake Using the initial seed contour, the snake algorithm automatically minimizes the energy function to find the strongest low-intensity edge that is assumed to be the fascia lata. However, this method fails on some occasions where there is poorer contrast between fascia and muscle. Orgiu noted that in older obese women, this method required over 1000 iterations to converge on a sufficiently fitted snake; and moreover, its performance compared to manual segmentation of the fascia was worse (average symmetric distance (ASD) = 1.14 ± 0.41mm versus 0.81 ± 0.37 mm overall) in this subgroup. In contrast, young normal weight women showed an ASD of only half as large (0.50 ± 0.06 mm) as for older obese women (44). While the performance of the snake algorithm is mostly effective for identifying the fascial boundaries around the whole limb, it would be less applicable to segmenting individual muscle groups, particularly in the calf muscles of leaner individuals where the distinction between muscle groups is less clear.

5.6 Muscle region and muscle groups – Atlas-based segmentation To target the goal of segmenting individual muscle groups, Karlsson (45) and Thomas (46) assigned voxel membership based on probability maps generated using multi-image atlases of the full body. Sample atlases of a range of individuals were non-rigidly co-registered using a phase-based morphon method on target individuals to quantify muscle groups. Individual muscle groups throughout the body were defined manually and the training set of image atlases were scaled and applied to target image sets, and muscle groups were assigned based on probability of a voxel registering successfully to each atlas’ muscle group. The disadvantage of this technique

22

Page 22 of 44

is that probability thresholds for classification needs to be defined using empirical data, which could change depending on the population examined. However, by using a machine learning approach, this method has the potential to improve over time with larger sets of training images for diseased sub-populations. 6. Advanced Segmentation Techniques 6.1

Segmentation to measure changes in bone over time

One of the advantages of the higher spatial resolution of HR-pQCT compared to pQCT is the ability to assess subtle differences and changes in bone over time with repeated measurements. A tool to match baseline and follow-up scans to minimize error due to repositioning between scans, based on cross-sectional area, was implemented by the manufacturer with the first release of the HR-pQCT (8). However, 3D image registration, using rigid body transformations, has been shown to improve reproducibility of density, trabecular and cortical outcomes (47,48). As well as matching the global ROI between baseline and follow-up scans, 3D image registration can be used to more specifically detect local changes. To examine bone resorption at the endocortical surface, Nishiyama et al. (49) compared three methods of segmenting the cortex after registration of baseline and follow-up scans: 1) independent segmentation of cortex at baseline and follow-up, 2) overlapping cortex on both baseline and follow-up, and 3) baseline cortex indexed to follow-up (Figure 8). They found that change in cortical porosity using the baseline-indexed cortex was most strongly related to bone turnover markers during significant bone loss 1-year after kidney transplant. This study indicates that the baseline-indexed cortical segmentation method is highly sensitive to changes in bone resorption at the endocortical surface.

23

Page 23 of 44

To examine local regions of bone resorption and formation, Christen et al. (50) registered baseline and follow-up images using rigid body registration as described above, and subtracted the densities of individual voxels in the common region (Figure 9). Any change in density less than or greater than 225 mg HA/cm3 was classified as bone resorption or formation, respectively. To minimize the effects of noise on the outcomes, clusters with fewer than 30 voxels were excluded. They then used the images to run finite element models to demonstrate that regions of bone resorption corresponded with regions of low strain energy density, and regions of bone formation corresponded with regions of high strain energy density. This approach could be applied in other contexts where measures of bone turnover are desired. A similar approach using micro-CT has been validated in animal models (51-53), however validation against a reference standard approach such as dynamic histomorphometry in humans is difficult because of the need for tissue samples.

6.2

Segmentation by applying local adaptive techniques

Segmentation approaches for HR-pQCT presented to this point have applied fixed global thresholds to distinguish between bone and non-bone voxels. This approach may be limited in cases where tissue mineral density is affected (e.g., chronic kidney disease (54) and type 2 diabetes (55). Several adaptive threshold techniques have been proposed to overcome this limitation. Burghardt et al (56), applied a 3D Sobel-based algorithm that includes criteria to distinguish voxels at edges from non-edges. This algorithm produced results equivalent or better

24

Page 24 of 44

than the standard HR-pQCT method for trabecular microarchitecture measurements when compared with micro-CT. Valentinitsch et al (57) used a fully threshold-independent segmentation tool based on local texture features to segment trabecular and cortical bone. The algorithm found textural features for each voxel using a three-dimensional gray level co-occurrence matrix (58,59), the local structure tensor, and the Hessian matrix (60). They then used a machine learning technique, specifically a random forest classifier, to more accurately assign membership of each voxel to the trabecular or cortical region. The same investigators later developed a tool to extract trabecular bone features (61) by forming clusters of sets of voxels with similar textural patterns, using a Gaussian mixture model (62). Through training, each voxel was assigned to one of three different textural clusters, which were attributed “high”, “medium” and “low” quality bone. Using this technique, they found that post-menopausal women with fragility fractures had a significantly lower fraction of bone in the “high quality” cluster. However, it is unknown whether or not motion grade 4 (11%) and 3 (17%) images in this cohort were excluded from analyses, which could have contributed to variability in textural patterns measured. While these advanced techniques my provide more accurate segmentation results in some situations, they are not readily available or implemented in the commercially available software packages for HR-pQCT.

7. Challenges and Incentives 7.1 Challenges of customized vs. standardized segmentation approaches

25

Page 25 of 44

One of the challenges of interpreting results from studies using non-clinical analysis applications is that the segmentation approaches are not always standardized. pQCT and HRpQCT systems are each produced by an individual manufacturer, and thus the default settings and segmentation approaches are available to all users of these systems. For HR-pQCT, the manufacturer of each system recommends segmentation approaches to be used in standard scan settings, and are implemented by default in the manufacturer’s provided software. For pQCT, although there is no standard recommendation by the manufacturer, many users converge on applying the same contour and peel mode segmentation approaches and thresholds for computing overall density values. In such cases, consensus is largely driven by the research community. In contrast, MRI systems are produced by a number of manufacturers and field strengths vary across academic centres. More importantly, there are no standardized approaches for segmentation – neither from the manufacturers nor from the research community. Compared to the peripheral CT systems, the consensus on the best analysis practices for musculoskeletal MR image quantification remains poor. In fact, different implementations of similar algorithms, on different software platforms, for example, may still yield differing results. Standardization can aid in comparability of results across studies. However, as all of these systems are research tools or clinical tools with only research-driven analysis methods, the segmentation approaches are continuously under development and the ‘standard’ approach may change over time. Further, new research questions and/or application to different populations may require adaptations of segmentation approaches to ensure accurate results. Authors should carefully describe their segmentation algorithms, citing appropriate references when necessary. Readers should note possible differences that may result when different segmentation approaches are applied. In a challenging environment where research drives standardization efforts for analysis,

26

Page 26 of 44

governing bodies would be useful. For example, the International Society for Clinical Densitometry has set out guidelines for bone mineral density testing using DXA. They have also created position statements for use of peripheral CTs. However, there are not guidelines published by any agency regarding the analysis of these modalities. No governing body currently oversees any musculoskeletal MR imaging applications. 7.2 Compatibility of segmentation approaches across peripheral systems Due to differences in spatial and contrast resolution across systems, segmentation approaches are not always compatible across systems, although the principles often remain similar. The optimal approach always depends on the desired outcome; for example, measurement of trabecular thickness from an HR-pQCT image requires more accurate segmentation of individual trabeculae, whereas this level of precision is not needed when measuring apparent bone mineral density of the whole trabecular region. Differences in resolution also affect the segmentation approach chosen – in some cases at the benefit of improving accuracy; for example, XCTII uses a smaller voxel size (62 µm vs. 82 µm) and can capitalize on this improvement to perform direct measurement of trabecular thickness, separation, and volume fraction rather than deriving these features using modeled equations (11) based on trabecular bone mineral density and trabecular number. Further, due to the fundamental differences in MRI and CT, segmentation tools applied on CT images that largely depend on calibrated bone density values cannot be transferred to compute bone microstructure on MR images. Because of the major disadvantage that MRI cannot directly measure bone mineral density, and that structural properties are often confounded by chemical shift artifact and magnetic susceptibility differences, efforts in further developing techniques in MR bone imaging have been developing only slowly. Despite its poorer application to metabolic bone diseases, there remain a number of

27

Page 27 of 44

knee osteoarthritis investigations that continue to use MR imaging of subchondral bone since knee cartilage and meniscal properties are already being obtained from the same modality. Just like CT has adopted a dual-application approach to bone and muscle, MRI could similarly benefit from continuing a multi-tissue approach, particularly for joint anatomy investigations that could also encompass bone imaging. 7.3 Future Directions Segmentation approaches are continuously evolving and improving. Upgrades in computational power mean that more complex algorithms may be utilized with less concern for computational time. Fully-automated segmentation tools that do not require user input are increasingly being developed and are leading the field towards improved reliability and reduced analysis time. Although the trade-off in efficiency and precision is often a loss in accuracy, for the purpose of longitudinal studies quantifying change over time, this limitation may be acceptable. In fact, even inaccurate measures have the potential to predict outcomes reasonably well in the absence of any other alternative. A prime example of this notion is areal bone mineral density obtained by DXA – although the results are derived from essentially a two-dimensional shadow of bone, thus lacking volumetric accuracy, it can still predict fragility fractures (63,64). One aspect of segmentation approaches that will continually be a challenge is the availability of standardized software packages to conduct analysis. While some groups develop new algorithms and/or software packages for commercialization, providing an open-access platform has its advantage in more widespread application and greater likelihood in uptake, dissemination, and accelerating advancements. Implementation and usability for the non-technical user however remains a challenge. With the advent of machine learning, atlas-based algorithms could be combined with other segmentation algorithms to build a more efficient engine that could

28

Page 28 of 44

address multiple limitations from individual algorithms. One can imagine a system that can segment all possible musculoskeletal tissues within a given image set. However, to lead an effort in combining a number of these algorithms as described above, a pool of expertise will be required and the individual efforts in developing proprietary software may slow down the ability to attain such a goal. While segmentation approaches will need to evolve with improvements in image quality due to improvements in equipment and image acquisition procedures, collaborative efforts in streamlining algorithms for widespread musculoskeletal imaging applications will need to be maintained and renewed. In the new era of big data analysis, big image analysis will become of major priority. Having the ability to integrate these combinations of algorithms into health systems image databases will contribute towards feeding large data initiatives that could lead to improvements in systems-based health outcome prognostication and decision making.

29

Page 29 of 44

References 1.

Sharma N, Aggarwal LM 2010 Automated medical image segmentation techniques. Journal of medical physics / Association of Medical Physicists of India 35(1):3-14.

2.

Folkesson J, Carballido-Gamio J, Eckstein F, Link TM, Majumdar S 2010 Local bone enhancement fuzzy clustering for segmentation of MR trabecular bone images. Medical physics 37(1):295-302.

3.

Falcao AX, Stolfi J, de Alencar Lotufo R 2004 The image foresting transform: theory, algorithms, and applications. IEEE transactions on pattern analysis and machine intelligence 26(1):19-29.

4.

Kass M, Witkin A, Terzopoulos D 1988 Snakes: active contour models. Int J Comput Vis 1(4):10.

5.

Buie HR, Campbell GM, Klinck RJ, MacNeil JA, Boyd SK 2007 Automatic segmentation of cortical and trabecular compartments based on a dual threshold technique for in vivo micro-CT bone analysis. Bone 41(4):505-515.

6.

Burghardt AJ, Buie HR, Laib A, Majumdar S, Boyd SK 2010 Reproducibility of direct quantitative measures of cortical bone microarchitecture of the distal radius and tibia by HR-pQCT. Bone 47(3):519-528.

7.

Zebaze R, Ghasem-Zadeh A, Mbala A, Seeman E 2013 A new method of segmentation of compact-appearing, transitional and trabecular compartments and quantification of cortical porosity from high resolution peripheral quantitative computed tomographic images. Bone 54(1):8-20.

8.

Laib A, Hauselmann HJ, Ruegsegger P 1998 In vivo high resolution 3D-QCT of the human forearm. Technology and health care : official journal of the European Society for Engineering and Medicine 6(5-6):329-337.

9.

Manske SL, Zhu Y, Sandino C, Boyd SK 2015 Human trabecular bone microarchitecture can be assessed independently of density with second generation HR-pQCT. Bone 79:213-221.

10.

Jorgenson BL, Buie HR, McErlain DD, Sandino C, Boyd SK 2015 A comparison of methods for in vivo assessment of cortical porosity in the human appendicular skeleton. Bone 73:167-175.

11.

Laib A, Ruegsegger P 1999 Comparison of structure extraction methods for in vivo trabecular bone measurements. Computerized medical imaging and graphics : the official journal of the Computerized Medical Imaging Society 23(2):69-74.

12.

Erlandson MC, Wong AKO, Szabo E, Vilayphiou N, Zulliger MA, Adachi JD, Cheung AM 2017 Muscle and Myotendinous Tissue Properties at the Distal Tibia as Assessed by High-Resolution Peripheral Quantitative Computed Tomography. Journal of clinical 30

Page 30 of 44

densitometry : the official journal of the International Society for Clinical Densitometry 20(2):226-232. 13.

Ferretti JL 1999 Peripheral quantitative computed tomography (pQCT) for evaluating structural and mechanical properties of small bone. In: An YH, Draughn RA (eds.) Practical guide for mechanical testing of bone. CRC Press, Boca Raton, FL.

14.

Sievanen H, Koskue V, Rauhio A, Kannus P, Heinonen A, Vuori I 1998 Peripheral quantitative computed tomography in human long bones: evaluation of in vitro and in vivo precision. J Bone Miner Res 13(5):871-882.

15.

Stratec 2004 XCT Research Series Manual Software, vol. Version 5.50.

16.

Hangartner TN, Short DF 2007 Accurate quantification of width and density of bone structures by computed tomography. Medical physics 34(10):3777-3784.

17.

Calder K, Inglis D, MacIntyre NJ 2010 Comparison of pQCT-based measures of radial bone geometry and apparent trabecular bone structure using manufacturer and in-house developed algorithms. J Clin Densitom 13(4):433-440.

18.

Wong AK 2016 A comparison of peripheral imaging technologies for bone and muscle quantification: a technical review of image acquisition. Journal of musculoskeletal & neuronal interactions 16(4):265-282.

19.

Gordon CL, Webber CE, Adachi JD, Christoforou N 1996 In vivo assessment of trabecular bone structure at the distal radius from high-resolution computed tomography images. Phys Med Biol 41(3):495-508.

20.

Zhang T, Suen C 1984 A fast parallel algorithm for thinning digital patterns. Commun Assoc Comput Mach 27:236-239.

21.

Wong AK, Hummel K, Moore C, Beattie KA, Shaker S, Craven BC, Adachi JD, Papaioannou A, Giangregorio L 2014 Improving Reliability of pQCT-Derived Muscle Area and Density Measures Using a Watershed Algorithm for Muscle and Fat Segmentation. J Clin Densitom.

22.

Majumdar S, Genant HK 1992 In vivo relationship between marrow T2* and trabecular bone density determined with a chemical shift-selective asymmetric spin-echo sequence. Journal of magnetic resonance imaging : JMRI 2(2):209-219.

23.

Allein S, Mihalopoulou E, Luypaert R, Louis O, Panayiotakis G, Eisendrath H 2000 MR phase imaging to quantify bone volume fraction: computer simulations and in vivo measurements. Magnetic resonance imaging 18(3):275-279.

24.

Machann J, Schnatterbeck P, Raible A, Lutz O, Claussen CD, Schick F 2000 Magnetic resonance osteodensitometry in human heel bones: correlation with quantitative computed tomography using different measuring parameters. Investigative radiology 35(7):393-400.

31

Page 31 of 44

25.

Majumdar S, Genant HK, Grampp S, Newitt DC, Truong VH, Lin JC, Mathur A 1997 Correlation of trabecular bone structure with age, bone mineral density, and osteoporotic status: in vivo studies in the distal radius using high resolution magnetic resonance imaging. J Bone Miner Res 12(1):111-118.

26.

Newitt DC, van Rietbergen B, Majumdar S 2002 Processing and analysis of in vivo highresolution MR images of trabecular bone for longitudinal studies: reproducibility of structural measures and micro-finite element analysis derived mechanical properties. Osteoporos Int 13(4):278-287.

27.

Magland JF, Wehrli FW 2008 Trabecular bone structure analysis in the limited spatial resolution regime of in vivo MRI. Academic radiology 15(12):1482-1493.

28.

Gomberg BR, Saha PK, Song HK, Hwang SN, Wehrli FW 2000 Topological analysis of trabecular bone MR images. IEEE transactions on medical imaging 19(3):166-174.

29.

Bhagat YA, Rajapakse CS, Magland JF, Wald MJ, Song HK, Leonard MB, Wehrli FW 2011 On the significance of motion degradation in high-resolution 3D muMRI of trabecular bone. Academic radiology 18(10):1205-1216.

30.

Folkesson J, Krug R, Goldenstein J, Issever AS, Fang C, Link TM, Majumdar S 2009 Evaluation of correction methods for coil-induced intensity inhomogeneities and their influence on trabecular bone structure parameters from MR images. Medical physics 36(4):1267-1274.

31.

Ouyang X, Selby K, Lang P, Engelke K, Klifa C, Fan B, Zucconi F, Hottya G, Chen M, Majumdar S, Genant HK 1997 High resolution magnetic resonance imaging of the calcaneus: age-related changes in trabecular structure and comparison with dual X-ray absorptiometry measurements. Calcified tissue international 60(2):139-147.

32.

Han M, Chiba K, Banerjee S, Carballido-Gamio J, Krug R 2015 Variable flip angle threedimensional fast spin-echo sequence combined with outer volume suppression for imaging trabecular bone structure of the proximal femur. Journal of magnetic resonance imaging : JMRI 41(5):1300-1310.

33.

Folkesson J, Goldenstein J, Carballido-Gamio J, Kazakia G, Burghardt AJ, Rodriguez A, Krug R, de Papp AE, Link TM, Majumdar S 2011 Longitudinal evaluation of the effects of alendronate on MRI bone microarchitecture in postmenopausal osteopenic women. Bone 48(3):611-621.

34.

Falcao AX, Bergo FP 2004 Interactive volume segmentation with differential image foresting transforms. IEEE transactions on medical imaging 23(9):1100-1108.

35.

Wong AK, Beattie KA, Min KK, Webber CE, Gordon CL, Papaioannou A, Cheung AM, Adachi JD 2014 A Trimodality Comparison of Volumetric Bone Imaging Technologies. Part I: Short-term Precision and Validity. Journal of clinical densitometry : the official journal of the International Society for Clinical Densitometry.

32

Page 32 of 44

36.

Kullberg J, Karlsson AK, Stokland E, Svensson PA, Dahlgren J 2010 Adipose tissue distribution in children: automated quantification using water and fat MRI. Journal of magnetic resonance imaging : JMRI 32(1):204-210.

37.

Wald D, Teucher B, Dinkel J, Kaaks R, Delorme S, Boeing H, Seidensaal K, Meinzer HP, Heimann T 2012 Automatic quantification of subcutaneous and visceral adipose tissue from whole-body magnetic resonance images suitable for large cohort studies. Journal of magnetic resonance imaging : JMRI 36(6):1421-1434.

38.

Wang D, Shi L, Chu WC, Hu M, Tomlinson B, Huang WH, Wang T, Heng PA, Yeung DK, Ahuja AT 2015 Fully automatic and nonparametric quantification of adipose tissue in fat-water separation MR imaging. Medical & biological engineering & computing 53(11):1247-1254.

39.

Csapo R, Malis V, Sinha U, Du J, Sinha S 2014 Age-associated differences in triceps surae muscle composition and strength - an MRI-based cross-sectional comparison of contractile, adipose and connective tissue. BMC musculoskeletal disorders 15:209.

40.

Joshi AA, Hu HH, Leahy RM, Goran MI, Nayak KS 2013 Automatic intra-subject registration-based segmentation of abdominal fat from water-fat MRI. Journal of magnetic resonance imaging : JMRI 37(2):423-430.

41.

Davison MJ, Maly MR, Adachi JD, Noseworthy MD, Beattie KA 2016 Relationships between fatty infiltration in the thigh and calf in women with knee osteoarthritis. Aging clinical and experimental research.

42.

Beattie K, Davison MJ, Noseworthy M, Adachi JD, Maly MR 2016 Quantifying fat and lean muscle in the lower legs of women with knee osteoarthritis using two different MRI systems. Rheumatology international.

43.

Positano V, Christiansen T, Santarelli MF, Ringgaard S, Landini L, Gastaldelli A 2009 Accurate segmentation of subcutaneous and intermuscular adipose tissue from MR images of the thigh. Journal of magnetic resonance imaging : JMRI 29(3):677-684.

44.

Orgiu S, Lafortuna CL, Rastelli F, Cadioli M, Falini A, Rizzo G 2016 Automatic muscle and fat segmentation in the thigh from T1-Weighted MRI. Journal of magnetic resonance imaging : JMRI 43(3):601-610.

45.

Karlsson A, Rosander J, Romu T, Tallberg J, Gronqvist A, Borga M, Dahlqvist Leinhard O 2015 Automatic and quantitative assessment of regional muscle volume by multi-atlas segmentation using whole-body water-fat MRI. Journal of magnetic resonance imaging : JMRI 41(6):1558-1569.

46.

Thomas MS, Newman D, Leinhard OD, Kasmai B, Greenwood R, Malcolm PN, Karlsson A, Rosander J, Borga M, Toms AP 2014 Test-retest reliability of automated whole body and compartmental muscle volume measurements on a wide bore 3T MR system. European radiology 24(9):2279-2291.

33

Page 33 of 44

47.

MacNeil JA, Boyd SK 2008 Improved reproducibility of high-resolution peripheral quantitative computed tomography for measurement of bone quality. Medical engineering & physics 30(6):792-799.

48.

Ellouz R, Chapurlat R, van Rietbergen B, Christen P, Pialat JB, Boutroy S 2014 Challenges in longitudinal measurements with HR-pQCT: evaluation of a 3D registration method to improve bone microarchitecture and strength measurement reproducibility. Bone 63:147-157.

49.

Nishiyama KK, Pauchard Y, Nikkel LE, Iyer S, Zhang C, McMahon DJ, Cohen D, Boyd SK, Shane E, Nickolas TL 2015 Longitudinal HR-pQCT and image registration detects endocortical bone loss in kidney transplantation patients. Journal of bone and mineral research : the official journal of the American Society for Bone and Mineral Research 30(3):554-561.

50.

Christen P, Ito K, Ellouz R, Boutroy S, Sornay-Rendu E, Chapurlat RD, van Rietbergen B 2014 Bone remodelling in humans is load-driven but not lazy. Nature communications 5:4855.

51.

Schulte FA, Lambers FM, Kuhn G, Muller R 2011 In vivo micro-computed tomography allows direct three-dimensional quantification of both bone formation and bone resorption parameters using time-lapsed imaging. Bone 48(3):433-442.

52.

Birkhold AI, Razi H, Weinkamer R, Duda GN, Checa S, Willie BM 2015 Monitoring in vivo (re)modeling: a computational approach using 4D microCT data to quantify bone surface movements. Bone 75:210-221.

53.

de Bakker CM, Altman AR, Tseng WJ, Tribble MB, Li C, Chandra A, Qin L, Liu XS 2015 muCT-based, in vivo dynamic bone histomorphometry allows 3D evaluation of the early responses of bone resorption and formation to PTH and alendronate combination therapy. Bone 73:198-207.

54.

Malluche HH, Porter DS, Pienkowski D 2013 Evaluating bone quality in patients with chronic kidney disease. Nature reviews Nephrology 9(11):671-680.

55.

Leslie WD, Rubin MR, Schwartz AV, Kanis JA 2012 Type 2 diabetes and bone. Journal of bone and mineral research : the official journal of the American Society for Bone and Mineral Research 27(11):2231-2237.

56.

Burghardt AJ, Kazakia GJ, Majumdar S 2007 A local adaptive threshold strategy for high resolution peripheral quantitative computed tomography of trabecular bone. Annals of biomedical engineering 35(10):1678-1686.

57.

Valentinitsch A, Patsch JM, Deutschmann J, Schueller-Weidekamm C, Resch H, Kainberger F, Langs G 2012 Automated threshold-independent cortex segmentation by 3D-texture analysis of HR-pQCT scans. Bone 51(3):480-487.

34

Page 34 of 44

58.

Kovalev VA, Kruggel F, Gertz HJ, von Cramon DY 2001 Three-dimensional texture analysis of MRI brain datasets. IEEE transactions on medical imaging 20(5):424-433.

59.

Haralick RM, Shanmuagam K, Dinstein I 1973 Textural Features for Image Classification. IEEE Trans Sys Man Cyber SMC-3(6):610-621.

60.

Sato Y, Nakajima S, Shiraga N, Atsumi H, Yoshida S, Koller T, Gerig G, Kikinis R 1998 Three-dimensional multi-scale line filter for segmentation and visualization of curvilinear structures in medical images. Medical image analysis 2(2):143-168.

61.

Valentinitsch A, Patsch JM, Burghardt AJ, Link TM, Majumdar S, Fischer L, SchuellerWeidekamm C, Resch H, Kainberger F, Langs G 2013 Computational identification and quantification of trabecular microarchitecture classes by 3-D texture analysis-based clustering. Bone 54(1):133-140.

62.

Reynolds D 2009 Gaussian Mixture Models. In: Li SZ, Jain AK (eds.) Encyclopedia of Biometrics, 2nd Edition. Springer US, MA, USA, pp 827-832.

63.

Cummings SR, Black DM, Nevitt MC, Browner WS, Cauley JA, Genant HK, Mascioli SR, Scott JC, Seeley DG, Steiger P, et al. 1990 Appendicular bone density and age predict hip fracture in women. The Study of Osteoporotic Fractures Research Group. Jama 263(5):665-668.

64.

Nevitt MC, Johnell O, Black DM, Ensrud K, Genant HK, Cummings SR 1994 Bone mineral density predicts non-spine fractures in very elderly women. Study of Osteoporotic Fractures Research Group. Osteoporosis international : a journal established as result of cooperation between the European Foundation for Osteoporosis and the National Osteoporosis Foundation of the USA 4(6):325-331.

35

Page 35 of 44

Figure 1. A 3D reconstruction of HR-pQCT volume. Image was acquired at the distal tibia showing the cortical (grey), trabecular (green) and cortical porosity (red) segmentation generated using the “autocontour” approach implemented by Scanco.

36

Page 36 of 44

Figure 2. Image of the cortex taken from a scanning electron microscope. Image shows the definition of the compact-appearing cortex, the transitional zone and the trabecular compartments as implemented by StrAx1.0 software (7). Reproduced from Zebaze et al. Bone 2013 with permission from Elsevier.

37

Page 37 of 44

Figure 3. Images of a human cadveric radius scanned on the first (XCTI) and second (XCTII) generation HR-pQCT scanners. The top images show the acquired grey-scale images while the bottom figures show the 3D bone structure after generation-specific filters and thresholds are applied. Reproduced and adapted from Manske et al. 2015 (9) with permission from Elsevier.

38

Page 38 of 44

Figure 4. pQCT trabecular segmentation and apparent microstructure computation. Density threshold 2SD above soft tissue peak from signal histogram was used to guide regiongrowing algorithm to segment bone from soft tissue. Inner threshold to separate cortical from trabecular bone was interactively selected to define medullary area with manual edit capabilities. Trabecular bone was segmented from marrow and skeletonized using a parallel thinning algorithm reported by Zhang and Suen (20). Figure provided by Wong.

39

Page 39 of 44

A

B

Younger Adults

C

Older Adults

D

SCI Adults *

Figure 5. Threshold and watershed-guided segmentations of muscle from pQCT images. Top row: Example of failed threshold-based segmentation due to fatty infiltration (A & B). Alternatively contoured muscle region providing reasonable segmentation (C & D). Contour mode 3 and peel mode 2 with threshold of 40 mg/cm3 is used to separate subcutaneous fat from muscle. Bottom row: Comparison of watershed segmentations among younger, older adults and those with spinal cord injury (SCI), adapted with permission from Wong et al 2012 (21).

40

Page 40 of 44

Figure 6. pMRI trabecular bone segmentation and skeletonization. Cost-minimizing functions expanded medullary region from a seed point to define medullary area, the boundaries of which were adjusted manually where necessary. The image foresting transform was applied on a multi-scale dimension to search for tissue boundaries, contour pixels and expanded seeds based on Eucleadian distance maps. Figure provided by Wong.

41

Page 41 of 44

Figure 7. MRI segmentation of muscle and adiposity guided by threshold selection. Comparison of inter- and intramuscular fat segmentation by adjusting the upper threshold for fat by just a single unit. Larger amounts of intramuscular fat are apparent when threshold was lowered. Example of intermuscular fat is highlighted within dashed blue rectangle and intramuscular fat in green solid rectangle. Arrow points to fascia surrounding muscles. Analyses were performed using the region-growing algorithm within Sliceomatic software (Tomovision, Magog, QC). Figure provided by Wong.

42

Page 42 of 44

Figure 8. Three options for segmenting the cortex using baseline and follow-up images on HRpQCT. Baseline and follow-up images were co-registered and aligned. The cortex can be defined independently in both images (left) using only the overlapping cortical region on both images (centre), or the cortex defined at baseline transformed to the follow-up image (indexed, right) to capture cortical porosity changes at the endocortical surface. Figure courtesy of Nishiyama.

43

Page 43 of 44

Figure 9. Bone remodeling dynamics depicted in XCTI image volume. a) 3D reconstruction of the distal tibia from an XCTI image after filtering and segmentation. b) Sites of bone remodelling as determined after registration of baseline and follow-up images. Reproduced and adapted from Christen et al. 2014 with permission from Nature Publishing Group.

44

Page 44 of 44