Introduction

Computer-aided diagnosis (CADx) or detection (CADe) systems are readily employed in life sciences, biomedicine and forensic science [1], to process medical image data for analytical purposes. For example, the computer-aided segmentation and classification of differences in organ morphology for patients with type 2 diabetes mellitus [2], polycystic liver disease [3], renal disease and obesity [4] has played a key role in supporting biomedical research. Moreover, recent literature reports that the measure of abnormal curvature growth can be used to stratify the severity and stages of Peyronie’s disease [5, 6]. The mean curvature of liver is used as an estimator for predicting the probability of suffering from and the severity of steatohepatitis [7]. Data about the structure and volume of psoas muscle can act as a predictor of outcomes for patients treated by chemotherapy for bladder cancer [8] and ovarian cancer [9].

Recent literature reports that the 3D brachial plexus reconstruction [10] determines the individual brachial plexus anatomy with maximum accuracy. Similar techniques are used to investigate craniocerebral trauma [11]. Furthermore, thin slice 3D reconstruction improves the detection of tumour margins in breast cancer patients compared to 2D CT images [12].

Considering the above successful applications of 3D rendering, the enhanced accuracy of 3D organ reconstruction in combination with morphological feature-based classification, may reveal previously unknown correlations between factors such as volume, curvature, spatial dimensions, anthropometry and health status. Although there are publications about varying correlations between organ-related morphological features and associated medical conditions, there are limitations to such studies, including a small patient dataset.

Recently, the relationship between nonalcoholic fatty liver disease (NAFLD) in type 2 diabetes has been investigated [13,14,15] together which increases the likelihood of developing complications of diabetes as well as augmenting the risk of more severe NAFLD [16]. Given the rising prevalence of type 2 diabetes [17, 18] and the broad clinical spectrum of the condition, there is a driving need for an improved understanding between morphology and function: this could provide objective measurements that help to establish progression, severity and remission.

Contributions

In order to investigate the possible relationship between organ morphology and anthropometry, this paper proposes a framework for automatic morphological feature extraction in 3D organ reconstructions in 3D radiological scans.

Current research literature, to the best of our knowledge, has not reported a generalisable framework for morphological feature extraction and correlative analysis. Thus, the proposed framework expands on the methodology presented in [19] by integrating an automatic organ segmentation approach [20], and furthermore, employs multiple case studies to analyse the relationship between computed morphological features. Moreover, the proposed tool presented in this paper can easily integrate into medical image analysis systems, and provide a potential “second opinion” before a final diagnosis. The analysis of organ curvature could provide quantitative guidance towards the assessment and stratification of a medical condition, thus enabling an improved therapy plan.

This novel approach produces detailed boundary tracing (contouring) for every protrusion and indentation as opposed to an approximate rating of the organ, as verified by an independent senior radiologist and radiographer. Two diverse MRI datasets are evaluated for the pancreas and liver, demonstrating the effective generalisability of the proposed segmentation approach. Computing morphological features following 3D organ reconstruction highlight a negative correlation between 3D curvature and volume with a statistical difference (p < 0.0001). The implementation will be available at https://github.com/med-seg/morph-feat-extract.

Section “Methods” details the methodology for automatic organ segmentation, reconstruction and morphological feature computation. Section “Results and discussion” presents and discusses two sets of computations. Section “Conclusions” provides a conclusion for the proposed framework, including probable future work.

Methods

The proposed framework consists of four stages: (1) an organ of interest is automatically segmented in an MRI volume; (2) an accurate 3D reconstruction of the organ of interest is generated; (3) rich organ features relating to the shape, surface (texture) and dimension are extracted; (4) these extracted features represent key imaging biomarkers that can be utilised for enhanced classification and study of the organ [21, 22]. A schematic representation of the sequence of stages in the proposed framework is provided in Fig. 1.

Fig. 1
figure 1

Overview of stages described in proposed framework. The organ is segmented, reconstructed and analysed to extract, classify and stratify morphological information

Full size image

Every stage in the framework can be treated as an independent process, potentially incorporated into another framework or pipeline. Replacing any stage in the framework will not significantly impact the complexity of a previous or later stage. The high modulatory enables a smooth substitution of different segmentation approaches without affecting the connections of subsequent stages. Moreover, the automatic segmentation stage is scalable depending on training data availability and capable of combining data augmentation. Both the 3D reconstruction and morphological-feature extraction stage can employ either manual or automatic segmentation using the target objects of interest.

