Similarity models for atlas-based segmentation of whole-body MRI volumes

Linköping University | Department of Biomedical Engineering

Master’s thesis, 30 ECTS | Biomedical Engineering

2020 | LIU-IMT-TFK-A–20/580--SE

Similarity models for atlas-based

segmentation of whole-body MRI

volumes

Elin Axberg

Ida Klerstad

Similarity models for atlas-based segmentation of whole-body MRI volumes

Master of Science Thesis in Biomedical Engineering

Elin Axberg, Ida Klerstad

Supervisor : Magnus Borga, IMT Linköping University

Supervisor : Sven-Patrik Hallsjö, HMV Linköping University

Examiner : Anders Eklund, IMT Linköping University

External supervisor : Hannes Järrendahl, AMRA Medical AB

Linköpings universitet SE–581 83 Linköping +46 13 28 10 00 , www.liu.se

Abstract

In order to analyse body composition of MRI (Magnetic Resonance Imaging) volumes, atlas-based segmentation is often used to retrieve information from specific organs or anatomical regions. The method behind this technique is to use an already segmented image volume, an atlas, to segment a target image volume by registering the volumes to each other. During this registration a deformation field will be calculated, which is applied to a segmented part of the atlas, resulting in the same anatomical segmentation in the target. The drawback with this method is that the quality of the segmentation is highly dependent on the similarity between the target and the atlas, which means that many atlases are needed to obtain good segmentation results in large sets of MRI volumes. One potential solution to overcome this problem is to create the deformation field between a target and an atlas as a sequence of small deformations between more similar bodies.

In this master thesis a new method for atlas-based segmentation has been developed, with the anticipation of obtaining good segmentation results regardless of the level of similarity between the target and the atlas. In order to do so, 4000 MRI volumes were used to create a manifold of human bodies, which represented a large variety of different body types. These MRI volumes were compared to each other and the calculated similarities were saved in matrices called similarity models. Three different similarity measures were used to create the models which resulted in three different versions of the model. In order to test the hypothesis of achieving good segmentation results when the deformation field was constructed as a sequence of small deformations, the similarity models were used to find the shortest path (the path with the least dissimilarity) between a target and an atlas in the manifold.

In order to evaluate the constructed similarity models, three MRI volumes were chosen as atlases and 100 MRI volumes were randomly picked to be used as targets. The shortest paths between these volumes were used to create the deformation fields as a sequence of small de-formations. The created fields were then used to segment the anatomical regions ASAT (ab-dominal subcutaneous adipose tissue), LPT (left posterior thigh) and VAT (visceral adipose tissue). The segmentation performance was measured with Dice Index, where segmentations constructed at AMRA Medical AB were used as ground truth. In order to put the results in relation to another segmentation method, direct deformation fields between the targets and the atlases were also created and the segmentation results were compared to the ground truth with the Dice Index. Two different types of transformation methods, one non-parametric and one affine transformation, were used to create the deformation fields in this master thesis. The evaluation showed that good segmentation results can be achieved for the segmentation of VAT for one of the constructed similarity models. These results were obtained when a non-parametric registration method was used to create the deformation fields. In order to achieve similar results for an affine registration and to improve the segmentation of other anatomical regions, further investigations are needed.

Acknowledgments

We would like to direct a big thank you to our supervisors at Linköping University, Magnus Borga and Sven-Patrik Hallsjö. Magnus Borga for all your support and great engagement to drive this project forward and Sven-Patrik Hallsjö for all the external support and inputs, which made it possible to finalize this thesis work.

We would also like to thank our colleagues at AMRA Medical for making us feel welcome at your company and for the opportunity to write this master thesis. An extra thanks goes to our supervisor Hannes Järrendahl for helping us with all code specific questions and to our thesis coworker Simon Kantedal for showing a great interest in our project.

Linköping, May 2020 Elin Axberg & Ida Klerstad

Contents

Abstract iii

Acknowledgments iv

Contents v

List of Figures vii

List of Tables x Notations xi 1 Introduction 1 1.1 Motivation . . . 1 1.2 Aims . . . 2 1.2.1 Problem statements . . . 3 1.3 Limitations . . . 3 2 Theory 5 2.1 Magnetic Resonance Imaging . . . 5

2.1.1 Fundamental theory . . . 6

2.1.2 Water and fat separation . . . 6

2.2 Image registration . . . 8

2.2.1 Interpolation . . . 9

2.2.2 Parametric image registration . . . 9

2.2.3 Non-parametric image registration . . . 10

2.3 Segmentation . . . 11

2.3.1 Atlas-based segmentation . . . 11

2.3.2 The Morphon . . . 11

2.3.3 Segmentation method at AMRA . . . 12

2.4 Image similarity measures . . . 13

2.4.1 Euclidean distance . . . 13

2.4.2 Dice Index . . . 13

2.4.3 Normalized Cross Correlation . . . 14

2.5 Graph theory . . . 14 2.5.1 Fundamental theory . . . 14 2.5.2 Matrix representation . . . 16 2.5.3 k-Nearest-Neighbours . . . 17 2.5.4 Dijkstra’s algorithm . . . 17 2.6 Related work . . . 18 3 Method 21 3.1 Data set . . . 21

3.2 Implementation overview . . . 23

3.3 Base model construction . . . 24

3.3.1 Similarity model 1 . . . 26

3.3.2 Similarity model 2 . . . 27

3.3.3 Similarity model 3 . . . 28

3.3.4 Summary of similarity models . . . 28

3.4 Transformations . . . 29

3.4.1 Morphon transformation . . . 29

3.4.2 Affine transformation . . . 30

3.5 Accumulation of deformation fields . . . 33

3.5.1 Implementation control . . . 33

3.6 Evaluation . . . 34

4 Results 38 4.1 Morphon transformation . . . 38

4.1.1 Segmentation results for VAT . . . 39

4.1.2 Segmentation results for ASAT . . . 41

4.1.3 Segmentation results for LPT . . . 42

4.2 Affine transformation . . . 44

4.2.1 Segmentation results for VAT . . . 45

4.2.2 Segmentation results for ASAT . . . 45

4.2.3 Segmentation results for LPT . . . 46

5 Discussion 48 5.1 Results . . . 48 5.1.1 Similarity measures . . . 49 5.1.2 Transformation methods . . . 52 5.2 Future work . . . 55 6 Conclusion 57 Bibliography 58

List of Figures

1.1 A representation of segmentation in multiple steps compared to direct segmenta-tion between target and atlas. . . 2 2.1 Illustration of the three orthogonal planes that are used to visualize an MRI volume. 5 2.2 A transverse oscillating magnetic field is applied to the subject (RF-excitation) in

order to construct an image of a slice of tissue. B0is the applied external field and ωLis the frequency of the proton spin. . . 6 2.3 Illustration of the resulting images when the fat and water signals are in phase (to

the left) and out of phase (to the right). . . 7 2.5 The leftmost illustration shows forward warping, where the displacement

