1 Introduction

A large number of studies have investigated the importance of abdominal fat to health, or used it as a co-variate to other factors. A correlation between abdominal fat, insulin resistance and other metabolic risk factors [5, 14, 16] has been shown. The longitudinal changes has also been studied, in e.g. children [18] and pre-/post menopausal women [4].

The gold standard for quantifying abdominal fat is Computed Tomography (CT) scanning [17]. The modality yields absolute values on the Hounsfiled Unit (HU) scale, which makes quantification relatively easy and the image acquisition is fast, which minimises effects of organ movement. However, CT has the drawback of ionising radiation which limits its use in healthy subjects. Another modality is Dual X-ray Absorptiometry (DEXA) [6], which is fast and cheap, but does not allow for differentiation of subcutaneous and visceral fat. Finally Magnetic Resonance Imaging (MRI) can be used. MRI allows not only for the differentiation of the subcutaneous and visceral fat, but also for splitting the subcutaneous fat into a superficial and deep compartment. The main drawbacks are that MRI is costly, scanning is time consuming and segmentation of the fat is difficult. The latter point is mainly due to the bias field image artefact, a consequence of magnetic field inhomogeneities, and the fact that the values in an anatomic MRI scan are not quantitative like CT.

For the task of segmentation several options exist. There is the option of manual segmentation, which is known to produce good results. However, it suffers from being time consuming and is subject to inter- and intra-observer variance. As an alternative automatic methods are attractive. Commercial software exist in the form of sliceOmatic (Tomovision, Inc., Magog, Canada) which handles the segmentation in a semi-supervised fashion. For the fully automatic case, several approaches have been tried. Zhou et al. [19], used chain-coding, a heuristic method that requires the scan to be centred on the L4-L5 vertebra. Leinhard et al. [11] proposed binary operations with some heuristics for the final segmentation. Knutsson et al. [9] used a morphon based approach, that – while yeilding good results – is very sparsely described with regards to implementation. Finally, Mosbech et al. [15] used active contours for automatic segmentation, and – as the only method – described a way to find the fascia of Scarpa, which divides the subcutaneous fat, using dynamic programming.

We present here a flexible package, available for download at GitHubFootnote 1, that can be used to segment a variety of MRI data. The initial step - bias field correction - uses the method presented by [10]. The segmentation is based on the graph-cut method proposed by Li et al. [12]. We have found a number of image based energies that can be combined to suit the problem at hand. A preliminary study of the method can be found here [2]. We compare the automatic segmentation against the gold standard CT [17], and obtain correlations and bias between the two modalities comparable to manual segmentation [8].

2 Materials and Methods

All code was written in MATLAB (The MathWorks Inc., Natick, Massachusetts, USA).

2.1 Bias Field Correction

The bias field is an image artefact present in all MRI-scans, and can be described as a low frequent noise over the image. To correct for the bias field, the method by Larsen et al. [10] was used, which assumes that the observed image originates from a generative model. The underlying ‘true’ and uncorrupted image is described using a mixture of Gaussians, and the bias field artifact, assumed to be multiplicative and smooth, is modeled using a linear combination of cubic b-splines with regularization on the bending energy. Model parameters are estimated using generalized expectation maximization (GEM). When model parameters have been estimated, the bias field can be computed and divided into the observed image in order to obtain its correction. The hyper parameters used as default in our package have been determined by testing on T1, T1 with water suppression and DIXON [13] sequences. For the first two on both the abdomen and the thigh, while the latter has only been tested on the abdomen.

2.2 Segmentation

Li et al. [12] formulated the segmentation of layers in volumetric data as a graph-problem. Representing the data as (x,y,z) where x and y are the horizontal plane and z the height, each terrain-like layer can only pass through each (x,y)-column once. For connectivity and to incorporate smoothness constraints, intra-column displacement of the layers is limited. The framework generalises easily for multiple layers, which can the be constrained with a min/max distance between them. The energy used (i.e. the vertex weights) are set depending on the problem. To solve the graph-cut problem, we use the solver implemented by Boykov and Kolmogorov [1].

To represent the MRI-scan in this fashion we start by unrolling the images using radial sampling, centred on the centroid of the slice. The abdomen is not round but elongated, and the subcutaneous layer of fat changes distance from the centroid rapidly in the lateral direction. We thus sample more densely in these regions, see Fig. 1a. This unrolling yields an image - Fig. 1b that can segmented by the framework of Li et al. [12].

Fig. 1.
figure 1

Radial sampling and the resulting unrolled image

Full size image

Subcutaneous Fat. The first step in the actual segmentation is the subcutaneous fat, defined by the inner and the outer surface. Those surfaces are found simultaneously, using a constraint imposed on the distance between the two layer. For thighs or other structures these settings can be varied. The segmentation employs a variety of image derived energies, that can be grouped in two cost-classes: gradient-based surface cost, and surface cost based on cumulative sums.