Automatic organ segmentation

The automatic organ segmentation methodology exploits a training dataset of expert-led manually annotated MRI volumes, that is, a dataset of contours for each organ of interest. The methodology of this technique [20, 23] progresses through three main stages. Firstly, a digital contrast is applied to a test MRI volume to enhance the organ tissue against background classes of non-organ tissue; furthermore, a bounding region is generated to identify (and isolate) the major organ area for every 2D slice in the test volume. Secondly, 3D image segmentation is performed via continuous max-flow and min-cuts approach [24, 25] to generate distinct segments of detailed boundary tracing. In order to retain organ tissue and eliminate maximum surrounding tissue, the boundaries or contours of segmented tissue are detected using structured forests [26,27,28]. Measures of similarity between segmented tissue contours and the manually annotated contours are computed using modified Hausdorff distance [29] and structural similarity index [30, 31] to yield a rough segmentation outcome as a volume. Thirdly, any remaining tissue deemed as non-organ of interest is eliminated via morphological operations including area, curvature and position between distinct tissue in the segmented organ volume.

Reconstructing the organ

A 3D binary volume is generated from the segmented MRI volume. In attempt to reduce image noise, the reconstruction process employs a Gaussian smoothing algorithm that is applied to the 3D interpolated data. Using such a smoothing technique enhances image structures at varying scales of visualisation [32] as a pre-processing stage in computer vision. The equation of a 3D Gaussian function is

$$ G(x,y,z)=\frac{1}{2\pi\sigma^{2}}{\textstyle \exp\left( -\frac{x^{2}+y^{2}+z^{2}}{2\sigma^{2}}\right)} $$
(1)

where x, y and z are the distances from the origin in the horizontal axis, vertical axis and depth axis, respectively, and σ is the standard deviation of the Gaussian distribution. The values from the Gaussian distribution are used to build a convolution array that is applied to the original organ data. The Gaussian-based smoothed organ data is described as an isosurface by a set of points, for which the function represented by that data takes on a standard isovalue. In order to construct the isosurface, the marching cubes algorithm initially takes eight neighbour locations at a time, generates an imaginary cube and then determines the polygons needed to embody the part of the isosurface passing through this cube. Next, the individual polygons are fused into the desired surface by creating an index to an array of 256 possible polygon configurations in the imaginary cube. Each one of the eight scalar values represents 1-bit in an 8-bit integer: this bit is set to 1 if the value of the scalar is higher than the isovalue; otherwise, the bit is set to 0. Once all eight scalars have been analysed, the final value is the actual index to the polygon indices array. Each vertex of the created polygons is placed on the appropriate position along the edge of the imaginary cube by performing linear interpolation on the two scalar values connected by that edge.

Laplacian smoothing is then applied, producing the Laplacian of the rectangular 3D grid-based mesh. The implementation of a Laplacian algorithm returns an array of vertices coordinates: a new position is chosen for each vertex in the mesh using local information that relates to the position of neighbouring vertices. For each vertex, the smoothing operation is

$$ \overline{v}_{i}=\frac{1}{N}\sum\limits_{j=1}^{N}\overline{v}_{j} $$
(2)

where N is the number of adjacent vertices to node i, \(\overline {v}_{j}\) is the position of the j-th adjacent vertex and \(\overline {v}_{i}\) is the new position for node i. Once the 3D reconstruction process is complete, the 3D organ model can be visualised against a slice (2D image) from the original non-annotated MRI volume.

Computing the volume and curvature

The total organ surface area of each slice is calculated as its respective segmented pixel area, and the organ volume per section is calculated as the product of each organ slice area and the MRI section thickness.

The 3D curvature of the organ is calculated on a triangular mesh [34], in which the curvature of a surface describes the local shape of that surface. A regular surface S is represented by d(x,y) in Fig. 2, in which the point q lies on this surface. The orientation of S at q is the unit length normal, N. Also, R describes a regular curve on S, parameterised by β(a) = d(x(a),y(a)) and where a is the arc length of R and β(0) = q. Every curve that lies on the surface S at a given point of q ∈S has the same tangent line and the same normal curvatures. Consequently, the normal curvature is referred to along a given direction at q.

Fig. 2
figure 2

A point q is at a surface S with unit length normal N [33]. A regular curve R on the surface S passes through point q