de-scribe how each pixel in the atlas image should move to fit the target image. The rightmost illustration shows backward warping, where the displacement describe where in the atlas image each pixel in the target image should be sampled. . . 9 2.6 Example of three different segmented regions. The red segment contain

subcuta-neous fat tissue, the blue segment contain visceral fat tissue and the green segment contain the muscle tissue of the left posterior thigh. Image source: AMRA, UK Biobank. . . 11 2.7 A representation of a graph G with five vertices and six edges. . . 15 2.8 The two upper graphs show the difference between an undirected graph (left) and

a directed graph (right). The two lower graphs show that an undirected graph (left) is equal to a graph where the two arcs are going in the opposite directions (right). . . 15 2.9 An example of an undirected weighted graph with three nodes and two edges.

The edges are labeled with a number which is called a weight. . . 16 2.10 An example of an undirected weighted graph with six nodes and eight edges. The

weights represents the distance between one node to another. . . 17 2.11 The figure illustrates the graph from Figure 2.10 after the shortest paths between

node x1and all other nodes have been found with Dijkstra’s algorithm. . . 18 3.1 Example of an MRI volume from the data set in the two planes; coronal (left) and

sagittal (right). . . 22 3.2 Example of an MRI volume from the data set in the axial plane. . . 22 3.3 Example of an MRI volume from the data set with the three masks ASAT (left),

LPT (middle) and VAT (right) are marked in red. . . 23 3.4 Flowchart over the steps performed to construct and use the similarity models and

evaluate the segmentation result. . . 23 3.5 Illustration of the variation of body types in the constructed similarity model of

3693 MRI volumes. The dots show the projections onto the two first eigenvectors of the normalized Laplacian matrix and the two branches that are seen represent the two sexes. . . 25

3.6 Illustration of the variation of body types in the constructed similarity model of 3693 MRI volumes. The dots show the projections onto the two first eigenvectors of the normalized Laplacian matrix and the figure demonstrates where specific MRI volumes are located in the graph. The two branches that are seen represent the two sexes. . . 26 3.7 Summary of the similarity measures used to construct the different similarity

models. . . 29 3.8 Selected points, marked in red, in a target image volume shown in the coronal plane. 31 3.9 Selected points, marked in red, in a target image volume shown in two different

axial planes. . . 31 3.10 Schematic illustration of how to find the secondary points for three primary

points, using a deformation field. . . 32 3.11 Schematic illustration of the accumulating process, where n defines the number of

steps between target and atlas, aithe continuously updated accumulated field and dievery single deformation field between two MRI volumes, where i=1, 2, ..., n. . 33 3.12 The final accumulated field constructed of a series of deformation fields between

nine cuboids with a ten degree displacement from each other. . . 34 3.13 The physical appearance of the three atlases selected for evaluation and their

ap-proximate position in the manifold. The dots show the projections onto the two first eigenvectors of the normalized Laplacian matrix. . . 35 4.1 Comparison of the mean Dice Index for the three masks segmented with the

dif-ferent segmentation methods, plotted together with the standard deviation. . . 39 4.2 The figure shows the difference in Dice Index between similarity model 2-VAT and

the direct method for all masks constructed from the paths in 3A, plotted together with the number of steps for each path. . . 40 4.3 The figure shows the difference in volume between similarity model 2-VAT and

ground truth for each segmentation in 1A, plotted together with the mean of each segmentation volume. . . 40 4.4 The figure shows the difference in Dice Index between similarity model 2-VAT and

the direct method for all masks constructed from the paths in 3A, plotted together with the number of steps for each path. . . 41 4.5 The figure shows the difference in volume between similarity model 2-VAT and

ground truth for each segmentation in 1A, plotted together with the mean of each segmentation volume. . . 42 4.6 The figure shows the difference in Dice Index between similarity model 2-LPT and

the direct method for all masks constructed from the paths in 3A, plotted together with the number of steps for each path. . . 43 4.7 The figure shows the difference in volume between similarity model 2-LPT and

ground truth for each segmentation in 1A, plotted together with the mean of each segmentation volume. . . 43 4.8 Comparison of the mean Dice Index for the three masks segmented with the

dif-ferent segmentation methods, plotted together with the standard deviation. . . 44 5.1 Segmented masks from similarity model 1 (upper row), similarity model 2-VAT

(middle row) and ground truth masks (bottom row). The images to the left show the ASAT mask, the images in the middle show the VAT mask and the images to the right show the LPT mask. . . 50 5.2 Illustration of a path between a target and an atlas in similarity model 1 (red) and

5.3 A zoomed-in version of Figure 5.2. As seen, the two paths sometimes change direction in their way from the target to the atlas and it does not look like they are taking the shortest path. The reason for this appearance is that the plot illustrates a path in a high dimensional space in a 2D plot. . . 51 5.4 The appearance of a deformed ASAT mask (upper left), a deformed VAT mask

(upper middle) and a deformed LPT mask (upper right) produced by performing a sequence of small deformations along similarity model 3 and their corresponding ground truth masks (below each deformed mask), plotted in the target body. . . . 54

List of Tables

4.1 DI statistics for Morphon transformations applied on the mask VAT. 3A represents the case with 300 paths and 1A represents the case with 100 paths. . . 39 4.2 DI statistics for Morphon transformations applied on the mask ASAT. 3A

repre-sents the case with 300 paths and 1A reprerepre-sents the case with 100 paths. . . 41 4.3 DI statistics for Morphon transformations applied on the mask LPT. 3A represents

the case with 300 paths and 1A represents the case with 100 paths. . . 42 4.4 DI statistics for affine transformations applied on the mask VAT. 3A represents the

case with 300 paths and 1A represents the case with 100 paths. . . 45 4.5 DI statistics for affine transformations applied on the mask ASAT. 3A represents

the case with 300 paths and 1A represents the case with 100 paths. . . 45 4.6 DI statistics for affine transformations applied on the mask LPT. 3A represents the

Notations

Abbreviation Meaning

ASAT Abdominal Subcutaneous Adipose Tissue

DI Dice Index

kNN k Nearest Neighbours

LPT Left Posterior Thigh

MRI Magnetic Resonance Imaging NCC Normalized Cross Correlation VAT Visceral Adipose Tissue

1

Introduction

Obesity is an increasing problem that is highly related to type-2 diabetes, cardio-vascular dis-eases and neurovascular disdis-eases, which all result in decreased quality of life [1]. To under-stand the correlation between obesity and metabolic risks, a more individual measure than BMI (Body Mass Index) is needed. For this, the analysis of body composition using Magnetic Resonance Imaging (MRI) has been accepted as gold standard [2]. In order to analyse body composition, the MRI volumes need to be segmented. It is very costly to manually segment MRI volumes and automatic segmentation methods have been developed in order to study large amounts of images [3]. Although, these methods are not without drawbacks and lim-itations. To optimise and further develop the autonomic segmentation methods, different approaches have been investigated [4].