With the slices along the z-direction, the surfaces are lying in the angular direction \(\theta \) while the radial direction r takes the role of the height. The gradient-based surface energy is defined as

$$\begin{aligned} G = \text {sign}\left( \frac{\partial I}{\partial r}\right) \cdot \left( \left|\frac{\partial I}{\partial r} \right|+ \left|\frac{\partial I}{\partial \theta } \right|\right) \end{aligned}$$
(1)

where \(\frac{\partial I}{\partial r}\) and \(\frac{\partial I}{\partial \theta }\) are the gradients in the radial and angular direction.

We further define \(G_m\) the gradient applied to the median filtered image, \(G_g\) as applied to the Gaussian filtered image. For the outer surface we have \(G_o\) as the negative gradient and \(G_{os}\) as the negative gradient applied to the median and then Gaussian smoothed image. The costs for the inner surface are illustrated in Fig. 2, and the outer in Fig. 3.

Fig. 2.
figure 2

The three gradient-based contributions to the cost of the inner surface. The cost is constructed such that it has low (dark) values where the surface is to be detected

Full size image
Fig. 3.
figure 3

The two gradient-based contributions to the cost of the outer surface

Full size image

The surface cost \(S_1\) is derived by taking the G image and setting all values below zero to zero. A cumulative sum is then taken in the radial direction. For \(S_2\) a cumulative sum is taken on the unfolded image. We define \(S_3\) as the element wise product of \(S_1\) and \(S_2\) filtered with a median filter. For the outer layer we define a single energy \(S_o\), which is the same as \(S_1\) but applied to \(G_o\). As illustrated by Fig. 4 these costs regulates how deep the layer can go in the segmentation.

Fig. 4.
figure 4

The four surface costs based on cumulative sums

Full size image

The two surfaces can then be found by weighting these energies depending on the anatomical region and sequence used to record the images. A median filter is applied to the inner surface to get a smooth segmentation, see Fig. 5.

Fig. 5.
figure 5

Segmentation of the subcutaneous fat

Full size image

Fascia of Scarpa. The fascia of Scarpa divides the subcutaneous fat into a deep and superficial part. In metabolic studies it is of interest, as the deep compartment may behave as visceral fat [3]. An advantage of MRI compared to CT is that the fascia is often visible and a segmentation thus possible. The identification of the fascia of Scarpa is optional, and would e.g. not be used on the thigh. The fascia is very fickle and is not always visible, in which case the segmentation is not meaningful.

When the subcutaneous layer has been found, we need only consider that part of the image. We use 3 energies: \(FS_r\) the raw image normalised, \(FS_l\) the raw image filtered with a Laplacian Of Gaussian (LOG) filter, and the surface cost \(S_1\). These energies are then weighted depending on the anatomical position. Examples can be seen in Fig. 6 and more details can be found in the code.

Fig. 6.
figure 6

Examples of the segmentation of the fascia of Scarpa

Full size image

Spine Extraction. The spinal compartment needs to be removed for accurate results. The marrow in the vertebra has comparable intensities to fat and is thus classified as such if not removed. The identification and removal of the spinal compartment is optional and would e.g. not be used on the thigh. The method is based on a new unfolding centred on the spine, with heuristic constraints added on the image derived energies. See the code for more details. An example is shown in Fig. 7.

Fig. 7.
figure 7

Spinal compartment

Full size image

Visceral Fat. The visceral compartment is defined by the inner layer, except for the compartment around the spinal column. The spinal compartment is removed - if available - before the segmentation of the fat. A k-means clustering using 5 groups is run and the median value is used as a threshold. Everything with a higher intensity is defined as fat.

2.3 Test Subjects

Data were acquired in connection with Copenhagen Women Study http://cws.ku.dk. All participants received written and oral information about the study, including risks and discomforts associated with participation, before they gave their written consent to participate. The study was conducted according to the Helsinki Declaration and approved by the ethical committee in the capital region of Denmark, protocol nr. H-1-2012-150.

6 subjects - 3 overweight (ID 1 to 3) and 3 normal weight (ID 4 to 6) - were chosen randomly among their weight groups. All subjects had a CT-scan performed on the abdomen. Within a mean/maximum of 11.7/24 days the subjects also had a MRI scan of the abdomen. The CT were acquired on a Siemens SOMATOM Definition scanner (Siemens, Erlangen, Germany) with settings 100 kVp and 60 mAs - the slice thickness 2 mm was and the pixel size was 1.5234 mm  \(\times \)  1.5234 mm. The MRI scans were acquired on a Siemens Avanto 1.5 T scanner (Siemens, Erlangen, Germany) using a T1 sequence with water suppression - the slice thickness was 7.2 mm and the pixel size was 1.1719 mm  \(\times \)  1.1719 mm. Before the MRI scan 20 mg of Hyoscinbutylbromide (Buscopan) were administered intramuscular to minimise organ movements.