Full size image

If q1 with unit length surface normal N1 is another different point on the surface very close to q, and t is the normalised projection of the vector \(\left (\boldsymbol {q}_{1}-\boldsymbol {q}\right )\) onto the tangent plane at q, the normal curvature, cn(t), along the tangent direction t can be approximated as

$$ c_{n}(\boldsymbol{t})=-\frac{<\boldsymbol{q}_{1}-\boldsymbol{q}, \boldsymbol{N}_{1}-\boldsymbol{N}>}{\Vert\boldsymbol{q}_{1}-\boldsymbol{q}\Vert^{2}} $$
(3)

where,

$$ \boldsymbol{t}=\frac{\left( \boldsymbol{q}_{1}-\boldsymbol{q}\right)-<\boldsymbol{q}_{1}-\boldsymbol{q}, \boldsymbol{N}>\boldsymbol{N}}{\Vert\left( \boldsymbol{q}_{1}-\boldsymbol{q}\right)-<\boldsymbol{q}_{1}-\boldsymbol{q}, \boldsymbol{N}>\boldsymbol{N}\Vert} $$
(4)

A triangular mesh, M = (P,C), is employed in the next stage and viewed as an approximation of an unknown smooth surface. Consider P as representing a set of data points, and therefore, C can be described as the connection of P to construct edges and faces in M. The motivation is to estimate the principal directions and curvatures at the vertices of M. Firstly, the normal vectors at the vertices of the triangular mesh are estimated. Let the triangular face in M be f. Secondly, each of these faces fi has a corresponding unit length normal vector \(\boldsymbol {N}_{f_{i}}\), and the triangular mesh is positioned so that the normal vectors point to the same side of the surface. The normal, N, at vertex q of M is a weighted average normal to the triangular faces adjacent to q, as

$$ N=\frac{{\sum}_{i=1}^{m}w_{i}\boldsymbol{N}_{f_{i}}}{\Vert{\sum}_{i=1}^{m}w_{i}\boldsymbol{N}_{f_{i}}\Vert} $$
(5)

where \(\boldsymbol {N}_{f_{i}}\) are the unit length normal to the triangles in the “one-ring” neighbourhood of q and wi is the weight, which is chosen based on the centre of the triangle face, fi. The number of points in a set of one-ring neighbour vertices of q is represented by m. For each neighbour qi of q, ti can be defined as the unit length projection of the vector \(\left (\boldsymbol {q}_{i}-\boldsymbol {q}\right )\) onto the tangent plane as

$$ \boldsymbol{t}_{i}=\frac{\left( \boldsymbol{q}_{i}-\boldsymbol{q}\right)-<\boldsymbol{q}_{i}-\boldsymbol{q}, \boldsymbol{N}>\boldsymbol{N}}{\Vert\left( \boldsymbol{q}_{i}-\boldsymbol{q}\right)-<\boldsymbol{q}_{i}-\boldsymbol{q}, \boldsymbol{N}>\boldsymbol{N}\Vert}\quad(i=1,\text{...}, m) $$
(6)

From here, it is possible to approximate the normal curvature cn(ti) as

$$ c_{n}(\boldsymbol{t}_{i})=-\frac{<\boldsymbol{q}_{i}-\boldsymbol{q}, \boldsymbol{N}_{i}-\boldsymbol{N}>}{<\boldsymbol{q}_{i}-\boldsymbol{q}, \boldsymbol{q}_{i}-\boldsymbol{q}>}\quad(i=1,\text{...},m) $$
(7)

Thus, the curvature of the organ is calculated on a triangular mesh that is based on local neighbourhood elements and vertices. The size of the neighbourhood can heavily affect results, and therefore increasing the neighbourhood size provides less sensitivity to noise [35], whereas a smaller neighbourhood size delivers better curvature estimates for less noisy data. From here, a global curvature Cg is computed as the mean value of all local curvatures, j = [1,...,Ng], such that

$$ C_{g} = \frac{1}{N_{g}}\sum\limits_{j=1}^{N_{g}}\boldsymbol{c_{j}} $$
(8)

The value of Cg can be used as a guide to estimate the roughness and smoothness of the organ surface.

Results and discussion