This master thesis was performed at AMRA Medical AB [5] and examined at the Department of Biomedical Engineering, IMT, at Linköping University. The main focus of the thesis has been to construct and evaluate a new model for segmentation of whole-body MRI volumes. In this chapter, a motivation to why it was relevant to investigate a new segmentation model is presented, followed by the problem statements and the limitations of the thesis.

1.1

Motivation

Medical images need to be segmented in order to perform calculations on specific organs or anatomical regions. This can either be done manually, where a person with knowledge of anatomy and medical images will mark the wanted segment, or it can be done automatically. When the latter approach is used, a number of problematic areas will become significant such as image artifacts, noise, variability between different patients and low contrast between different regions. One way to overcome these problems is to use atlas-based segmentation [3]. The theory behind this method is that an already segmented image volume, an atlas, is used to segment the target image volume by rigidly or non-rigidly registering the target to fit the atlas. When performing this registration a deformation field will be calculated, which can be seen as a description of how each point in the target should be moved to fit the atlas. This deformation field is applied to a segmented part of the atlas which will result in a segmentation of the corresponding anatomic region in the target. Since this segmentation is based on the segmentation of the atlas, the quality of the segmentation depends on how well the atlas image volume can represent the target image volume. This means that the resulting segmentation will not be successful if the volumes are too different [6].

To further develop the existing segmentation methods, different approaches have been in-vestigated, such as GRAM [7], LEAP [8] and GIF [9]. Another method, described in [4], uses empirically learned multi-scaled anatomical manifolds. These manifolds are representations of different body types in a high-dimensional volume space, which can be seen as a graph where each body is an image volume represented as a node in the graph. The distances be-tween the nodes represent the simiarities bebe-tween the bodies. All nodes can be looked at as two different sets: labeled (atlases) and unlabeled (targets) volumes. The graph can be used

1.2. Aims to find the shortest path, the path with the least dissimilarity, from target to atlas. By using this path, the segmentation can be done in several small deformations instead of one large deformation, as shown in Figure 1.1 [4].

Figure 1.1: A representation of segmentation in multiple steps compared to direct segmentation

be-tween target and atlas [6].

The graph representing the manifold was, in the study described above, constructed by cal-culating the similarity between all image volumes and with these calculations constructing a similarity matrix, which can be seen as an undirected representation of the graph. The similarities were calculated after the image volumes had been registered to each other. The method was performed on 36 whole-body MRI volumes and the used similarity measure for comparing the image volumes was normalized cross correlation, NCC. The results from this article showed that the segmentations were improved when using the shortest path com-pared to a direct segmentation, especially in the worst-case performance. This indicates that when the target and the atlas are similar the proposed method and the direct method per-formed equally, but when the target and the atlas are very dissimilar the proposed method is better. The conclusion of this study implies that multiple small deformations between MRI volumes results in a better segmentation than a direct deformation between target and atlas [4].

1.2

Aims

A large limitation of the study described above was the small amount of data sets used [4]. The aim of this master thesis is to construct and evaluate a similarity model, earlier men-tioned as manifold/graph, of whole-body MRI volumes. In the context of this thesis, a simi-larity model is a description of the similarities in a set of MRI volumes, which is presented in a matrix. The construction of the similarity model is done similar to the mentioned study but

1.3. Limitations with a larger scale of data sets, using approximately 4000 whole-body MRI volumes. Consid-ering the large amount of data sets the used similarity measure, for construction of the model, needs to work without having to register the image volumes. This is due to the fact that it is too time consuming to register all image volumes to all other image volumes to construct the model. The idea is that the constructed model will represent a large variety of human bodies which potentially can be used in the future to guarantee a good atlas-based segmentation for all body types.

1.2.1

Problem statements

The master thesis aims to answer the following questions.

1. How can a similarity model be constructed for a large set of whole-body MRI volumes? 2. Can good segmentation results be achieved when performing a sequence of small

de-formations along the constructed model from target to atlas?

3. Can equally good segmentation result be achieved if an affine registration method is used instead of a non-parametric registration method, when performing a sequence of small deformations along the constructed model?

1.3

Limitations

To answer the stated questions, the master thesis is limited by the following: • The allocated time is limited to 20 weeks.

• The provided MRI volumes are collected from the UK-biobank, which only represents the variety of human bodies of British people between the ages 40 to 80 year. This is further explained in section 3.1.

2

Theory

In this chapter the theory behind the master thesis is presented. This includes an overview of magnetic resonance imaging, image registration, segmentation, similarity measures, graph theory and related work. The reader is assumed to be studying on a master degree level in engineering and fundamental concepts are therefore not explained.

2.1

Magnetic Resonance Imaging

Magnetic resonance imaging (MRI) is a noninvasive technique that provides information about the internal anatomy and physiology of living subjects [10]. The technique has no known side effects which makes it suitable for human studies and the MR images contain information about both the rough structural properties of the tissue and its biochemistry. This means that they are extremely versatile and can be used for many purposes, such as analysing tissue metabolites and analysing specific anatomic structures. Depending on the area of interest, MRI can be used to construct images in 2D or 3D and the visualization of the images are done in the three orthogonal planes; sagittal, coronal and transverse/axial [10]. The orientation of these three planes are illustrated in Figure 2.1.

Figure 2.1: Illustration of the three orthogonal planes that are used to visualize an MRI volume.

2.1. Magnetic Resonance Imaging

2.1.1

Fundamental theory

The MRI technique is based on the interaction of an external magnetic field with the spin of a nucleus. In the case of MRI, the dominant nucleus is the proton in the hydrogen atom which is present in water, fat and other organic molecules [12]. During an MRI examination an external magnetic field, B0, is applied to the subject which causes the spin of the hydrogen atoms in the tissue to precess parallel or anti-parallel about the direction of the applied field [10]. In order to detect a measurable signal, a transverse oscillating magnetic field is applied to a specific region of the subject. This will cause the hydrogen atoms in that region to tip away from the external field and create a detectable magnetization component in the trans-verse plane. Since all nuclei in the selected tissue are excited at the same time, the detected signal contains information from all points of the selected tissue. In order to reconstruct the signal into an image, each part of the signal needs to be correctly identified and mapped onto the correct pixel within the image [10].

B0 ωL> ωRF ωL= ωRF ωL< ωRF RF excitation M0 M0 M

Figure 2.2: A transverse oscillating magnetic field is applied to the subject (RF-excitation) in order to

construct an image of a slice of tissue. B0is the applied external field and ωLis the frequency of the proton spin.

2.1.2

Water and fat separation

The MR signal contains information from both water and fat tissue [13]. In order to estimate the amount of fat and muscle volume in different organs or anatomical regions, different methods can be used to separate the signal into one fat signal and one water signal [12]. One of these methods is called Dixon imaging, which provides one image showing the content of fat and one image showing the content of water [13]. This is done by utilising the chemical shift difference between water and fat [13], which is a small difference in the proton’s preces-sion frequency due to molecular differences in the tissues [10]. The chemical shift is utilised by the construction of one image when the signal from fat and water are in phase, equation 2.1, and one image when the signals are 180° out of phase, equation 2.2 [13]. w and f repre-sent the magnitudes of the magnetizations at all given pixels from the water signal and the fat signal, respectively. I1and I2are the two constructed images. After some calculations, one fat image and one water image can be obtained as the sum and difference of the in-phase and off-phase image. An example of the results of Dixon imaging can be seen in Figure 2.3 and 2.4.