To compare consistent regions, slices were included from the first caudal slice were the iliac crest were not visible, moving in the cranial direction until the last slice before the liver becomes visible.

CT Segmentation. The CT scans were segmented using the method by Kim et al. [7] with the modification that diagonal directions were included in the step to find the visceral mask. All voxels with values between −150 HU and −50 HU were classified as fat. All slices were inspected for major errors.

MRI Segmentation. The images were bias field corrected using the following settings: ‘stepSize’ was set to the voxel size, the ‘desiredStepSize’ was set to 4, ‘numberOfGaussians’ 8, ‘smoothingDistance’ 25, ‘smoothingRegularization’ 5e-1, and ‘mask’ was an image of all positive voxels.

The subcutaneous fat was found with

$$\begin{aligned} \text {Energy}_\text {inner}&= 0.5\cdot G + 0.1\cdot G_g + 0.1\cdot G_m + 0.05\cdot S_1 + 0.25 \cdot S_3 \end{aligned}$$
(2)
$$\begin{aligned} \text {Energy}_\text {outer}&= 0.5\cdot G_{os} + 0.5\cdot S_o \end{aligned}$$
(3)

for the inner and outer layer respectively.

The spinal compartment was removed, and 5 clusters used for determination of the visceral fat threshold. All slices were inspected for major errors.

3 Results

The obtained volumes are given in Table 1 for the total volume i.e. the entire body, for the subcutaneous fat and for the visceral fat.

Table 1. Segmentation Volumes
Full size table

Further, as the total volume deviates between the two modalities we normalise the subcutaneous and visceral with the total volume, to get the two fat measures as a fractional fat content. These are given in Table 2. The mean difference±standard deviation between CT and MRI was −0.005±0.0122 for the subcutaneous fat and −0.013±0.0088 for the visceral fat.

Table 2. Relative fat Volumes
Full size table

A correlation and a Bland-Altman (Tukey mean-difference) plot for the relative volumes are shown in Fig. 8.

Fig. 8.
figure 8

Correlation and Bland-Altman plot for subcutaneous and visceral fat

Full size image

Subject 3 has the largest deviation from CT. All the slices from the subject are shown in Fig. 9.

Fig. 9.
figure 9

All slices for subject 3. Blue is the outer layer, green the fascia of Scarpa, red the inner layer of the subcutaneous fat, and yellow the visceral fat (Color figure online)

Full size image

4 Discussion

A software package for segmentation of abdominal fat in MRI has been implemented and validated by comparison with CT.

Even though the same anatomical regions have been selected, we do not obtain the same volume, with a difference of 5.7%. This might be due to partial volume effects, as the voxel size in the MRI is more than twice the size of what it is in the CT images. Further, there exist an uncertainty on the automatic segmentation of the CT data, and because the CT does not have the spinal compartment excluded we cannot expect the visceral data to match exactly. The main uncertainty, however, stems from subject placement and movement. Breath cycle affects the position of the abdominal organs and thus visceral fat, and to a lesser degree the subcutaneous fat. The placement on the table may push the buttocks more caudal or cranial, which affects the measurement of subcutaneous fat. Both resulting in a high variance, which can be minimised by following a strict protocol for subject placement. If not reduced in this way, a larger sample size is needed to detect any changes. The data used here were acquired with different aims, and do not follow the same protocol. With these caveats the mean difference of 1.9% for subcutaneous and −0.5% for visceral fat is acceptable.

When comparing the relative amount of fat present in each modality, we get very good correlations, better than with manual segmentation. For visceral fat manual/automatic r=0.89/r=0.97, subcutaneous fat manual/automatic r=0.92/r=0.99 [8]. This may partly be caused by the fewer number of subjects (6 vs. 27) included in our study. Our bias is also comparable to or better than for manual segmentation. For visceral fat we get −0.013 ± 0.0088 compared to −0.029 ± 0.019, and for subcutaneous fat −0.005±0.0122 compared to 0.004 ± 0.026 [8]. This is further confirmed by the Bland-Altman plots where all points are within the limits of agreement. The subject with the largest absolute error is shown. When inspected by an experienced medical doctor only minor errors were found.

Although not tried here, examples of the method applied to other sequences can be found in [2], where it also performed good as judged by inspections. The package is thus flexible across T1, T1 with water suppression and DIXON sequences. It should be noted that the data quality is very important to the quality of the segmentations. Artefacts from organ or breathing motion can easily lead to incorrect segmentations. Further, parameter tuning is often necessary, as the images can vary widely depending on the field strength of the MR-scanner and the sequence employed.

5 Conclusion

A flexible software package for automatic segmentation of abdominal fat has been implemented. Both the bias and variance is comparable to or smaller than for manual segmentation.