The proposed approach employs two datasets of T2-weighted (fat-suppressed) abdominal MRI image volumes, obtained using a Siemens Trio 3T scanner. The organ of interest, either the pancreas or liver, has been manually annotated by an expert-operator in every image volume. It is noted that the volunteers who underwent MRI scanning were aged over 18 and displayed early signs of type 2 diabetes. The pancreas dataset consists of 185 image volumes, split into 100/85 for training/test evaluation. Every image volume in the dataset consists of 80 slices with 2.0mm spacing, with each slice of spatial size 320 × 260 and 1.3128mm pixel interval in the axial and sagittal direction. The liver dataset consists of 50 image volumes, split into 20/30 for training and testing. Every image volume in the dataset consists of 370 slices with 3.0mm spacing, with each slice of spatial size 224 × 173 and 2.2321mm pixel interval in the axial and sagittal direction.

For the pancreas, the automatic segmentation method achieves a mean Dice Similarity coefficient (DSC) ± standard deviation (SD) of 81.14 ± 6.54%, and a mean Jaccard Index (JI) ± SD of 68.60 ± 7.68%. Analysing the liver segmentations produces a mean DSC ± SD of 92.17 ± 4.42%, and the mean JI ± SD is 85.73 ± 6.63%. The automatically segmented contours are elaborated to generate an image volume and an approximated surface using the techniques described in Sections “Introduction” and “Methods”. A number of 50 neighbouring points is chosen to calculate the curvature.

Organ contouring

Figures 3 and 4 display the segmentation contouring (green) against the ground-truth (red) contouring in eight distinct slices in eight different pancreata and livers, respectively. The top row in Fig. 3 displays relatively smoother contouring in comparison to the bottom row. Furthermore, as verified by a senior radiologist and a senior radiographer, the segmentation output captures detailed protrusions and dents of the pancreas’ boundaries in comparison to the ground-truth.

Fig. 3
figure 3

a, b, c and d display the segmentation (green) against ground-truth (red) contouring for eight slices in eight different pancreata, with a “smooth” surface. e, f, g and h the segmentation (green) against ground-truth (red) contouring for two slices in another two different pancreata, with a “ragged” surface

Full size image
Fig. 4
figure 4

a–h display the segmentation (green) against ground-truth (red) contouring in eight slices from eight different livers

Full size image

Aiming to investigate and thus derive a possible organ classification model, the volume and 3D global curvature are calculated for each MRI scan (case). The extraction and potential classification of parameters, relating to shape and texture, can help to tackle fundamental clinical questions about organ changes associated with obesity and type 2 diabetes mellitus, which are coexisting conditions frequently associated with NAFLD disease. Furthermore, since NAFLD is an increasingly recognised condition that could progress to end-stage liver disease [36] improved image guidance systems could potentially help to assist final diagnosis.

Figures 5 and 6 illustrate six diverse volumetric pancreas and liver reconstructions, respectively, with the segmentation (green) overlapping the ground-truth (red). Notice the visibly high variation in pancreatic structure and curvature, which is also reflected in a higher curvature standard deviation (6.83 × 10− 3) compared to the standard curvature deviation for the liver dataset (3.62 × 10− 3) (Figs. 7 and 8).

Fig. 5
figure 5

A list of six pancreas 3D reconstructions in six different MRI scans (volumes). The segmentation outcome (green) overlaps the ground-truth (red)

Full size image
Fig. 6
figure 6

A list of six liver 3D reconstructions in six different MRI scans (volumes). The segmentation outcome (green) overlaps the ground-truth (red)

Full size image
Fig. 7
figure 7

Box plots for (a) volume and (b) 3D curvature for 85 pancreas image volumes (cases)

Full size image
Fig. 8
figure 8

Box plots for (a) volume and (b) 3D curvature for 30 liver image volumes (cases)

Full size image

Analysis of automatic segmenation results

The following evaluation is based on automatic segmentation outcomes. The 3D reconstruction presented in Fig. 9 illustrates the depth of the pancreas curvature, with the “ragged” surface concentrated at the higher threshold with varying shades of yellow. Differences in pancreas volume may also be associated with age, gender and health status of the participants [37], and a wide variation in organ shape, size and curvature are established from this sample of 85 image volumes: the mean volume is 69.30 ± 32.50cm3 and the mean 3D global curvature is (35.23 ± 6.83) × 10− 3. Figure 7 provides a statistical representation of these computations. Interestingly, Fig. 10 illustrates relatively less variation in curvature compared to the curvature of the pancreas. Analysing the results from 30 image volumes achieves a mean volume of 1547.48 ± 204.19cm3 and a mean 3D global curvature of (19.87 ± 3.62) × 10− 3. A statistical representation of the data is provided in Fig. 8.