2.1. Magnetic Resonance Imaging

I2=w ´ f (2.2)

Figure 2.3: Illustration of the resulting images when the fat and water signals are in phase (to the left)

and out of phase (to the right). Image source: AMRA, UK Biobank.

Figure 2.4: Illustration of the resulting fat and water images obtained from the images in Figure 2.3.

2.2. Image registration

2.2

Image registration

Medical images can be used in numerous fields of health care, such as diagnosis, treatment planning or disease monitoring [14]. This can for example be used to monitor the develop-ment or growth of a tumour, where images are taken at different times to enable comparison between them. To analyse the images, image registration plays an essential role. Because of its importance, many different image registration algorithms have been proposed, developed and used over the years [14].

One can describe image registration as a preprocessing step in medical imaging, where it is used to align two images, atlas and target, of the same anatomical structure [14]. The two images are bound to be misaligned since they have been taken at different times, with different modalities [14] or contain two different subjects [15]. The alignment is done by estimating a spatial transformation, T, between the images, such that when transforming the atlas image, IA˝T, it will be as similar as possible to the target image, IT [15]. The transformation can in turn be defined as T(x) = x+u(x), where x describes each pixel or voxel (3D equivalence to pixel) and u is the displacement. The main goal is to minimize the distance between all anatomically corresponding points in the atlas and the target [15]. There are three components that can describe all registration techniques; transformation, sim-ilarity measure and optimization [14]. The transformation is described above, but the simi-larity measure is used to measure the simisimi-larity between atlas and target. The optimization is in turn used to optimize the transformation estimation as a function of the similarity mea-sure [14]. There exist different methods/approaches for all these three components that can be combined in numerous ways. The image registration algorithms can be divided into two types of registration methods, parametric and non-parametric [15]. The parametric registra-tion algorithms include those transformaregistra-tions with reduced number of degrees of freedom, e.g. rigid or affine transformations. Non-parametric methods are not limited to a reduced space and estimates a displacement vector for each pixel or voxel [15]. These two methods are described further in the subsequent sections.

The transformation can be done in two different directions, either by forward warping or by backward warping [15]. Figure 2.5 illustrates the two approaches, where the displacement either go from atlas to target (forward warping) or from target to atlas (backward warping). Since it is the atlas image that will be deformed to match the target image, each pixel or voxel in the target needs to be defined with a label from the atlas image. Therefore, it is important to know, for each pixel or voxel in the target image, where to retrieve the information in the atlas image. When using the forward warping approach, this might not be ensured and holes can appear in the image. Because of this, the backward warping approach is normally used in image registration [15].

2.2. Image registration

IA

IT

IA

IT

Figure 2.5: The leftmost illustration shows forward warping, where the displacement describe how

each pixel in the atlas image should move to fit the target image. The rightmost illustration shows backward warping, where the displacement describe where in the atlas image each pixel in the target image should be sampled.

2.2.1

Interpolation

Something that can be noted in Figure 2.5 is that, regardless of the warping method, the displacement might not point to an exact pixel centre. To find the value in the correspond-ing pixel, some kind of interpolation is needed. There exist different interpolation methods, where nearest neighbour and linear interpolation are two of them.

Nearest neighbour

The nearest neighbour interpolation method assigns the intensity value of the spatially closest neighbour to the target point [16]. It is seldom used in image registration because of its poorly suitable properties in an optimization framework, in combination of also generating visually unpleasant results [15]. Regardless of this, the nearest neighbour interpolation is preferably used when transforming binary masks [15].

Linear interpolation

The linear interpolation method assigns the combined intensity values of the neighbouring points to the target point [16]. The linear distance between the different neighbouring points and the target point dedicate how much each neighbour adds to the final intensity [16]. Lin-ear interpolation is a well used interpolation method in image registration that benefits low computational cost and ensures that no spurious oscillations occur [15].

2.2.2

Parametric image registration

For parametric image registration methods, a parameterization of the transformation has been done [15]. There are different kinds of parametric registration methods, which are de-fined according to the degrees of freedom [16]. For example a rigid transformation only includes translation and rotation, where an affine transformation also include scaling and shearing [16]. This means that a three-dimensional rigid transformation would contain six independent parameters, three translations along the x-, y- and z-axis and three rotations around the x-, y- and z-axis [17]. Since the affine transformation method is of essence in this thesis, it is described further below.

2.2. Image registration

Affine transform

A three-dimensional affine transformation contains twelve independent parameters, three translations, three rotations, three scalings and three shearings [17]. One way to compute such a transformation can be done by selecting corresponding landmarks, or points, in the atlas and target image [15]. The landmarks can be marked both manually or automatically, either way there is no straight-forward approach how to select these landmarks. The idea is to compute the transformation by finding the minimum distance between the landmarks in the atlas image with the corresponding landmarks in the target image. This can be done by solving a least squares method, assuming that the number of landmarks in the atlas image are larger than the amount of independent parameters and that the landmarks are not selected on a straight line [15]. The major drawback with this landmark-based approach is that the transformation only relies on the selected landmarks [15].

A three-dimensional affine transformation can be represented by a three-by-four matrix [17]. This matrix applied to any coordinate (x, y, z) describes the new coordinate (x’, y’, z’), as is shown in equation 2.3. As can be seen, the original coordinates are assigned to a vector with the final element equal to 1. This element is needed to make the mathematics work properly [17], i.e. it ensures that the translation components (e14, e24, e34) are added to the transformation.   x1 y1 z1  =   e11 e12 e13 e14 e21 e22 e23 e24 e31 e32 e33 e34       x y z 1     (2.3)

2.2.3

Non-parametric image registration

As it appears, most of the human body does not conform to a rigid or affine approximation, which has lead to much research in the construction of non-parametric image registration methods [14]. Non-parametric transformation is a challenging task, partly because of the high degrees of freedom [14]. Numerous algorithms have been constructed and they all require long computation times, which is one of the major drawbacks with this registration method [14]. On the other hand, a non-parametric transformation normally returns a much better representation of the deformed image than parametric methods [18].

Instead of parameterizing the transformation, a deformation field is constructed in non-parametric image registration. This deformation field contains a displacement vector for each pixel in the image, describing how to transform the pixel to imitate the corresponding atlas image pixel [19]. An example of a non-parametric registration method is described further in section 2.3.2, The Morphon.

For non-parametric image registration, regularization is needed, which can be seen as a low pass filtering on the registration process [15]. This is because the deformation field estimation neither models spatial dependencies among the deformation vectors nor enforces smoothness [20]. The regularizer is incorporated into this registration method to ensure smoothness of the deformation field [16]. The regularizer is also used to favor any specific properties in the solution that the user requires [19].

2.3. Segmentation

2.3

Segmentation

