ANALYSIS FOR JUVENILE IDIOPATHIC ARTHRITIS
MARTIN T.J. VAN KUIK
A THESIS PRESENTED FOR THE DEGREE OF MASTER OF SCIENCE
MEDICAL INFORMATICS
DEPARTMENT OF RADIOLOGY AND NUCLEAR MEDICINE UNIVERSITY OF AMSTERDAM
This thesis “Improvement of enhancement point detection in Dynamic Contrast-Enhanced MRI analysis for juvenile idiopathic arthritis” is the result of my scientific research project for graduation of the master Medical Informatics at the department of Radiology and Nuclear Medicine at Amsterdam UMC - location AMC.
I would like to thank the following people, as without their assistance, guidance and partici-pation this work would not have been possible.
First, I want to thank Cristina Lavini for her role as mentor and daily supervisor providing guidance and advice during this project. Also, I want to thank Josien van den Noort for her assistance and participation in this project. I want to thank Frans Voorbraak for his role as tutor and the assistance he provided during my project. Next, I would like to thank Robert Hemke, the radiologist who participated in the evaluation of the software and who made this research possible.
Also, I would like to thank all my colleagues at the department, including Anne-Sophie, Anouk, Ayla, Cato, Dana, Joost, Kyra, Leonard, Manon, Marijn, Miou, Rik, Ruud, and Sofieke for their presence during my project. I would like to thank my classmate Alma and my friend Giovanni for all the cups of coffee we drank together.
And last but certainly not least, I would like to thank my family and friends for their continuous support.
Introduction - Juvenile idiopathic arthritis (JIA) is a type of rheumatic disease in children characterized by prolonged synovial inflammation of a joint. The Dynamic Contrast Enhanced Magnetic Resonance Imaging (DCE-MRI) imaging technique, in which contrast fluid is administered during the MRI scan, is the preferred imaging modality by radiologists to quantify disease activity and therapy response in JIA. The radiology department of the Academic Medical Center (AMC) has developed research software (Dynamo) to analyse DCE-MRI images. Dynamo classifies so-called time-intensity curves (TICs) into different types depending on the TIC shape. This classification in its clinically oriented version (Auto-ME-TIC) is used by the radiologist for disease assessment in JIA. The aim of this research is to improve the accuracy of the classification by developing a new enhancement point detection algorithm in which the moment in time when contrast fluid enters each location in the knee is assessed automatically.
Methods - In this research the new enhancement point detection algorithm is evaluated by two evaluators and compared with the algorithms already implemented in Dynamo and Auto-ME-TIC. New classification software (Varius) is developed in which the new enhancement point detection algorithm is implemented and TIC classification is performed. Dynamo, Auto-ME-TIC and Varius are then evaluated by three independent evaluators including a radiologist to determine and compare the accuracy of the classification of these three classification algorithms. In this assessment, the radiologist represents the gold standard, and all algorithms are assessed against his judgement.
Results - The new enhancement point detection algorithm in Varius compared to Dynamo and Auto-ME-TIC is more accurate (with accuracy improvement of up to 91.0% vs. 71.3% with respect to Auto-ME-TIC) whereas the accuracy of the classification in Varius is not higher than the existing classification algorithm (with values of 40.3% vs. 47.9% when compared to the best performing classification, Auto-ME-TIC).
Discussion - The low accuracy of all classification software reveals there might be limits to the current approach in the classification, which in its current implementation makes use of thresholds to segment the TICs. This suggests it is worthwhile to consider other techniques for classification of TICs. A question that remains to be answered is the validity of the radiologist as a reliable reference standard.
Conclusion - The new enhancement point detection algorithm is more accurate than the existing enhancement point detection algorithms. The successful implementation of the new enhancement point detection algorithm in Varius did not improve the accuracy of the classification. The new enhancement point detection algorithm provides an opportunity for further experimentation in the development of software for DCE-MRI analysis.
Keywords: JIA, DCE-MRI, TIC, classification
Introductie - Juveniele idiopatische arthritis (JIA) is een reumatische aandoen-ing bij kinderen gekenmerkt door onstekaandoen-ing van het synovium van een gewricht. De beeldvormingstechniek Dynamische MRI (DCE-MRI), waarbij contrastvloeistof wordt toegediend tijdens de MRI scan, heeft de voorkeur onder radiologen om de ziekteactiviteit en therapierespons in JIA vast te stellen. De radiologieafdeling van het AMC heeft voor onderzoeksdoeleinden software ontwikkeld (Dynamo) om DCE-MRI beelden te analyseren. Dynamo classificeert zogenaamde signaal-intensiteitscurves (TICs) in verschillende typen afhankelijk van de vorm van de TIC. Deze classificatie in een klinisch georiënteerde setting (Auto-ME-TIC) wordt gebruikt door radiologen om ziekteactiviteit van JIA vast te stellen. Het doel van dit onderzoek is om de accuraatheid van de classificatie te verhogen door een nieuw enhancement punt detectie algoritme te ontwikkelen dat automatisch voor iedere locatie in de knie vaststelt wanneer contrastvloeistof zich in het kniegewricht bevindt.
Methoden - In dit onderzoek is het nieuwe enhancement punt detectie algoritme geevalueerd door twee beoordelaars en vergeleken met de algoritmen die al waren geïmplementeerd in Dynamo en Auto-ME-TIC. Nieuwe classificatie software (Varius) is ontwikkeld waarin het nieuwe enhancement punt detectie algoritme is geïmple-menteerd en TIC classificatie heeft plaatsgevonden. Dynamo, Auto-ME-TIC en Varius zijn vervolgens geëvalueerd door drie onafhankelijke beoordelaars waaronder een radioloog om de accuraatheid van de classificatie te beoordelen voor deze drie classificatie algoritmen. Bij deze beoordeling representeert de radioloog de gouden standaard, en zijn alle algoritmen geëvalueerd naar zijn oordeel.
Resultaten - Het nieuwe enhancement punt detectie algoritme in Varius vergeleken met Dynamo en Auto-ME-TIC is accurater (met een maximale verbetering in de accuraatheid tot 91.0% vs. 71.3% volgens Auto-ME-TIC) terwijl de accuraatheid van de classificatie in Varius niet hoger is dan het huidige classificatie algoritme (accuraatheid van 40.3% vs. 47.9% vergeleken met de best presterende classificatie Auto-ME-TIC).
Discussie - De lage accuraatheid van alle classificatie algoritmen laat zien dat er mogelijk beperkingen zijn in de huidige aanpak van classificatie, waarin in de huidige implementatie gebruik wordt gemaakt van thresholds om TICs te segmenteren. Dit suggereert dat het de moeite waard is om andere technieken voor classificatie te overwegen. Een vraag die onbeantwoord blijft is of de radioloog als referentie standaard geschikt is.
Conclusie - Het nieuwe enhancement punt detectie algoritme is accurater dan de bestaande enhancement punt detectie algoritmen. De successvolle implementatie van het nieuwe enhancement punt detectie algoritme in Varius heeft de accuraatheid van de classificatie niet verhoogd. Dit nieuwe algoritme vormt een mogelijkheid om verder te experimenteren in de ontwikkeling van software voor DCE-MRI analyse. Trefwoorden: JIA, DCE-MRI, TIC, classificatie
Chapter
1
Introduction
Since 2007 the radiology department of the Academic Medical Center (AMC) makes use of software to assess disease activity for juvenile idiopathic arthritis (JIA). The images that are being used by the radiologist are created by a Magnetic Resonance Imaging (MRI) technique, also called DCE-MRI, in which contrast fluid is administered to the patient. The aim of this research is to improve the algorithm used in the DCE-MRI analysis software which classifies so-called time-intensity curves (TICs) in different types. In this chapter the necessary background information is provided regarding JIA, (DCE)-MRI, enhancement point detection, and TIC types. The rest of this chapter describes our research questions.
1.1
A description of Juvenile Idiopathic Arthritis
Explanation of the diseaseJIA is a type of rheumatic disease that has an onset before the age of 16 for a minimal duration of 6 weeks. JIA is the most prevalent cause of chronic joint inflammation in childhood (C. M. Nusman et al.(2013)). The prevalence in the western world is estimated to be 1 in 1000 children (Ravelli & Martini(2007)). JIA is characterized by prolonged synovial inflammation of the joint, depicted in Fig. 1.1. The left image shows a healthy joint. The right image shows where the joint is affected by arthritis resulting in bone growth, inflammation of the synovial membrane, thinning of the cartilage and presence of excess synovial fluid.
Figure 1.1: Anatomical differences of a normal joint next to a joint affected by JIA. Reprinted and modified from (Weiss & Ilowite(2005)).
Diagnosis of JIA
JIA diagnosis is performed as a diagnosis-per-exclusionem (i.e. a diagnosis after having excluded a set of other diseases (C. M. Nusman et al. (2013))). Traditionally a diagnosis of JIA is given after anamnesis and physical examination is conducted. However, the possibility of subtle presence of arthritis, signs of joint swelling and limitation of motion (used for diagnosis) appear to be difficult to assess in a physical examination (Guzman et al.(1995)). Therefore imaging methods play an increasingly important role in the diagnosis of JIA.
In addition to anamnesis and physical examination, both ultrasound and MRI are being used. However, MRI is the only imaging modality that is able to show objectively and without the use of ionizing radiation the manifestations of JIA with all its relevant structures (C. M. Nusman et al. (2013)). Conventional radiography (e.g. X-rays) as additional diagnostic evidence was used in the past but nowadays has low value in contributing to the diagnosis of JIA.
Treatment of JIA
The treatment of patients with JIA aims at clinical remission by prescribing medication (often in combination with occupational therapy, physical therapy and psychosocial support) (C. M. Nus-man et al.(2013)). To achieve clinical remission the condition of inactive disease activity must be met. A number of different types of medication exist to achieve inactive disease activity, varying from disease-modifying anti-rheumatic drugs (DMARDs) to biologicals that inhibit activity of the immune system, often in combination with IAGCs (Intra-articular glucocorticoids) or non-steroidal anti-inflammatory drugs (NSAIDs) that decrease inflammation caused by the immune system (Benjamin & Lappin (2019);Saad & Pellegrini (2019)).
Monitoring disease activity and therapy response in JIA
To evaluate therapy response various validated outcome measures are available such as the Visual Analog Scale (VAS) score, the Juvenile Arthritis Multidimensional Assessment Report (JAMAR) (Martin et al.(2018)) and the Juvenile Arthritis Disease Activity Score (JADAS) (Nordal et al. (2011)). There is also the possibility of measurement of clinically related biomarkers for prediction of disease activity and severity.
Another objective method to evaluate therapy response is by making use of Dynamic contrast-enhanced Magnetic Resonance Imaging (DCE-MRI). DCE-MRI is a method used to quantitatively evaluate disease activity and therapy response in children with Juvenile Idiopathic Arthritis. (Hemke et al. (2014,2017);Boesen et al. (2018)).
1.2
Magnetic Resonance Imaging
Magnetic resonance imaging (MRI) is a medical imaging technique widely used in radiology. An MRI scanner makes use of magnetic fields and radio waves to form pictures of the anatomy and physiological processes of the body. The strength of a magnetic field is described in Tesla (T). A common magnetic field strength is 1.5T and 3T. In Fig. 1.4 a typical MRI scanner of 1.5T is shown. Nowadays MRI scanners exist with a magnetic field strength of 7.0T. MRI does not involve ionizing radiation which distinguishes it from other imaging modalities such as Computed Tomography (CT) or Positron-emission tomography (PET). MRI is the preferred method for imaging examination on pediatric patients since risks by exposure to radiation are higher for children compared to adults (Portelli et al. (2016)).
Figure 1.2: The alignment of hydrogen protons. Reprinted and modified from (Berger(2002)).
Figure 1.3: MRI image of one slice of the knee in axial orientation of a patient with JIA.
In an external magnetic field the magnetic movement, also described as spin, of the protons in the nuclei of hydrogen atoms align to the main magnetic field1 B as depicted in Fig.~ 1.2
(Berger (2002)). MRI is an imaging modality that is able to measure a signal from the hydrogen proton spins (Philips Medical Systems(1999)). An example of an MRI image of a slice of a knee in axial orientation of a patient with JIA is shown in Fig. 1.3.
1.3
Dynamic Contrast Enhanced Magnetic Resonance Imaging
The article from Walsh & Pearson(2001) states that the more inflammation is present at a joint, the more vascularized this tissue area can become. Angiogenesis and chronic inflammation are codependent (Jackson et al.(1997)). By measuring the degree of vascularization (also known as angiogenesis) by examining signal changes in this tissue over time one indirectly measures the degree of inflammation (Costa et al.(2007)).
Figure 1.4: Overview of the imaging modality DCE-MRI.
1
“The hydrogen proton can be likened to the planet earth, spinning on its axis, with a north-south pole. In this respect it behaves like a small bar magnet. Under normal circumstances, these hydrogen proton “bar magnets” spin in the body with their axes randomly aligned. When the body is placed in a strong magnetic field, such as an MRI scanner, the protons’ spin axes tend to line up. This uniform alignment creates a magnetic vector oriented along the axis of the MRI scanner” modified fromBerger(2002).
Figure 1.5: Time-Intensity Curves (TICs) from DCE-MRI data with various shapes. The x-axes represents the number of dynamic scans (e.g. acquisition points or time points) that are measurements of signal intensity through time, the y-axes represents signal intensity (SI).
The main difference of DCE-MRI compared to conventional MRI is that the scan is repeated in time during and after the administration of a contrast agent. DCE-MRI is schematically represented in Fig. 1.4. Contrast agent influences the T1 magnetic properties of hydrogen proton spins in their vicinity. In DCE-MRI it is common to use gadolinium (Gd) as contrast agent. After administration the contrast agent enters the bloodstream and after a while the contrast agent reaches the joint. Before wash-in, during wash-in and during wash-out stages2
signal changes in the tissue are measured over time resulting in a so-called time-intensity curve (TIC) (Lavini et al. (2007,2013); Hemke et al.(2014, 2017)) as shown in Fig. 1.5. The TIC (or signal-time curve) is defined as the measured Magnetic Resonance (MR) signal in a pixel or region of pixels over time. Several shapes of TICs that occur frequently have been identified, where some types correspond to a particular degree of disease activity.
Pixel-by-pixel analysis
In the beginning research in DCE-MRI analysis has mostly focused on a Region Of Interrest (ROI) approach in which the TIC-curve was created by averaging all ROI pixel values assessed by the radiologist by looking at it and then making a selection of pixels (Reddick et al.(1999); Kiessling et al. (2005)). Since 2007 a new approach has been introduced in which each pixel is analysed individually by means of an algorithm. Each pixel is classified in one of seven pre-defined shape types. The benefits of the pixel-by-pixel approach have been described in (Lavini et al.(2007); O. A. Kubassova et al.(2007)). For these two reasons the pixel-by-pixel
approach has gained more interest among researchers.
1. Lesions are typically characterized by more than one TIC type. The ROI approach is not able to discover multiple types within a lesion.
2. The pixel-by-pixel approach lends itself better for an analysis of the statistical distributions of the TIC-types that are present in the image.
At this moment there are many techniques available to analyze DCE-MRI. DCE-MRI analysis may be performed in a fully quantitative manner using pharmacokinetic modeling (PKM) which describes the changes in contrast agent concentration derived from observed changes in the signal intensity (SI) (Jackson et al. (1997); Hodgson et al. (2007)) based on a mathematical
2
Before wash-in: before the contrast agent has reached the target tissue; during wash-in: when contrast agent reached the target tissue; during wash-out: when contrast agent is leaving the tissue.
model. However, the most accessible and efficient method to analyze DCE-MRI currently used in clinical practice is a semi-quantitative analysis of the TIC-shape by analyzing parameters derived from the TIC curve. In this research the main technique is a semi-quantitative analysis that involves the calculation of so-called TIC parameters and curve-fitting.
Dynamo and Auto-ME-TIC
Dynamo is the name of the first version of the TIC-classification software developed by C. Lavini (Lavini et al.(2007)). An example of the user interface of the Dynamo software is found at AppendixA.1.
Many parameters can be calculated and estimated. An example of a parameter is the maximum enhancement (ME). Other parameters such as the Time-To-Peak (TTP) and mean slope of increase (MSI) are calculated for every TIC. A full list of parameters is found at AppendixA.8. Dynamo was used in a research setting from 2004 onward. During this period important changes have taken place such as an increase of the spatial and temporal resolution3 in the scanning protocol of the MRI-images.
The implementation of Dynamo in a clinical setting with a modification in the enhancement point detection has been given the name Auto-ME-TIC4. Auto-ME-TIC was implemented by the Medical Imaging Quantification Center (MIQC) and operational as of 2017. An example of the user interface is found at Appendix A.2.
Enhancement point detection
In this research the enhancement point is defined as the last point in time before the administered contrast agent enters each location in the knee as described by point A in Fig. 1.6. Note that Dynamo and Auto-ME-TIC are developed to identify the point located one to the right of the enhancement point as described by point B. As there is only one point difference between point A and point B, it is possible to directly compare the two enhancement point detection algorithms.
Figure 1.6: Visualisation of the enhancement point detection of TIC. The TIC is classified as Type 6. Type 6 is a special Type of TIC which represents a location in the artery. The x-axis represents the number of time points, the y-axis represents the signal intensity.
3
The number of dynamic scans. 4
Throughout this thesis Dynamo will be referred to as the original software used in research and Auto-ME-TIC as the modified version currently in use by radiologists in the clinic of the AMC.
Clinical relevance of TIC-shape classification
The article Lavini et al.(2007) states that for JIA certain types of TICs correspond to a certain degree of disease activity. At this moment seven different TIC shapes are used by Auto-ME-TIC as shown in Fig. 1.7. After performing the pixel-by-pixel approach and classifying each TIC into a specific type, a classification map is generated as shown in Fig. 1.8. The clinical relevance of shape types may be different between different anatomical locations or between different types of disease.
Figure 1.7: Classification of TICs. Type 1 (gray): no enhancement; Type 2(green): slow enhancement, maximum of the curve is reached after half scan; Type 3 (blue): quick enhancement, followed by a signal plateau; Type 4 (magenta): fast enhancement and quick washout; Type 5 (yellow): quick enhancement, followed by a slow constant
enhancement; Type 6 (red): artery; Type 7 (white/light gray): all others.
The clinical meaning of shape types for the knee for JIA is as follows.
Type 1: Represents tissue where no enhancement takes place. This means contrast agent has not reached this tissue or the signal contained too much noise.
Type 2: The clinical relevance of a slow and steadily increasing wash-in is unknown. Type 3: The clinical relevance of a quick wash-in with a horizontal wash-out stage is unknown.
Type 4: A quick wash-in and a quick wash-out indicates an increase in disease activity (Hemke et al. (2014)).
Type 5: A quick wash-in followed by a steady rise in signal after the peak indicates lowered disease activity (C. Nusman et al.(2017)).
Type 6: The artery containing a steep wash-in and steep wash-out. The artery has a quick wash-in and a quick wash-out because the artery is the first pass before the dilution of contrast agent in blood takes place.
Type 7: All other TICs with no particular shape.
The Maximum Enhancement (ME) parameter derived from the TIC shape proves a reliable marker for predicting joint inflammation for both the knee and the wrist (C. M. Nusman et al. (2017)). Both the classification map as well as the ME parameter are consulted by the radiologist
Figure 1.8: Pixel by pixel analysis and color rendering of the time intensity curve shape analysis. Reprinted fromLavini et al.(2013).
1.4
Research questions and aim of study
The Medical Imaging and Quantification Center (MIQC) at the AMC already makes use of an algorithm in Auto-ME-TIC (based on Dynamo) to provide an automatic analysis of the DCE-MRI scans to produce maps of ME and classification maps of time-intensity curves (TICs)5. The classification algorithm6 calculates multiple TIC curve parameters. However, there are some limitations associated with the algorithm in the DCE-MRI analysis software currently in use. A limitation is that Dynamo and Auto-ME-TIC perform the identification of the enhancement point only for a single artery pixel of one slice located in the middle. The enhancement point that is found is then applied to all the pixels in all slices. This is not good enough since it assumes there is no difference in value between pixels between different locations in the knee.
The enhancement point is important as it is used to calculate TIC curve parameters that are used in the classification of TICs. We expect that by improving the detection of the enhancement point, some TIC curve parameters will be estimated more accurately leading to a higher accuracy in classifying TIC shapes into types.
An algorithm already had been developed by C. Lavini to calculate the enhancement point for all pixels but is currently not in use in clinic because the duration of execution is too long. In this thesis the name Dynamo+ will be used as a reference to the Dynamo software with the implementation of an enhancement point algorithm calculating the enhancement point for all types. The enhancement point detection in Auto-ME-TIC originated from Dynamo but has been modified. The radiology department of the AMC has the desire to gain insight in the accuracy of the enhancement point detection of Dynamo, Dynamo+ and Auto-ME-TIC.
At this moment the classification of the seven different types of TICs needs improving because
5Throughout this thesis the term shape will be used as a synonym for a TIC. In this thesis TICs may also be referred to as curves.
6
The classification algorithm is the part of the DCE-MRI analysis software where TIC parameters are calculated and where the Lavini classification takes place using conditional statements.
there are new scanning protocols available with higher temporal resolutions compared to the temporal resolutions used when the software Dynamo was developed. This work will focus on both the improvement of the current method to automatically recognize the point in time where the enhancement begins as well as the required changes and improvements of the calculation of multiple TIC parameters in the classification of TICs. In this research the tissue region of interest is the knee joint as presented in Fig. 1.1.
The accuracy of Dynamo on lower temporal resolution DCE-MRI series is known from previous research (Lavini et al.(2007)). The accuracy of Dynamo on higher resolution DCE-MRI series is unknown. The accuracy of Auto-ME-TIC has never been evaluated. A new evaluation of Dynamo and a first evaluation of Auto-ME-TIC is therefore desired. The TIC-curves will be assessed by a blind evaluation of at least one experienced radiologist. This enables a comparison between the accuracy of Auto-ME-TIC, Dynamo and Varius.
Varius is the new classification software developed in this thesis. Varius consists of the original Lavini classification in Dynamo and Auto-ME-TIC with a novel implementation of a variable enhancement point detection algorithm along with any required modifications in the calculation of TIC curve parameters. These modifications involve a different method to calculate the baseline in a TIC, minor modifications to calculate the mean slope of increase (MSI) (e.g initial slope in a TIC), a different method to calculate the time-to-peak (TTP) and a different method to fit a line to calculate the slope at the end of a TIC. It is desired to find out the accuracy of Varius as to compare its accuracy to Dynamo and Auto-ME-TIC. Our overarching goal is to increase the accuracy of the classification of the TICs in order to provide a more accurate description of disease activity and therapy response of JIA.
Problem Definition: There is a desire by researchers and radiologists to more accurately measure disease activity to monitor the therapy response of patients suffering from JIA. Aim: The aim of this research is to develop an enhancement point detection algorithm for every pixel and every shape type in order to improve the DCE-MRI analysis. TIC curve parameters dependent on the point of enhancement will be calculated based upon the new enhancement point detection algorithm. The enhancement point detection algorithm should determine the position of the last acquisition point before the presence of contrast agent.
Research Questions
1. How can the detection of the enhancement point be improved in such a way that is fast enough for clinical use?
2. What is the accuracy and performance of the new enhancement point detection algorithm and the current enhancement point detection algorithms that detect the enhancement point?
3. What changes are required to TIC curve parameters in order to successfully implement the new enhancement point detection algorithm in the classification algorithm?
4. What is the accuracy and performance of the DCE-MRI analysis with higher temporal resolutions for Auto-ME-TIC, Dynamo and Varius?
1.5
Outline of thesis
In Chapter 2 various algorithms and software modifications developed for this thesis are presented. The algorithm for the detection of an enhancement point and its comparison to the existing enhancement point detection algorithms, the evaluation of the new enhancement point detection algorithm, the evaluation of the enhancement point detection algorithms in Dynamo, Dynamo+ and Auto-ME-TIC, the modifications required to implement these changes resulting in the classification system Varius, and lastly the evaluation of the accuracy and performance of Dynamo, Dynamo+, Auto-ME-TIC and Varius.
In Chapter 3 all results are presented. These include the accuracy of the new enhancement point detection algorithm, the accuracy of the existing enhancement point detection algorithms in Dynamo, Dynamo+ and Auto-ME-TIC, as well as the accuracy in classification of Dynamo, Auto-ME-TIC and Varius. The calculation time is estimated for enhancement point detection and classification by Dynamo, Dynamo+, Auto-ME-TIC and Varius.
In Chapter 4 a discussion is found on the main findings of this research with answers to the research questions, a discussion of our results, the limitations and strengths in this research, the implications of this research for clinic, as well as to provide an agenda for future research.
Chapter
2
Methods
In this chapter the following activities are discussed.
1. Pre-processing steps.
2. The development of the new enhancement point detection algorithm. 3. The evaluation of all enhancement point detection algorithms.
4. The modifications that were applied to TIC curve parameters and required to implement the new enhancement point detection algorithm.
5. The evaluation of the classification of Dynamo, Auto-ME-TIC and Varius.
6. The methods used for calculation of execution time for the enhancement point detection algorithms in Dynamo, Dynamo+, Auto-ME-TIC and Varius.
7. The methods used for calculation of the execution time of the classification by Dynamo, Auto-ME-TIC and Varius.
All algorithms are written in MATLAB R2016a. DCE-MRI images of the knee from eight patients with JIA were used as training data for development. DCE-MRI images of the training data were used for evaluation of the enhancement point detection algorithms. DCE-MRI images of the knee from nine other patients were used for the evaluation of the classification software. In total data of 17 patients were used for development and evaluation. The scan protocol that was used for the DCE-MRI scans is identical to the scan protocol currently used in clinic. Images were acquired with a 3.0T MR imaging equipment (Ingenia from Philips Medical Systems). The scan protocol is found at AppendixA.16. All algorithms with their names and descriptions are listed in Appendix B.12. A list of all thresholds used is found in Appendix B.1.
2.1
Pre-processing
Noise and background filter
A noise and background filter is used to remove noisy pixels in close proximity and further away from the outer border regions of the knee as depicted in Fig. 2.1a. In our set of nine patients the noise filter threshold is set to a signal value of 20. This threshold is a default parameter defined in Dynamo and Auto-ME-TIC. This leads to measured signal values below the threshold of 20 labeled as background and ignored.
Determining noise at a fixed oscillation threshold
An algorithm already existed which determines the allowed oscillation (%) of the TIC in time. Setting a threshold on this parameter results in a segmentation as depicted in Fig. 2.1b. An oscillation threshold is a decimal value between 0.00 and 1.00. An oscillation threshold of 0.20
means an oscillation of up to 20% is allowed for each pixel. The range of an oscillation threshold remains within 0% and 100%. The lower this oscillation threshold, the more pixels will be labeled as noisy since less oscillation is allowed. The higher the oscillation threshold less pixels will be labeled as noisy since more oscillation is allowed.
(a) Noise mask. (b) Oscillation mask. (c) Oscillation map.
Figure 2.1: A noise mask, an oscillation mask and an oscillation map. In Fig. 2.1athe background is red and the area of interest is black. In Fig. 2.1bthe red area represents noise after allowing oscillation up to 20%. In Fig
2.1cthe brighter the color, the more oscillation is present. The darker the color, the less oscillation is present as shown in the color bar.
Determining different noise level oscillation groups below an oscillation threshold for sigma calculation
A new algorithm is written as found in Appendix B.12.1. It is possible to create an oscillation mask repeatedly for many oscillation thresholds. Each oscillation mask has its own oscillation threshold. Combining all oscillation masks each with its own oscillation threshold results in an oscillation map as depicted in Fig. 2.1c. The brighter the color, the more oscillation is present. The darker the color, the less oscillation is present as shown in the color bar. The resulting noise oscillation map, containing pixel groups with different (%) oscillation, is used in the sigma calculation for smoothing.
2.2
Development of an algorithm to automatically recognize the
first point in time where the enhancement begins
The new enhancement point detection algorithm performs the enhancement point detection for all pixels and shape types 2 to 7 which are present in the image. During development two key aspects were taken into consideration. The new enhancement point detection algorithm must be as accurate as possible while at the same time its processing speed must be preserved or even faster compared to the current enhancement point detection algorithms. The existing enhancement point detection algorithms in Dynamo and Auto-ME-TIC are only suited for enhancement point detection of Type 6 and not suited to be applied to pixels of all types present in a slice. The result of the enhancement point calculation for one pixel of Type 6 is applied to all pixels in all slices. The new enhancement point detection algorithm calculates the enhancement point for each and every pixel of the knee joint. Two examples of enhancement point detection are shown in Fig. 2.3. The enhancement point detection algorithm is found at B.12.10.
An algorithm already was created by C.Lavini with the purpose of calculating the enhancement points for all types in a slice as implemented in Dynamo+. This algorithm is found at Appendix
B.12.13. An example of Dynamo+ calculating all the enhancement points in a slice is found in AppendixA.4.
Defining the search space for candidate enhancement points
A new algorithm is written that determines the desired search space. The search space is defined as the range of time points or acquisition points in a TIC which are candidate enhancement points. One of the candidate enhancement points is chosen as enhancement point. In this patient candidate enhancement points may be chosen in the search space with a range between time point two and time point fifteen. This range might slightly differ between different patients. The search space is identical for all slices belonging to the same patient. The implications of this decision on the timing is discussed later in the discussion section. The search space is fixed around a central point in time of the search space window and is defined in seconds. The central point in time is defined as the time point in a TIC where change in signal value f0(x) is highest. Most often this is the center in the initial slope of a TIC.
The center of the search space is found as follows. The first step is to extract the signal values of all slices belonging to the same patient. The next step is to calculate the average TIC of all pixels of each slice. After the mean TIC of every slice is calculated, a mean TIC of all slices of the patient is calculated. The resulting TIC which includes the search space for candidate enhancement points is depicted in Fig. 2.2. The time point is selected where the first order derivative f0(x) is the highest and referred to as the central point in time. The algorithm used for finding the center of the search space is listed in Appendix B.12.2.
In the range of 49 seconds before and 49 seconds after this central point in time the enhance-ment point search is conducted which is visualized as the gray area in Fig. 2.2. This corresponds to about eight time points before and eight time points after the central point in time in our dataset. A constraint is that the search is conducted starting from time point two and above since a baseline is always required. The baseline is defined as a series of points in a TIC during which the injected contrast agent has not yet reached the target tissue. The baseline can in principle have different lengths in each pixel present in the image.
Figure 2.2: An average TIC of the whole slice of all slices belonging to one patient. An example of search space selection of candidate enhancement points. The blue-dotted line is the central point in time. The search space of candidate enhancement points is visualized as a gray area between time point 2 and time point 15.
(a) Enhancement point detection in a TIC Type 3. (b) Enhancement point detection in a TIC Type 6. Figure 2.3: Two examples of the enhancement point detection. The red line represents the signal intensity. The pink line represents the first order derivative f0(x) of the signal intensity and the black line represents the second order derivative f00(x) of the signal intensity. The blue dot and the blue-dotted line represents the location of the enhancement point. Classification was performed by Auto-ME-TIC.
Smoothing and calculating a variable sigma
A new algorithm is written that calculates a variable sigma value as described in equation
2.1 in order to apply a 1D horizontal Gaussian smoothing filter on a TIC as depicted in Fig.
2.4. Using a Gaussian filter it is possible to blur images and remove detail and noise. The Gaussian smoothing filter is performed by the function imgaussfilt of the Image Processing Toolbox available for MATLAB. This function is also suitable for the smoothing of TICs. The smoothed TIC is a TIC upon which a gaussian smoothing filter is applied. The original signal is filtered and sharp edges present as noise are suppressed.
(a) Low smoothing (σ = 0.5). (b) Medium smoothing (σ = 1.0). (c) High smoothing (σ = 1.5). Figure 2.4: Sigma values and their smoothing effect on a single TIC. The x-axes represents the number of time points, the y-axes represents signal intensity. The higher the sigma in equation2.1the higher the smoothing effect that is applied on a TIC.
G(x) =√1 2πσe
−x2
2σ2 (2.1)
where G = the new signal value, x = the distance from the origin in the horizontal axis, σ = the standard deviation of the Gaussian distribution and G(x) = Gaussian smoothing.
The parameter σ described in equation 2.1 is important because the degree of smoothing applied to each individual pixel in the knee is determined by the value of this parameter. The default parameters are a minimum sigma of 0.5 and a maximum sigma of 1.5. In total four ranges of sigma values have been tested. A sigma range of 0.5 to 1.5, a sigma range of 0.75 to 1.75, a sigma range of 1.00 to 2.00 and a sigma range of 1.25 to 2.25.
The parameter σ is automatically adjusted depending on the noise oscillation level of the pixel group of each individual pixel as depicted in Fig. 2.1c. In total 29 noise oscillation pixel groups are currently present given the set of parameters listed in AppendixB.2. An algorithm is written which contains equation 2.2to calculate the resulting sigma. This algorithm is listed in AppendixB.12.7.
S(k) = σmin+ k ·
σmax− σmin
oscthr (2.2)
where k = oscillation values (%), oscthr = oscillation threshold, σmax = the upper bound of the sigma range, σmin the lower bound of the sigma range, S(k) = resulting sigma value. The values in k are derived from the oscillation map. The equation in2.2 is applied to groups of pixels with identical (%) oscillation. The number of noise level groups is calculated from the number of unique noise oscillation values present in the oscillation map (see Fig. 2.1c). A higher sigma leads to more smoothing applied to a group of pixels (containing identical oscillation values (%)) while a lower sigma leads to less smoothing applied to a group of pixels (containing identical oscillation values (%)) as depicted in Fig. 2.5.
By adaptive smoothing of the TICs noise in the knee is suppressed. This leads to a more distinct shape. This shape is less prone to errors. Note that the smoothed TIC is only used in the enhancement point detection and not further on in the classification algorithm.
The author performed a manual visual inspection by-eye of at least five pixels that differed between the resulting enhancement point slices each created from different sigma ranges. After comparing the chosen points by the algorithm with the enhancement point chosen upon by-eye as human observer, the most reliable sigma range was established.
Figure 2.5: Oscillation map and effect of adaptive smoothing. The x-axes of the TICs represents the number of time points that are measurements of signal intensity through time, the y-axes represents signal intensity (SI). Note the different values of oscillation (%) and sigma that were applied to reduce noise. The color bar shows the (%) oscillation whereas a darker color (e.g. blueish) represents lower noise oscillation (%) whereas a brighter color (e.g. reddish) represents higher noise oscillation (%).
Approximating the baseline
After the adaptive smoothing has taken place but before the actual selection of the enhancement point the baseline is approximated as depicted in Fig. 2.6. Only the first half of all acquisition points in a TIC are candidate points over which to calculate a mean signal baseline value as the baseline is expected to be present in the first half of the TIC. A baseline is a signal value that is calculated over a range of time points in a TIC also described as the baseline range. The baseline range is defined as a range of signal values within one standard deviation above and below the mean baseline signal value.
Baseline calculation is performed over a range of time points, between time point one and the time point that belongs to the maximum change in signal value of the first order derivative f0(x). A mean baseline value is calculated over this range of time points. The resulting baseline value is used later on in a conditional statement to determine whether to select the enhancement point determined by the first derivative f ’(x) or the second derivative f ”(x).
Figure 2.6: A description of the baseline in a TIC.
Final selection of the enhancement point
The enhancement point is chosen where the first order derivative of the signal of the smoothed TIC f ’(x) is highest. However, if the signal value of the enhancement point chosen by the first order derivative f ’(x) does not remain within the baseline range as depicted in Fig. 2.6, the enhancement point is chosen where the second order derivative f ”(x) of the signal of the smoothed TIC is highest. This chosen point can be either inside or outside the baseline range. In many situations the enhancement point of the second order derivative f ”(x) is located one point to the left from the enhancement point chosen by the first order derivative f ’(x).
Overview of the enhancement point detection
The whole process of enhancement point detection is summarized below.
1. Apply a noise filter to remove noise around the knee and the outer edges as depicted in Fig. 2.1a. This mask is the background and will be ignored.
2. Apply a noise oscillation filter to remove severe noise levels equal to or above an oscillation threshold as depicted in Fig. 2.1b. This step is applied many times in order to successfully perform step3.
3. Generate a noise oscillation map to quantify the noise under the oscillation threshold as depicted in Fig. 2.1c.
4. Define the search space of candidate enhancement points as depicted in Fig. 2.2.
5. Calculate a sigma that is dependent on the (%) noise oscillation for different groups of pixels.
6. Apply the correct level of smoothing for the pixel groups depending on their noise levels as depicted in Fig. 2.5.
7. Calculate the first order f ’(x) and second order f ”(x) derivative of the TICs and find the maximum of both derivatives as depicted in Fig. 2.3. This results in two time points in a TIC where change in signal intensity is highest.
8. Calculate the baseline over the baseline range. Calculate an upper and lower boundary of one standard deviation above and below the baseline representing the expected maximum and minimum baseline signal values as depicted in Fig. 2.6.
9. Set the enhancement point as the point belonging to the maximum change in signal value of the first order derivative f ’(x). In the situation that the chosen point has a signal value above or below the baseline range, the enhancement point of the second order f ”(x) derivative is chosen.
2.3
Evaluation of the enhancement point detection algorithm
in Dynamo, Dynamo+, Auto-ME-TIC and Varius
The enhancement point detection algorithm is evaluated by two methods. After enhancement point detection calculation is completed for all pixels of one slice, all TICs having been assigned and identified by the algorithm with the same enhancement point were grouped together and plotted in a graph showing individual TICs divided by type and a graph showing an average TIC of each candidate enhancement point for all types. The individual TIC is defined as one of many TICs that has been assigned the same enhancement point as all the other individual TICs in a slice. The average TIC is defined as a TIC that is created by calculating the mean of all individual TICs of all types in a slice that have been assigned with the same enhancement point. TICs were classified by Auto-ME-TIC in order to distinguish between types.
The first evaluation method is a visual inspection by-eye by the author of averaged TIC graphs and individual TIC graphs to inspect the accuracy of the new enhancement point detection algorithm. The location of the enhancement point detection on the averaged TIC graph and the individual TIC graph is checked against the chosen enhancement point (by-eye) by the author. When most enhancement points on the averaged TIC graph or individual TIC graph are correctly chosen, the enhancement point that was detected by the enhancement point detection algorithm should correspond to the enhancement point chosen (by-eye) by the author. The individual TIC graphs, showing individual TICs having been assigned by the algorithm with identical enhancement points, were plotted in order to detect individual cases of incorrect enhancement point detection.
The second evaluation method consists of manual labeling of the enhancement points of a set of TICs, by-eye, by C. Lavini and the author. A dataset of in total 384 enhancement points are evaluated for Type 2 to 7. The dataset was created from eight different patients from the training set. The set contains shape Type 2 to shape Type 7 in equal proportions. Classification was performed by Auto-ME-TIC for an even proportion of shape types. TICs were selected from slice 9 for all eight patients. Slice 9 is exactly in the middle or center of the image stack for 6 out of 8 patients. For 2 out of 8 patients slice 9 is located slightly below the middle.
After calculation of the enhancement points for Dynamo, Dynamo+ and Auto-ME-TIC the calculated point (as described by point B in Fig. 1.6) was subtracted by one as to describe point A (see point A in Fig. 1.6), so that results may be compared with Varius as Varius is developed to search for Point A. After removal of TICs that were considered too noisy 326 TICs were left for further analysis. Each evaluator made his or her own decision when a TIC was considered too noisy to be able to evaluate.
2.4
Classification procedure by Dynamo, Auto-ME-TIC and
Varius
The classification procedure in Varius is based on the original classification (also referred to as the Lavini classification) used in Dynamo and Auto-ME-TIC. The Lavini classification is the last part of the algorithm where TICs are classified to a Type by executing a set of conditional statements as found in Appendix A.17. The Lavini classification is used in Dynamo, Auto-ME-TIC and Varius. In contrast to Dynamo and Auto-ME-TIC where the enhancement point is fixed for all pixels in all slices, Varius has a parameter representing the enhancement point which is variable for each pixel. No changes to the Lavini classification procedure were made except for successful implementation of the new enhancement point detection algorithm. The following input parameters and important thresholds used in Dynamo, Auto-ME-TIC and Varius are set by default. A list of input parameters and thresholds is presented in AppendixB.1 andB.3. A multitude of TIC curve parameters are calculated. An overview and description of all TIC curve parameters used to classify a TIC shape into a Type is found at Fig. 2.7 (Lavini et al. (2007)). A complete overview of all parameters calculated for one slice is found at Appendix A.8.
Figure 2.7: A description of the curve parameters that are used in the classification. Parameters used for the features definitions: SB, signal baseline; MSD, maximum signal difference; TTP, time to peak; MSI, maximum slope of increase; SM, signal maximum; α and β are the intercept and tangent of the slope fitting the tail of the curve, respectively (Lavini et al.(2007)).
1. Baseline (SB)*: The average signal intensity at a given location in the knee before the enhancement takes place.
2. Enhancement point: The point in time when the contrast agent reaches each location in the knee. This point is described as the last baseline point for each pixel in a TIC curve.
3. Maximum Signal Difference (MSD): The difference in signal between the maximum (e.g. top) the curve and the baseline (SB).
4. Maximum Enhancement (ME): The ME is defined as the MSD divided by the baseline. 5. Maximum Slope of Increase (MSI)*: The maximum slope during enhancement.
6. Time-to-peak (TTP)*: The time difference between the moment where the ME occurs and the enhance-ment begins.
7. Tangent*: Tangent of the slope β fitting the tail of the curve. 8. Intercept*: Intercept of the line fitting the tail of the curve. 9. Relative tangent: Tangent divided by the MSD.
10. Relative ratio: MSD divided by (Intercept - Baseline).
* This parameter is modified to implement the new enhancement point detection algorithm in Varius. The modification of this parameter is discussed at section2.5.
(a) Baseline (SB). (b) Enhancement points. (c) Mean Signal Difference (MSD). Figure 2.8: Three curve parameters used in the classification. For Fig. 2.8a and Fig. 2.8cthe darker (e.g. blueish) the lower the signal intensity, the brighter (e.g. reddish) the higher the signal intensity. For Fig. 2.8bthe brighter (e.g. more red) the color the more early in time the enhancement took place at that pixel whereas the darker the color (e.g. more blue) the later in time the enhancement took place at that pixel.
A complete overview and visualization of all TIC curve parameters is found at AppendixA.8. The classification is executed by fulfilling a set of conditions for each shape type for each individual pixel. In order to classify a TIC shape in Type 6, a list of conditional statements have to be fulfilled as shown in Fig. 2.9. A full overview of the whole classification in Varius including all conditions derived from the Lavini classification is found at Appendix A.17.
Type 6 relenhancementa ≥ histo_thrd slopyb > MSD / 2.0 TTP ≤ enhpoints+expectedPeakWidthc TTP ≥ enhpoints - 2 reltangent < tan_thrd relratio ≥ rrat_thrd
a This is the ME parameter. This parameter describes the relative enhancement. b This is the MSI parameter.
c This is a parameter that describes the estimated width of the peak in a TIC as a number of time points. d A description of the thresholds is given in AppendixB.1.
Figure 2.9: Classification of Type 6. All conditions that have to be fulfilled to classify a TIC shape in Type 6.
After all conditions in Fig. 2.9are fulfilled the pixel will be labeled as Type 6. In a comparable manner other types are classified. We can classify and distinguish types by inclusion and exclusion of certain conditional statements.
2.5
Development and modifications of TIC curve parameters
Additional changes have been performed to the classification algorithm in order to successfully implement the new enhancement point detection algorithm in Varius. An overview of the required modifications applied to TIC curve parameters are listed below.
Implementation of a method to calculate the baseline (SB) parameter
A new algorithm is developed to calculate the signal baseline (SB) over a variable baseline distance instead of a fixed distance. Calculating a baseline variable for each pixel is possible
(a) Baseline in Auto-ME-TIC (b) Baseline in Varius
Figure 2.10: The average signal intensities of the baseline for every pixel from Auto-ME-TIC and Varius of slice 17 of a patient with JIA. The darker (e.g. blueish) the lower the average signal intensity of the baseline points, the brighter (e.g. reddish) the higher the average signal intensity of the baseline points. The baseline in2.10ais fixed and for2.10bvariable for every pixel.
since the enhancement point in the new algorithm is variable. The baseline parameter is used in the classification as depicted in Fig. 2.10.
Mean Slope of Increase (MSI) modifications
The current algorithm to calculate the MSI is modified to be able to work with a variable enhancement point. The function of the modified slope algorithm with implementation of a variable enhancement point is found at AppendixB.12.9. The slope calculation itself has not changed. The slope is calculated by measuring the difference in signal intensity between two time points in a TIC. The range over which to slide a window along the x-axis of a TIC is modified to accommodate a variable enhancement point starting at time point two. This is more elaborately visualized and described in AppendixA.10. The slope is calculated as the difference in signal values between two time points in a TIC. The higher the difference in signal value, the larger the slope. The MSI parameter is visualized as shown in Fig. 2.11.
(a) Slope in Auto-ME-TIC. (b) Slope in Varius.
Figure 2.11: Two different slopes from Auto-ME-TIC and Varius of slice 17 of a patient. The darker (e.g. blueish) the lower the MSI, the brighter (e.g. reddish) the higher the MSI.
Calculation of the Time-To-Peak parameter
A new algorithm has been developed to calculate the TTP parameter. The TTP parameter in Dynamo and in Auto-ME-TIC is calculated differently. The TTT is defined as the number of time points in a TIC starting at tstart = 1 and ending at tend where tend has the highest
signal value (e.g. the top or maximum). This is illustrated in Appendix A.9. Dynamo and Auto-ME-TIC calculate the TTT whereas Varius calculates the TTP. The TTT curve parameter is used in Dynamo and Auto-ME-TIC. The TTP curve parameter is used in the classification software Varius.
(a) TTT in Dynamo. (b) TTP in Varius.
Figure 2.12: TTP calculation from Auto-ME-TIC and Varius of slice 17 of a patient with JIA. The brighter (e.g. reddish) the earlier the TTP, the darker (e.g. blueish) the later the TTP.
Implementation of a method to fit a line to calculate a tangent over a variable distance
At this moment Dynamo and Auto-ME-TIC fit a line to calculate the tangent over a fixed distance. A more elaborate visualization is found at Appendix A.11. The TIC is divided into three parts. The third part is used for fitting a line. A tangent is calculated from this fitted line. The following linear equation is solved to calculate the intercept and tangent: y = a + b · x where y = the signal value(s), a = the intercept, b = the tangent of the slope and x = the slope β.
After solving the linear equation the tangent is projected by a multiplier of 5 minutes to calculate the whole tangent distance. An algorithm has been developed and implemented in Varius to fit a tangent over a variable distance instead of a fixed distance at the end of a TIC curve as depicted in Fig. 2.13. Tangents with different lengths are identified and grouped together. Every group of tangents is fitted separately. The first time point in a TIC where the fit of the line starts is calculated follows: sindex= enhpoint + 0.15 · tres where sindex = the time point in a TIC where the fit of the line starts and tres = the number of time points in a TIC. The last time point of the fit of the line is always the last time point of the TIC.
(a) Tangent fit (distance = 24). (b) Tangent fit (distance = 30). (c) Tangent fit (distance = 32). Figure 2.13: The fitted line (in red) for three different TICs for tangent calculation in Varius. Each fitted line describes the slope at the end of a TIC. The time point sindexis the time point where the fitted line begins.
2.6
Evaluation of Dynamo, Auto-ME-TIC and Varius
MATLAB was used to develop a script which visually displays a TIC in its original state. TICs from Type 2 to 7 were shown to the user in a pseudo-random fashion. A snapshot of the evaluation is shown in AppendixA.12. The user is asked for input and enters the Type as a number in the console. All details are saved as MATLAB files and written to a text file for later analysis. A snapshot of how the evaluation took place is found at Appendix A.12.
Three evaluators performed in the evaluation. The first participant R.Hemke is a radiologist. This radiologist is an end-user of Auto-ME-TIC and has performed additional research on DCE-MRI. The second participant is C. Lavini who has written the original classification software Dynamo. Lavini has been involved for a long period in this research. The third participant J. van den Noort is involved in the clinical implementation of Auto-ME-TIC.
All three evaluators classified the same 702 TIC shapes from Type 2 to 7 from nine different JIA patients of the evaluation set. In total 75 TICs were excluded because one or more evaluators classified these TICs as too noisy to be able to classify them. In total 627 TICs remained for further analysis for evaluation of the accuracy of the software. Pixels were chosen and types were evenly distributed according to the classification of Auto-ME-TIC. In total 50% of all slices in the center region of the joint of the image were selected. In the situation that an image contains 28 slices TICs would be pseudo-randomly selected in each slice from slice 7 to 21. The reason for this selection method is that most inflammation is expected to appear in the center region of the image.
Pixels in each slice were pseudo-randomly selected after having applied the noise filters as depicted in Fig. 2.1. The training set was kept separate from the evaluation set as to prevent the possibility of showing pixels that were used during development. TICs were shown with their original signal intensity. The algorithm used for pixel selection is found at AppendixB.12.18. Part of the analysis is a calculation of overlap (e.g. matches (%) and mismatches (%)) between evaluators as well as a calculation of a kappa inter rater agreement.
2.7
Run time of the enhancement point detection algorithms in
Dynamo, Dynamo+, Auto-ME-TIC and Varius
MATLAB contains several methods to evaluate run time. Three different methods were used. 1. The tic and toc commands are used to measure run time over a number of lines of code
located in a script or function.
2. The command timeit measures the time required to run a function. This function performs multiple runs of the input function and calculates a median time in seconds of the measurements. Timeit is considered more accurate than the tic-toc measurement. 3. The profile command is used to track execution time and provides a global overview of
Chapter
3
Results
In this thesis the following objectives were investigated: whether the development of a new enhancement point detection algorithm for all types was achievable, how the accuracy of the new enhancement point detection algorithm would compare with the enhancement point detection algorithms currently in use, what the impact of this new enhancement point detection algorithm has on accuracy of the DCE-MRI analysis (Varius) and the duration of the execution time of the enhancement point detection algorithms and the classification software.
3.1
Development of the new enhancement point detection
algo-rithm
The resulting enhancement point calculation of a whole slice of a patient with JIA is depicted in Fig. 3.1. Dynamo and Auto-ME-TIC chose enhancement point 5.
Figure 3.1: The enhancement points of slice 17 of a patient with JIA from the training set. The color of the color bar corresponds to the data point in time when the enhancement occurs. The brighter (e.g. more red) the color the more early in time the enhancement took place at that pixel whereas the darker the color (e.g. more blue) the later in time the enhancement took place at that pixel.
Calculation of a variable sigma for smoothing
In total four different ranges of sigma values as used in equations 2.1 and2.2 were applied in an effort to reduce noise before enhancement point detection. Two sigma ranges are depicted in Fig.
3.2. During visual inspection of the TICs of at least five pixels the author manually inspected each pixel by-eye and chose an enhancement point and checked whether or not the calculated enhancement point by Varius was identical. The chosen sigma range of σmin = 0.5 and σmax = 1.5 was established as most reliable.
(a) Sigma 0.50 to 1.50. (b) Sigma 1.25 to 2.25.
Figure 3.2: The effect of different sigma minima and maxima on the enhancement point detection. The brighter (e.g. more red) the color the more early in time the enhancement took place at that pixel whereas the darker the
color (e.g. more blue) the later in time the enhancement took place at that pixel.
Impact on classification of choosing one point to the left or one point to the right from the chosen enhancement point
The impact on the final classification by deviating one point is shown in Fig. 3.3. Purposely deviating one point from the chosen point in the new enhancement point detection algorithm, the overlap in Fig. 3.3bcompared to Fig. 3.3a, and Fig. 3.3bcompared to Fig. 3.3cwas still 93%. The highest impact of being one point off remains on Type 6 since the artery was more difficult to detect as shown in Fig. 3.3c.
(a) Point left. (b) Chosen enhancement point. (c) Point right.
Figure 3.3: Impact on the classification by Varius by deviating one point. The color bar displays the color that belongs to a TIC-shape Type.
3.2
Evaluation of the enhancement point detection algorithm
in Dynamo, Dynamo+, Auto-ME-TIC and Varius
Accuracy of Dynamo, Auto-ME-TIC and Varius
A verification by-eye by the author is performed for the new enhancement point detection algorithm in Varius. As Fig. 3.4 and Fig. 3.5 show in most situations the TICs show the same enhancement point as chosen by the author. A more thorough analysis of all candidate enhancement points for one slice is found at Appendix A.13,A.14.
(a) Enhancement point 4. (b) Enhancement point 8. (c) Enhancement point 12. Figure 3.4: Three examples (slice 14, the middle) of averaged TIC graphs of the evaluation of enhancement point detection of one slice of one patient. The x-axis represents the number of time points, the y-axis represents the signal intensity. Every TIC is a calculated average of all TICs showing the same enhancement point.
(a) Type 2. (b) Type 3. (c) Type 4.
Figure 3.5: Three examples (slice 14, the middle) of the evaluation of enhancement point detection of individual TICs for Type 2, Type 3 and Type 4 of one slice of one patient. The x-axis represents the number of time points, the y-axis represents the signal intensity. Every TIC is plotted separately. This evaluation is performed for every combination of chosen enhancement point and Type. This is the evaluation for enhancement point 6.
All enhancement point algorithms were evaluated by C. Lavini and the author. A complete table with the accuracy of enhancement point detection for each and every type by Dynamo, Auto-ME-TIC and Varius is found at Appendix B.4, Appendix B.6 and Appendix B.7. A summary of the results is found at Table 3.1and Table3.2. The accuracy of the enhancement point detection algorithm in Dynamo and Auto-ME-TIC reveal that the modifications that were applied in the enhancement point detection algorithm in Auto-ME-TIC proved to have been beneficial.
Table 3.1: Accuracy of Type 6 of the enhancement point detection algorithm in Dynamo, Auto-ME-TIC and Varius. Type 6 Exact match by C. Lavini (%) Exact match by M.van Kuik (%) One point off by C.Lavini One point off by M. van Kuik Mean exact match (%) Mean one point off (%) Dynamo (N = 326) 41.0 37.7 60.7 59.0 39.4 59.9 Auto-ME-TIC (N = 326) 44.3 41.0 72.1 70.5 42.7 71.3 Varius (N = 326) 34.4 34.4 91.8 90.2 34.4 91.0
As shown in Table 3.1the exact match between both evaluators and Dynamo is relatively poor with an exact match accuracy of 39.4% for Type 6. Allowing to deviate one point increased the match to 59.9%. The exact match for Type 6 for Auto-ME-TIC performs better with an exact accuracy of 42.7% which increases to 71.3% if allowed to deviate one point. The exact match for Type 6 for Varius performs poor with 34.4% but increases to 91.0% if allowed to deviate one point. The Table in AppendixB.7andB.5 show that the enhancement point detection for Type 6 stands out as relatively well-detected compared to the other types.
Table 3.2: Accuracy of all Types of the enhancement point detection algorithm in Dynamo+ and Varius.
All Types (2 to 7) Exact match by C. Lavini (%) Exact match by M.van Kuik (%) One point off by C.Lavini One point off by M. van Kuik Mean exact match (%) Mean one point off (%) Dynamo+ (N = 326) 25.1 24.2 56.4 53.0 24.7 54.7 Varius (N = 326) 40.5 33.5 81.4 74.6 37.0 78.0
The accuracy of Dynamo+ for an exact match over all Types is 24.7% which increases to 54.7% if allowed to deviate one point. Varius has an exact match accuracy over all Types of 37.0% which increases to 78.0% if allowed to deviate one point.
3.3
Classification by Dynamo, Auto-ME-TIC and Varius
Dynamo, Auto-ME-TIC and Varius have all been integrated in order to simultaneously classify the same selection of slices and compare their results. An example of the classification map of one slice created by Varius is shown in Fig. 3.6. The classification of one slice by Dynamo, Auto-ME-TIC and Varius is shown in Fig. 3.7. At first sight a larger variety of types are present in Auto-ME-TIC and Varius when compared to Dynamo. Also, the classification in Dynamo shows a relative large area of Type 2 compared to the the other classification software. Auto-ME-TIC seems to classify more areas in Type 5 compared to the other classification software. A difficulty that already existed in Dynamo and Auto-ME-TIC, and also persists on Varius, is that Type 4 and Type 6 are sometimes difficult to distinguish from each other.
Figure 3.6: Classification by Varius of slice 17 of a patient with JIA. The color bar displays the color that belongs to a TIC-shape Type.
(a) Dynamo. (b) Auto-ME-TIC. (c) Varius.
Figure 3.7: The classification of slice 11 of a patient with JIA performed by Dynamo, Auto-ME-TIC and Varius. The color bar displays the color that belongs to a TIC-shape Type.
3.4
Evaluation of Dynamo, Auto-ME-TIC and Varius
A second evaluation of Dynamo and a first evaluation of Auto-ME-TIC and Varius is performed to compare the overlap (e.g. matches (%) and mismatches (%)) in classification and gain insight in the accuracy of Dynamo, Auto-ME-TIC and Varius. In this research the radiologist is used as reference to be able to compare the accuracy in classification between Dynamo, Auto-ME-TIC and Varius. The accuracy of the classification algorithm in Dynamo had been previously evaluated.
Table 3.3: The accuracy of the classification algorithm in Dynamo from previous research (Lavini et al.(2007)). Dynamo is the reference. Important to note is that for the results in the evaluation from there was no separate training and test set used for evaluation of TICs.
Type Matches(%) 2 78.6 3 78.6 4 90.0 5 67.1 6 82.8 7 28.0 All types 70.5
The overlap of matches (%) as shown in Table
3.3 and Table 3.4 was calculated by identi-fying the classification of the algorithm and comparing this with the human evaluator. In other words, how well does the human evalu-ator agree with the classification algorithm as reference? This direction of overlap was cal-culated in order to compare the results with results from previous research. The confu-sion matrices in Fig. 3.8, 3.9 and 3.10show the relative overlap of matches (%) and mis-matches (%) and is calculated by identifying the classification of radiologist and comparing this with the classification of the algorithm. In other words, how well does the classifica-tion algorithm agree with the radiologist as reference? The results from this research of the re-evaluation of Dynamo with the software as reference is presented in Table 3.4(see
Ap-pendix B.8and B.9for Auto-ME-TIC and Varius). An inspection of the classification of the radiologist reveals the classification of Type 2 has a moderate accuracy and Type 3 a high accuracy. The accuracy in classification of Type 4 and Type 6 is moderate. Type 5 and Type 7 are very poorly classified which indicates these are most difficult to classify.
Table 3.4: Percentage of matches (%) of the evaluation between Dynamo and all blind evaluators. Dynamo is the reference standard. A description of how well the evaluators agree with Dynamo. The interpretation of accuracy (% of matches) is defined as follows: 1-20% very poor, 21-40% poor, 41-70% moderate, 71-80% high, 81-90% very high, 91-100% almost perfect or perfect.
Type R.Hemke (%) C.Lavini (%) J.van den Noort (%) Mean (%) 2 (N = 242) 64.9 61.6 45.5 57.3 3 (N = 88) 70.5 43.2 60.0 58.0 4 (N = 89) 49.4 65.2 29.2 47.9 5 (N = 37) 18.9 35.1 37.8 30.6 6 (N = 73) 47.9 74.0 50.7 57.5 7 (N = 92) 12.0 16.3 29.3 19.2 All types (N = 627) 43.9 49.3 42.1 45.1
The results of the re-evaluation for Dynamo with the radiologist as the reference is shown in Table 3.5 and Fig. 3.8. A visual inspection of Fig. 3.8 reveals a moderate accuracy in classification for Type 2 and a very high accuracy for Type 6. Classification accuracy of Type 3 and Type 4 are moderate. Accuracy in classification of Type 5 is very poor and of Type 7 is poor which indicates that these are the most difficult to classify for Dynamo.
Table 3.5: Percentage of matches (%) of the evaluation between all evaluators and Dynamo. Each human evaluator was considered the reference standard. A description of how well Dynamo agrees with the evaluator. The interpretation of accuracy (% of matches) is defined as follows: 1-20% very poor, 21-40% poor, 41-70% moderate, 71-80% high, 81-90% very high, 91-100% almost perfect or perfect.
Type Number of TICs R.Hemke (%) Number of TICs C.Lavini (%) Number of TICs J.van den Noort (%) Mean (%) 2 231 68.0 199 74.9 139 79.1 74.0 3 149 41.6 73 52.1 170 31.2 41.6 4 104 42.3 125 46.4 77 33.8 40.8 5 53 13.2 97 13.4 108 13.0 13.2 6 40 87.5 69 78.3 41 90.2 85.3 7 50 22.0 64 23.4 92 29.3 24.9 Mean (%) (N = 627) 45.8 48.1 46.1 46.6
Figure 3.8: Confusion matrix of matches (%) and mismatches (%) between the radiologist and Dynamo. The x-axis and y-axis represent a Type. The color bar represents the overlap in the number of matches between the radiologist and Dynamo. The more green the color the more overlap was found. The more white the color the less overlap was found. The diagonal of the confusion matrix starts in the upper left corner and ends in the lower right corner. Ideally, the overlap in this diagonal is 100% and in all other cells 0%.
The results of the evaluation for Auto-ME-TIC with the radiologist as the reference is shown in Table 3.6and Fig. 3.9. A visual inspection of Fig. 3.9 reveals Type 6 is most accurate with an accuracy that is almost perfect. Accuracy in classification of Type 3 and Type 4 is moderate. Accuracy in classification of Type 2 and 5 is poor, and of Type 7 is very poor which indicates that these are most difficult to classify for Auto-ME-TIC.
Table 3.6: Percentage of matches (%) of the evaluation between all evaluators and Auto-ME-TIC. Each human evaluator was considered the reference standard. A description of how well Auto-ME-TIC agrees with the evaluator. The interpretation of accuracy (% of matches) is defined as follows: 1-20% very poor, 21-40% poor, 41-70% moderate, 71-80% high, 81-90% very high, 91-100% almost perfect or perfect.
Type Number of TICs R.Hemke (%) Number of TICs C.Lavini (%) Number of TICs J.van den Noort (%) Mean (%) 2 231 37.2 199 42.2 139 47.5 42.3 3 149 53.0 73 60.3 170 38.8 50.7 4 104 48.1 125 56.8 77 33.8 46.2 5 53 35.8 97 32.0 108 31.5 33.1 6 40 97.5 69 88.4 41 95.1 93.7 7 50 16.0 64 18.8 92 23.9 19.6 Mean (%) (N = 627) 47.9 49.7 45.1 47.6
Figure 3.9: Confusion matrix of matches (%) and mismatches (%) between the radiologist and Auto-ME-TIC. The x-axis and y-axis represent a Type. The color bar represents the overlap in the number of matches between the radiologist and Auto-ME-TIC. The more green the color the more overlap was found. The more white the color the less overlap was found. The diagonal of the confusion matrix starts in the upper left corner and ends in the lower right corner. Ideally, the overlap in this diagonal is 100% and in all other cells 0%.