Fig. 9
figure 9

The depth of a pancreas curvature is highlighted. The varying shades of yellow indicate the highest degree of curvature dent (0.25 to 0.30) relative to global curvature of the organ

Full size image
Fig. 10
figure 10

The depth of a liver curvature is highlighted. The varying shades of yellow indicate the highest degree of curvature dent (0.20 to 0.25) relative to global curvature of the organ

Full size image

When comparing the automatic segmentation results to respective the ground-truth for evaluation purposes, Table 1 displays the MAE (mean absolute error), RMSE (root mean square error) and MAP (mean absolute percentage error) computations for volume and curvature in respect to corresponding ground-truth values of 85 pancreas and 30 liver image volumes.

Table 1 MAE (mean absolute error), RMSE (root mean square error) and MAP (mean absolute percentage error) computations for volume and curvature
Full size table

Analysis of morphological features

In order to demonstrate the reliability of the proposed framework, the relationship between ground-truth morphological features is analysed alongside the morphology of automatic segmentation results. Figures 11 and 12 illustrate the relationship between pancreas volume and curvature for both the ground-truth (red) and the automatic segmentation (blue), respectively. As highlighted in both plots, the pancreas curvature versus volume displays a weak, but significant, negative correlation. Similarly, Figs. 13 and 14 illustrate a slight negative correlation for liver curvature and volume. For each dataset, the Wilcoxon signed-rank is computed using the differences between paired ground-truth volume and curvature values, and the differences between paired segmentation volume and curvature values. All four tests reveal a correlation that is statistically significant from zero (p < 0.0001).

Fig. 11
figure 11

Ground-truth pancreas curvature in relation to volume a displays negative correlation

Full size image
Fig. 12
figure 12

Automatic segmented pancreas curvature in relation to volume displays a negative correlation

Full size image
Fig. 13
figure 13

Ground-truth liver curvature in relation to volume displays a slight negative correlation

Full size image
Fig. 14
figure 14

Automatic segmented liver curvature in relation to volume displays a slight negative correlation

Full size image

Analysing the results presented in this paper, if the morphology of the pancreas can be described by the global curvature of the organ shape, then increasing values of global curvature are to be anticipated with decreasing organ size and thus giving rise to shape irregularities. Previous literature [38] reports that, for some patients with type 2 diabetes, there was an observed change in pancreas volume and morphology, in which the organ displayed an involuted morphology with serrated borders. Given the increase in global incidences of type 2 diabetes [17], there is a significant need to better understand the relationship between morphology and function. Furthermore, since NAFLD and type 2 diabetes mellitus frequently coexist and share the abnormalities of excess adiposity and insulin resistance, it can be necessary to compute objective measurements to stratify condition progression and remission, especially where changes in organ morphology have been observed but not quantified. Consequently, there is an apparent demand for objective methodologies that produce reliable morphological characterisations of the pancreas, liver and other major organs.

The segmentation program ran via a workstation with i7-59-30k-CPU at 3.50 GHz, and the mean time for segmenting the pancreas and liver is 25 and 32 minutes, respectively. Also, the computation of morphological measurements of volume and curvature is approximately 3 minutes. Using a GeForce Titan X GPU could potentially result in a tenfold decrease in run-time.

Conclusions

Accurate 3D organ reconstruction, coupled with the analysis of size and structure, can support medical professionals perform clinical detection, diagnosis and planning of treatment. This paper describes a computing tool that automatically extracts and analyses organs based on morphological features using anonymised volunteer data from magnetic resonance imaging (MRI) volumes. In addition to offering radiologists a 3D organ reconstruction view and a “second opinion”, the proposed computational tool can easily integrate into forensic science and biomedical research. Although the case studies in this paper focus on MRI modality, the proposed framework is extendable to computer tomography (CT) and ultrasound (US). Due to differences in the imaging acquisition of MRI, CT and US, the automatic segmentation approach may differ, but since the framework is modular, it is also adaptable to other anatomical structures. Future work aims to utilise the proposed framework in a large cohort of abdominal image volumes, and, explore the potential correlations between corresponding anthropomorphic and environmental factors, and organ structure and shape.