In order to enable computerized analysis of a certain part of an image, the part needs to be segmented [21]. A segmentation is an extraction of a region in an image and in medical imaging this region can for example be organs, body tissue or even tumours [22]. Figure 2.6 shows an example of three different segmented regions marked in red, blue and green. Segmentation plays an essential role in medical imaging to enable computer-aided diagnosis systems. However, construction of a segmentation is challenging and the outcome might be affected by various problems such as noise, low contrast, illumination, and irregularity of the object boundaries [22].

Figure 2.6: Example of three different segmented regions. The red segment contain subcutaneous fat

tissue, the blue segment contain visceral fat tissue and the green segment contain the muscle tissue of the left posterior thigh. Image source: AMRA, UK Biobank.

2.3.1

Atlas-based segmentation

An atlas is typically an average representation of a group of patients [15]. Atlas-based seg-mentation is a segseg-mentation method where a labeled, i.e. segmented, atlas image is used to segment an unlabeled target image [15]. This is done by image registration, where the con-structed deformation field is applied to a segmented part in the atlas image to describe the corresponding part in the target image [15]. This segmentation method is normally used in medical imaging to overcome the various problems associated with segmentation [3].

2.3.2

The Morphon

The Morphon is a non-parametric atlas-based image registration and/or segmentation algo-rithm [15]. The algoalgo-rithm operates iteratively and updates the deformation field between target and atlas as the images go from coarse to fine scale [23]. The deformation is more rigid for the coarse scales and less rigid in the end of the process, at the fine scale. For each scale the deformation field is successively updated, or accumulated, by an incremental defor-mation field at the current scale. This continuously updated defordefor-mation field is called the accumulated deformation field. Each scale is iterated several times to ensure better precision of the incremental field, and for each iteration the certainty of the deformation fields are

con-2.3. Segmentation trolled [23]. The certainty is represented by a field that describes the certainty of a rightful deformation for each point in the target image.

To compute the incremental field, the Morphon uses a phase-based quadrature filter approach [23]. For further understanding of this approach the reader is referred to the original article of the Morphon presented by Knutsson and Andersson [23]. During the field computation the estimation of a corresponding certainty field is also done [15]. This incremental certainty field plays a significant role in the regularization and accumulation steps [15], which will be explained in more detail in the sections below.

Regularization

To provide a robust segmentation algorithm, the incremental deformation field is regularized in every iteration [23]. This is done by filtering with a Gaussian lowpass kernel and normal-ized convolution. The function for this can be seen in equation 2.4, where diis the incremental deformation field, ci is the incremental certainty field and g is the Gaussian kernel. As the scales go finer the variance of g decreases, resulting in a less rigid deformation [23].

di=

(cidi)˚g

ci˚g (2.4)

Field accumulation

Both the incremental deformation field and certainty field are updated for each iteration, resulting in a continuously accumulated deformation field and certainty field until the last iteration for the final scale [23]. Equation 2.5 and 2.6 show how the accumulated deformation field and certainty field are updated, respectively, wheres di is the incremental deformation field, ciis the incremental certainty field, ca is the accumulated certainty field and da is the accumulated deformation field [15]. Regarding the accumulated deformation field, it is clear that the dihas very low impact on da, unless the cihas similar magnitude as ca[23].

da= cada +ci(da+di) ca+ci (2.5) ca= c 2 a+ci(ca+ci) ca+ci (2.6)

When all iterations have been run, the final accumulated deformation field describes how the target volume can be deformed to represent the atlas volume. This means that when the accumulated deformation field is applied on a certain segmented area in the atlas, e.g. the liver, a muscle or fat region, it transforms that anatomical region from atlas to target [3].

2.3.3

Segmentation method at AMRA

AMRA uses atlas-based segmentation to segment whole-body MRI volumes. Instead of using one atlas to register the target a few atlases, similar to the body type of the target, are chosen

2.4. Image similarity measures on all selected atlases and the final segmentation is then created by combining the different segmentations by using a voting scheme. This means that if a certain fraction of the different segmentations include a specific region in the MRI volume this region will be a part of the final segmentation [3]. To choose atlases similar to the target, a signature vector of the target is constructed. This is done by firstly doing a centring of the fat image volume and then computing the projections of the centred image volume in three orthogonal planes (sagittal, coronal and axial). The pixel values of the three plane projections are then reformatted into a single vector. By dimensionality reduction of the constructed vector, the resulting signature vector consists of 20 floating-point numbers. This signature vector is then compared with the signature vectors in a database of many atlases and the seven most similar ones are chosen to segment the target. In the last step of this process, the segmentation is observed manually, where adjustment can be done if needed.

2.4

Image similarity measures

Image similarity usually refers to the fact that corresponding spatial positions in two im-ages have similar intensity values [15]. When the comparing imim-ages are medical imim-ages the concept commonly means that corresponding spatial positions refer to the same anatomical locations [15]. The concept of image similarity plays an important role in a number of image processing algorithms such as image restoration, image registration and image segmentation [24]. Different indices have been developed to measure image similarity, where some of them determines the similarity by comparing the corresponding pixel intensities, which relies on the assumption that the images are perfectly registered and are at the same scale. This will become problematic when there for instance exist small translations or rotations between the images, since this will result in a low similarity even if the objects in the images are very similar. Other similarity indices firstly find pixel correspondence based on intensity and then calculates the similarity by comparing the geometric transformations between corresponding pixels [24]. Some commonly used similarity measures will be described in the subsections be-low.

2.4.1

Euclidean distance

Consider two M by N images x and y with pixel intensity value vectors x= (x1, x2, ...., xMN) and y = (y1, y2, ...., yMN). The Euclidean distance between these two images is defined in equation 2.7, where a high Euclidean distance represents a low similarity and a low Euclidean distance represents a high similarity [25]. In image processing the Euclidean distance is often called mean squared error (MSE) or sum of squared differences (SSD).

dE(x, y) = g f f e MN ÿ k=1 (xk´yk)2 (2.7)

2.4.2

Dice Index

The Dice Index (DI) is a very simple calculation to determine similarity of two binary images by measuring the degree of spatial overlap between the images [24]. A Dice Index value of 1

2.5. Graph theory signifies a total spatial overlap and a Dice Index value of 0 signifies no overlap. The equation for the Dice Index is shown below, where X and Y represents the two binary images [24].

DI= 2|XŞ Y|

|X|+|Y| (2.8)

2.4.3

Normalized Cross Correlation

Another image similarity measure is Normalized Cross Correlation (NCC) [4]. The equation for this measure is shown in equation 2.9, where X and Y are vector representations of the images, X ¨ Y is the scalar product of the vectors and } ¨ } is the euclidean norm [4]. A NCC value close to 1 represents high similarity and a value close to 0 represents low similarity.

NCC= (X ¨ Y)

}X}}Y} (2.9)

2.5

Graph theory

A graph is, very simply put, a set of points or elements where pairs of the elements can be connected [26]. In graph theory, the points are called vertices and the connections are called edges. Graphs are very useful tools in numerous applications and fields of sciences to efficiently solve complex problems. A graph might, for example, represent a map that can describe the fastest way between one city to another, wheres the vertices are cities and the edges are roads. It can also represent a social network, where the vertices are people and the edges describe who is friends with whom [26].

Graph theory has its own terminology and denotations. There exist various types of graphs, different representations of them and related algorithms to extract useful information from them. The sections bellow only addresses those parts in the theory that are useful and rele-vant for this thesis.

2.5.1

Fundamental theory

As described above, the points in a graph are called vertices and the lines are called edges. If two vertices are connected by an edge, they are said to be adjacent [26]. One vertex can have numerous vertices adjacent [26]. In Figure 2.7 an example of a graph is shown, where the vertices are denoted x and the edges are denoted e.

2.5. Graph theory

Figure 2.7: A representation of a graph G with five vertices and six edges.

There are numerous ways to store graphs in computer memory [26]. One way is to use a so called adjacency matrix. This is a matrix where each row and column represent a vertex, resulting in a square matrix with the same size as the square of the number of vertices in the graph. If two vertices, e.g. xi and xj, are adjacent it would be represented with a 1 in the matrix at row i and column j and at row j and column i. If they however are not adjacent it would be represented with a 0 [26]. The adjacency matrix of the graph in Figure 2.7 would thus be as follows: A(G) =       0 1 1 0 0 1 0 1 0 0 1 1 0 1 1 0 0 1 0 1 0 0 1 1 0       (2.10)

As can be noted, this particular matrix is symmetrical as it is equal to its transpose. An adjacency matrix is always symmetrical if the graph is undirected. In an undirected graph the way from xi to xjis the same as going from xjto xi[26]. However, if the direction does matter, the edge only exists in one direction. This kind of graph is called a directed graph. A pair of vertices that are connected by an edge in one direction, are called an arc. If two arcs are going in the opposite direction of each other it would result in an equal edge as in an undirected graph. Figure 2.8 summarises this and shows the relation between arcs and edges, undirected and directed graphs.

Figure 2.8: The two upper graphs show the difference between an undirected graph (left) and a directed

graph (right). The two lower graphs show that an edge in an undirected graph (left) is equal to a graph where the two arcs are going in the opposite directions (right).

2.5. Graph theory The graph shown in Figure 2.7 is a connected graph. This is because all vertices are somewhat connected by edges to each other. However, if we were to remove the edges e2and e3, the graph would be a disconnected graph. This is because vertices x1and x2still belong to the graph G but are not connected to the rest of the graph by any edges [26].

2.5.2

Matrix representation

The graph presented in Figure 2.7 is an example of an undirected simple graph. This means that all edges are pointed in both directions and are represented with the value 1 in the ad-jacency matrix. Another example of a graph can be seen in Figure 2.9, where the edges are labeled with a number. This number indicates the weight of the edge, which for example can represent the similarity between two connected nodes [27]. A high weight would in this case represent high similarity and a low weight would represent low similarity. One matrix rep-resentation for this type of graph is the weighted adjacency matrix W [27], which for the graph in Figure 2.9 is illustrated in matrix 2.11.

Figure 2.9: An example of an undirected weighted graph with three nodes and two edges. The edges

are labeled with a number which is called a weight.

Wij=   0 0.2 0 0.2 0 0.8 0 0.8 0   (2.11)

Another matrix representation of a weighted graph is the weighted degree matrix D [27]. As can be seen in matrix 2.12 this matrix is a diagonal matrix where the diagonal entries Diiare the sum of the weights (of the edges) connected to the node xi.

Dij=   0.2 0 0 0 1 0 0 0 0.8   (2.12)

When the weighted adjacency matrix and the weighted degree matrix are known, one can construct a Laplacian matrix L by performing the calculation L=D ´ W, or a normalized Lapla-cian matrix LNby preforming the calculation LN=D´1/2(D ´ W)D´1/2[27]. The Laplacian matrix for the graph in Figure 2.9 is presented in matrix 2.13. These types of matrices contain a lot of information about the properties of the corresponding graph and can for example be used to investigate the connectivity in the graph or to visualise high dimensional data in a low dimensional space [27]. The connectivity can be investigated by analysis of the

eigenval-2.5. Graph theory of isolated subgraphs within the graph. This means that if the graph only has one eigenvalue equal to zero the graph is completely connected [27]. In order to visualise high dimensional data, the Laplacian matrix is used to find low-dimensional representations where nodes with high similarity are close to each other in the embedding [27].

Lij=   0.2 ´0.2 0 ´0.2 1 ´0.8 0 ´0.8 0.8   (2.13)

2.5.3

k-Nearest-Neighbours

In some applications, the possibility to control the amount of connected nodes in a graph is desirable [28]. This can be done by construction of a k-Nearest-Neighbours (kNN) graph, which is a graph in which each node is connected to its k most similar neighbours. By using an algorithm that for each row in a similarity matrix keeps the k edges with the highest weights and sets all other edges to zero, a kNN graph can be constructed [28].

2.5.4

Dijkstra’s algorithm

The graph shown in Figure 2.10 consists of six nodes and eight edges. Consider that this graph is a model of a map where the nodes represent cities and the edges represent roads. The weights represent the distance from one city to another. As can be seen in the figure there are several ways to go from one city to another and finding the shortest path between two cities can become problematic when dealing with larger matrices and graphs. This problem can be solved by using Dijkstra’s algorithm [29] which, when given a weighted adjacency matrix and a starting node, finds the shortest path from the starting node to all other nodes in the graph. When Dijkstra’s algorithm is applied to the graph in Figure 2.10 with starting node x1it would result in the paths shown in Figure 2.11.

Figure 2.10: An example of an undirected weighted graph with six nodes and eight edges. The weights

2.6. Related work

Figure 2.11: The figure illustrates the graph from Figure 2.10 after the shortest paths between node

x1and all other nodes have been found with Dijkstra’s algorithm.

2.6

Related work

The use of anatomical manifolds for atlas-based segmentation have been investigated in nu-merous studies. In this section, a few of these studies will be described in more detail to get a better understanding of the connection between them and this master thesis.

GIF [9], LEAP [8], GRAM [7], MUPPS [30] and a Markov Random Field Approach [31] are all examples of studies that relates to the subject of this master thesis. However, all of them are used for segmentation of the brain and have not been applied to whole-body imaging. GRAM [7], Geodesic Registration on Anatomical Manifolds, is an experiment performed to investigate the segmentation performance when using the shortest path (the path with the least image dissimilarity) in an anatomical manifold. The shortest path was used to construct the deformation field between a target and an atlas as a sequence of small deformations. The approach was tested on both simulated and real images of the human and mouse brain and the segmentation performance was calculated using Mean Square Error (MSE). The manifold in this experiment was constructed by registering every image to every other image on a course scale and then measuring the similarity as the distance between two images after registration. When the approach was used to segment parts of the human brain, 416 MRI volumes were used to construct the manifold. The results showed that the approach of using the shortest path lead to decreased MSEs with smoother deformation fields in comparison to a direct deformation between target and atlas.

Another study that utilizes anatomical manifolds for atlas-based segmentation is LEAP [8], Learning Embeddings for Atlas Propagation. This study differs from GRAM since the mani-fold in this study is used to construct new atlases in a multi-atlas segmentation process, which means that a few atlases are chosen to segment a target. In order to construct the manifold, all images are compared with the intensity-based similarity measure Normalized Mutual Infor-mation (NMI) over a region of interest (ROI) around the hippocampus in the brain. This ROI is defined by multi-atlas segmentation of the hippocampus prior to the process of comparing all images. The study uses 30 MR images as initial atlases and 796 MR images as initial tar-gets to construct the manifold. The initial atlases are then used to segment the tartar-gets that are closest to them in the manifold in an iterative process, meaning that when a target has been segmented it will act as an atlas in the next iteration. This is repeated until all targets in the manifold are segmented. The results from this study showed that the segmentation perfor-mance increased with the proposed method compared to standard multi-atlas segmentation since the proposed method avoids the need to estimate large deformations. This means that this approach could be a good way of constructing new atlases in an efficient way.

The study described in [4] was briefly explained in section 1.1 which showed a strong corre-lation to this master thesis. In this study, 36 whole-body MRI volumes were used to construct

2.6. Related work a manifold of human bodies. The similarity between two image volumes was calculated with the similarity measure NCC after one of the image volumes had non-parametrically been registered to the other image volume. By repeating this calculation for all 36 image volumes, a similarity matrix could be constructed. This similarity matrix was then used to find the shortest path between target and atlas. The segmentation performance, when the deformation field was created as a sequence of small deformations along the shortest path, was evaluated with the Dice Index on segmented livers by using manually segmented livers as ground truth. This was compared to the segmentation performance when the deformation field was created directly between target and atlas. The results from this study showed that the segmentation performance in the worst-case performance increased from 0.24 to 0.67, the standard deviation decreased from 0.10 to 0.04 and the mean increased from 0.82 to 0.87 with the proposed method. This strongly indicates that it is beneficial to create the deformation field between target and atlas as a sequence of small deformations along the shortest path. The largest difference between the studies described above and this master thesis is the amount of used image volumes and that none of the studies, apart from [4], has been ap-plied to whole-body imaging. This master thesis aims to fill this gap by using around 4000 whole-body MRI volumes to construct the manifold. By taking the larger data set into con-sideration, the approach for construction of the manifold needs to be adapted since it will be too time consuming to register all image volumes to all other image volumes. Instead, other methods for comparing images will be used. The hypothesis is that the larger set of MRI volumes will lead to a more robust model that will represent a larger variety of human bodies.

3

Method

In this chapter the chosen method for the master thesis project is explained. This includes an overview of the different parts of the method, followed by a more detailed explanation of the two parts Base model construction in section 3.3 and Transformations in section 3.4. In the last part of the chapter the chosen method for evaluation is explained.

3.1

Data set

The MRI volumes that were used in this master thesis were acquired from the UK Biobank [32], which is a very large, population-based prospective study with over 500 000 participants in the ages of 40 to 80 years. The study aims to be used for detailed investigations of diseases that affect people of middle and old age by continuous collection of data from physical mea-surements, medical imaging and questionnaires [32]. The used MRI volumes in this thesis are neck to knee 3D images of fat and water, with the approximated size of 224×174×370 vox-els and a resolution of 2×2×3 mm per voxel. All volumes have been segmented and binary masks for each segmented volume are also provided by AMRA. There were approximately as many men as women in the data set and the sexes are therefore considered to be equally represented. An example of an MRI volume from the data set is visualized in Figure 3.1 and 3.2. Another example is shown in Figure 3.3 where the segmented anatomical regions ASAT (abdominal subcutaneous adipose tissue), LPT (left posterior thigh) and VAT (visceral adipose tissue) are marked in red.

3.1. Data set

Figure 3.1: Example of an MRI volume from the data set in the two planes; coronal (left) and sagittal

(right).

3.2. Implementation overview

Figure 3.3: Example of an MRI volume from the data set with the three masks ASAT (left), LPT

(middle) and VAT (right) are marked in red.

3.2

Implementation overview

Figure 3.4 shows a fundamental overview of the chosen method for this master thesis. This in-cludes construction of similarity models (in the Introduction referred to as manifold or graph) and transformations which were used to create the deformation fields between a target and an atlas as a sequence of small deformations. The created fields were used to segment spe-cific parts of the MRI volumes which were compared with ground truth segmentations. The ground truth segmenations are the binary masks provided by AMRA. Each of the headings in the figure below are described in more details in the subsequent sections of this chapter, except for Results, which is presented in chapter 4. Furthermore, there were three different similarity models constructed to prosecute this master thesis, each with a different similarity measure. Also, two different image registration methods were used, one non-parametric and one parametric. The implementations were done in the programming language Python.

Base model Similarity models Transformations

Field accumulation Comparison with

ground truth Results

Figure 3.4: Flowchart over the steps performed to construct and use the similarity models and evaluate

3.3. Base model construction

3.3

Base model construction

The first step of the implementation was to construct the similarity model that would be used as a basis for similarity model 1, 2 and 3. Since the set of MRI volumes was going to be large, a simple similarity measure was needed to compare the image volumes. Therefore, the method used at AMRA today was chosen, where the Euclidean distance between two signature vector representations of the image volumes is calculated. After investigation of the database at AMRA, it turned out that the signature vectors for 3693 of the 3D images in the UK Biobank already had been calculated and these images were therefore chosen to construct the model. The similarity between two image volumes was calculated with equation 3.1 where w is the similarity and d is the Euclidean distance between the two signature vectors. The equation is a Gaussian kernel which gives a similarity measure that decreases when the Euclidean distance increases and ranges between zero and one. The value 1000 was chosen because it resulted in similarity values in a range that was easy to use for the master thesis students, without receiving too low values or a lot of values very close to zero. The process of calculating the similarity was repeated until all image volumes had been compared to each other and the calculated similarities were saved in a 3693x3693 similarity matrix. This similarity matrix is a weighted adjacency matrix where the weights represents the similarity between two image volumes.

w

=

e

In order to visualize the graph representation of the constructed similarity matrix in a 2D plot, a normalized Laplacian matrix was calculated. The rows of the constructed similarity matrix (which represents the similarity between one image volume and all other image volumes) were then projected onto the first two eigenvectors of the Laplacian matrix which resulted in one x-value and one y-value for each image volume. The resulting plot is visualized in Figure 3.5 where each dot represents an MRI volume. This plot was used to investigate if the constructed similarity model represented the MRI volumes in a reliable way by finding where specific volumes were located in the graph. The result from this investigation is illustrated in Figure 3.6 where it is clear that MRI volumes that are very dissimilar are located far away from each other in the graph. This indicates that the constructed model represents similarities in a reliable way. The representation of the sexes was also investigated which showed that the two branches that are seen in the plots represents the two sexes.

3.3. Base model construction

Figure 3.5: Illustration of the variation of body types in the constructed similarity model of 3693 MRI

volumes. The dots show the projections onto the two first eigenvectors of the normalized Laplacian matrix and the two branches that are seen represent the two sexes.

3.3. Base model construction

Figure 3.6: Illustration of the variation of body types in the constructed similarity model of 3693 MRI

volumes. The dots show the projections onto the two first eigenvectors of the normalized Laplacian matrix and the figure demonstrates where specific MRI volumes are located in the graph. The two branches that are seen represent the two sexes.

As mentioned earlier, the constructed similarity model is a weighted adjacency matrix which is completely connected since all edges are represented by a number larger than zero. Since the idea is to use the similarity model to construct the deformation field between target and atlas as a sequence of small deformations, the MRI volumes that are dissimilar do not need to be connected. In order to reduce the amount of connected nodes, a kNN matrix was com-puted which means that the k highest values in each row of the similarity matrix were kept and the rest of the values were set to zero. The value of k and the used similarity measures differ between similarity model 1, 2 and 3 and will be explained in more detail in the sub-sections below. In order to find the shortest path, which in this context refers to the path with the least total dissimilarity between target and atlas, Dijkstra’s algorithm was used on all similarity models.

3.3.1

Similarity model 1

The first version of the similarity model (similarity model 1) that was constructed was a kNN matrix, created from the base similarity model, with k= 3. The decision to use a low value for k was based on the fact that a low value would result in fewer connections in the graph and therefore fewer transformations would need to be constructed later. This was important since that process was very time consuming. The value k=3 was chosen because it was the lowest value that could be used without resulting in a disconnected graph. This was tested

3.3. Base model construction by investigation of the eigenvalues of the normalized Laplacian matrix. The result after kNN was approximately 16 000 connections in the graph. This is a larger number than 3693 ¨ 3 since one connection leads to a path in both directions and some nodes will therefore have more than three connections because they are forced to be connected to nodes that are connected to them. In order to evaluate the model and perform a sequence of small deformations from target to atlas, the transformations for the 16 000 connections were created as described in section 3.4.

3.3.2

Similarity model 2

The idea with this similarity model was to optimize the performance of similarity model 1. In order to accomplish this, a more accurate similarity measure than the Euclidean distance was needed and a decision was made to construct a model based on the Dice Index for specific masks. This could be done since the Morphon transformations had been created and saved for all connections in similarity model 1. The approach was to replace the values in similar-ity model 1 with the Dice Index between segmented masks (that were formed by using the created transformations) and the ground truth masks that were provided by AMRA. To get a better understanding of how this was done a simple example will be used. Consider the matrix S1below to represent similarity model 1 and S2to represent similarity model 2. In this example the similarity models only consist of five MRI volumes which results in a 5x5 matrix. The values in S1ans S2are only used to describe the used method and are not real calculated similarity measures.

S1=       0 0.7 0.8 0.6 0 0.7 0 0.6 0.7 0.2 0.8 0.6 0 0.3 0.4 0.6 0.7 0.3 0 0.8 0 0.2 0.4 0.8 0       (3.2)

To construct similarity model 2 in this example, the values in S1 will be replaced row by row. In the first row, one can see that image volume one is connected to image volume two, three and four. In order to find the new similarity value, between for example image volume one and two, the saved deformation field between these two bodies will be used to create a segmented mask of a region in volume one by applying the field to the ground truth mask in volume two. The Dice Index will then be calculated between the created mask and the ground truth mask for volume one. The same approach will be used to replace the rest of the values in row one and the process will continue in this way, row by row. As can be seen in the matrix S2this will lead to the same connections as in S1but with different values.

S2=       0 0.6 0.9 0.7 0 0.6 0 0.5 0.5 0.4 0.9 0.5 0 0.4 0.2 0.7 0.5 0.4 0 0.9 0 0.4 0.2 0.9 0       (3.3)

The different values in similarity model 2 compared to similarity model 1 will result in differ-ent paths between target and atlas when Dijkstra’s algorithm is applied. The masks that were

3.3. Base model construction used to construct this model were ASAT, LPT and VAT which means that there are three ver-sions of this model, one for each mask. These masks were only used to construct the models and all three versions of similarity model 2 can be used to segment other anatomical regions.

3.3.3

Similarity model 3

The third version of the similarity model (similarity model 3) that was constructed was a kNN matrix, created from the base similarity model, with k = 10. The reason for choosing a larger value for k in this model was because there was no guarantee that the similarity values in the base model were totally reliable and with a larger k the probability of losing good connections in the model decreased. After the kNN matrix had been constructed, the similarity values were replaced in a similar way as when constructing similarity model 2, but instead of using Dice Index, NCC was used. However, the volumes were not registered to each other prior to the similarity calculation as done in similarity model 2. Before calculating NCC between two fat-channel MRI volumes, the centre of mass for one of the volumes was changed to match the other. This was done to elude the problem of receiving a low NCC value if one of the bodies for example laid obliquely at the time of the MRI examination. When all values had been replaced, a new kNN matrix was constructed with k = 3. This was done in order to reduce unnecessary connections and therefore reducing the amount of transformations needed later. In the same way as for similarity model 2 the new values in this similarity model, and in some cases other neighbours, will result in different paths being chosen between target and atlas.

3.3.4

Summary of similarity models

Figure 3.7 shows a summary of the constructed similarity models and the used similarity measures. All models are kNN matrices of the base similarity model and the first similarity measure used for all models is therefore the Euclidean distance between signature vector representations of the image volumes. The additional similarity measures are shown in the figure.

3.4. Transformations Base model kNN (k=3) Similarity model 1 kNN (k=10) NCC kNN (k=3) Similarity model 3 kNN (k=3) Image registration Dice Index Similarity model 2 VAT Similarity model 2  LPT Similarity model 2 ASAT

Figure 3.7: Summary of the similarity measures used to construct the different similarity models.

3.4

Transformations

To be able to perform a sequence of small deformations along the constructed models, the information about how an MRI volume should be changed to match its neighbours must be known. In this thesis two kinds of transformation methods have been used. These two methods are further described below.

3.4.1

Morphon transformation

One of the used transformation methods was the non-parametric registration of dense de-formation fields, constructed by the Morphon (see section 2.3.2), where a dense dede-formation field contains a displacement vector for each voxel in the target image. These dense fields were constructed between all neighbours in similarity model 1, which resulted in approxi-mately 16 000 deformation fields. Since similarity model 1 and 2 have the same neighbours, the fields constructed for similarity model 1 could also be used for similarity model 2. As mentioned in section 2.2.3 the major drawback with a non-parametric transformation is the demanding computation time. With the Morphon the construction of one deformation field took approximately four minutes. This, together with the fact that 16 000 deformation fields were to be constructed, resulted in the computation time of approximately 44 days. Compared to the time allocated for this thesis work, that computation time was too long. This resulted in the decision of removing the highest resolution scale when constructing the deformations fields, which led to shorter computation time with the drawback of lower