Stochastic Gradient Langevin dynamics for joint parameterization of tracer kinetic models, input functions, and T1 relaxation-times from undersampled k-space DCE-MRI

https://doi.org/10.1016/j.media.2020.101690Get rights and content

Highlights

Abstract

Dynamic Contrast Enhanced (DCE) Magnetic Resonance Imaging (MRI) is an important diagnostic technique that can quantify the structure and function of microvasculature processes, using T1 relaxation times and tracer kinetic maps. However, a series of methodological limitations affect both the accuracy and standardisation of the quantified maps, and consequently their diagnostic ability. The main methodological challenge in the quantification of tracer kinetics is a multi-parameter optimization, with correlated parameters that have different scales, which results in local minima particularly when measurements are highly undersampled. This work suggests a novel data driven optimization scheme, based on a variation of the Stochastic Gradient Langevin dynamics (SGLD) Markov chain Monte Carlo algorithm, which combines stochastic gradient descent and Langevin dynamics. The proposed SGDL algorithm avoided local minima and accurately quantified proton density, T1 relaxation times and tracer kinetics. Joint direct parameterization significantly benefited the quantification of proton density, T1 relaxation times, and the selection of a suitable tracer kinetic model per tissue type. Model based arterial and portal vein input functions were automatically determined during the joint direct parameterization. Observations made on simulated fully and highly undersampled liver DCE MRI data were confirmed on acquired clinical data.

Introduction

Dynamic Contrast Enhanced (DCE) Magnetic Resonance Imaging (MRI) can quantitatively assess the characteristics of the blood vessels using mathematical models that describe how the tracer perfuses through the vessels (Dikaios et al., 2017). Malignant tumours depend upon neovascularization or angiogenesis to proliferate and survive. Tumour vessels produced by angiogenesis are structurally abnormal and hyper-permeable compared to normal vessels. Accurate, repeatable and reproducible quantification of tracer kinetics can help clinical oncologists detect and characterize tumours and monitor their response to treatment. However, there are limitations that affect the quantification of tracer kinetics (Schabel and Parker 2008; Subashi et al., 2014): (i) MRI measurements depend on field inhomogeneities, 3D gradient non-linearity, chemical shifts, etc., causing image distortions (Heye et al., 2013), (ii) spatiotemporal resolution (Yuan et al., 2012), (iii) quantification of the native T1 of the tissue before bolus injection (T10) (Heye et al., 2014), and (iv) estimation of the tracer's input functions (Woolf et al., 2016).

High temporal resolution is necessary to capture rapid microvascular exchanges, and higher spatial resolution can secure more coverage. Parallel imaging and/or Compress Sensing can improve the temporal and/or spatial resolution but result in images with lower SNR especially for high undersampling rates. Errors in T10 quantification also affect the estimation of tracer kinetics. To derive patient specific T1 maps, variable flip angle methods based on spoiled gradient echo sequences (SPGR) are commonly preferred because they are fast and reliable. However, the reconstructed images suffer from low SNR and the T1 mapping is not robust in regions with high T1 and low steady state magnetization (Sandmann et al., 2016). Patient-specific arterial input function (AIF) estimated from manual ROIs drawn by radiologists are commonly used in clinical studies, but they are affected by partial volume effects and by corrupted MR signal in the arteries due to flow artifacts (van Osch et al., 2005). Alternatively, population-based AIFs are utilized, but they ignore differences in injection rate and cardiac output that vary per patient and result in non-reproducible tracer kinetic quantification (Galbraith et al., 2002). Automated AIF identification from large vessels is a recent development for patient specific AIF and has been shown to be more reliable (Chen et al., 2008).

Direct reconstruction of T1 maps (Velikina et al., 2013) or tracer kinetics (Dikaios et al., 2014a, Dikaios et al., 2014b; Felsted et al., 2009; Guo et al., 2017) from highly undersampled DCE MRI k-space data were suggested to avoid optimization problems due to image distortions and low SNR. Such methods have shown promising results, however because they operate in the k-space instead of the enhanced images, certain key issues need to be resolved, i.e. (i) estimate a patient specific IFs during the optimization, without the need for a radiologist to draw ROIs in the image domain and (ii) select a tracer kinetic model per tissue type. Joint direct parameterization of tracer kinetics and arterial input functions has been suggested in the literature (Guo et al., 2018), but it is currently only applicable to organs such as the brain which can be described by the Toft model (and its extension) and requires the identification of an arterial ROI in the contrast enhanced images. To apply joint direct methods to other organs such as the liver, more tracer kinetic models need to be accounted for, which adds extra degrees of complexity making the optimization more difficult. Specifically, to robustly quantify the tracer kinetics the optimization needs to: (i) handle highly undersampled k-space data, (ii) be independent of the parameters initialization or any other trade off parameters in the likelihood function, (iii) avoid local minima when the number of recovered parameters is high, and (iv) avoid bias in areas with different tracer kinetics due to high correlation of the kinetic basis functions.

This work focusses on the optimization problem suggesting a novel model-based algorithm that jointly parameterizes tracer kinetics, proton density, T10, input functions and selects a suitable tracer kinetic model per tissue type. The previously reported Bayesian Inference Markov Chain Monte Carlo (MCMC) algorithm (Dikaios et al., 2014a, Dikaios et al., 2014b) parameterized 3 different tracer kinetics related to the extended Toft model from highly udersampled k-space, but results indicated bias in areas with different tracer kinetics due to high correlation of the kinetic basis functions. The direct parameterization of more parameters with different scales that are correlated, will require a different approach to ensure convergence. Specifically, a variation of the Stochastic Gradient Langevin dynamics (SGLD) algorithm (Welling and Teh, 2011, Patterson and Teh, 2013) is suggested, which accelerates the mixing of Langevin dynamics, while it still guarantees convergence. The main contributions of this work are:

  • (1)

    A novel stochastic optimization algorithm (SGLD) that can jointly parameterize proton density, T10 and tracer kinetics

  • (2)

    Automatically select a tracer kinetic model per tissue type and update during the optimization

  • (3)

    Mathematically model arterial and portal vein input functions and update the model parameters during the optimization.

The suggested joint direct parameterization algorithm is evaluated on a simulated respiratory liver DCE digital phantom, where contrast enhancement was modeled using different tracer kinetic models and arterial/portal vein input functions. The evaluation was performed for fully sampled and highly undersampled (k, t)-space data. Radially undersampled liver DCE MRI data were also used to evaluate the proposed SGLD algorithm on real data.

Section snippets

Mathematical model of the DCE MRI signal

Following the intravenous injection of a paramagnetic contrast agent its concentration in the extracellular extravascular space (EES) is usually described with the extended Tofts model (Tofts 1997) using the tracer kinetic parametersC(r,t;ψ,φ)=vp(r)·Ca(t)+Ktrans(r)·0tCa(r,τt0)e(Ktrans(r)ve(r)·(tτ))dτ,or the dual input function Orton model (Orton et al., 2009) for the liver, where a combination of an arterial and a portal vein input function is used.C(r,t;ψ,φ)=Ktrans(r)·0tCp(τt0)e(Ktrans

Simulated liver DCE MRI (k, t)-space data

T1-weighted abdominal images of a normal volunteer were acquired in multiple time frames without contrast injection using a fast gradient echo DCE-MRI protocol (flip angle 10°, repetition time 2.3 ms). The first time-frame was manually segmented into: liver, bowel, right and left heart, aorta, portal vein to simulate contrast enhancement using either the arterial input function (flow, (modified) Toft) or dual input functions (for the liver). Each organ was assigned different tracer kinetic

Results

The proposed SGDL algorithm is validated on simulated and acquired liver DCE MRI data. For the simulated dataset, results are shown for the 20-fold spiral undersampled (k, t)-space DCE MRI.

Discussion

Mathematical optimization determines the “best” possible solution to problems expressed in terms of parameters (also referred to as degrees of freedom). It has applications in scientific areas, where models are used to describe measured data. SGLD is an optimization method which combines stochastic gradient descent and Langevin dynamics. This work suggests a variation of the SGDL algorithm, to estimate vasculature (tracer kinetics) and structural (proton density, T1 relaxation times) parameters

CRediT authorship contribution statement

Nikolaos Dikaios: Conceptualization, Methodology, Software, Writing - original draft.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgement

This work has been supported by Royal Society (INF\R1\191030).

References (37)

  • L. Feng et al.

    Golden-angle radial sparse parallel MRI: combination of compressed sensing, parallel imaging, and golden-angle radial sampling for fast and flexible dynamic volumetric MRI

    Magn. Reson. Med.

    (2014)
  • D. Flouri et al.

    Fitting the two-compartment model in DCE-MRI by linear inversion

    Magn. Reson. Med.

    (2016)
  • S.M. Galbraith et al.

    Reproducibility of dynamic con-trast-enhanced MRI in human muscle and tumours: comparison of quantitative and semi-quantitative analysis

    NMR Biomed.

    (2002)
  • Y. Guo et al.

    Direct estimation of tracer-kinetic parameter maps from highly undersampled brain dynamic contrast enhanced MRI

    Magn. Reson. Med.

    (2017)
  • Y. Guo et al.

    Joint arterial input function and tracer kinetic parameter estimation from undersampled dynamic contrast-enhanced MRI using a model consistency constraint

    Magn. Reson. Med.

    (2018)
  • T. Heye et al.

    Reproducibility of dynamic contrast-enhanced MR imaging. Part II. Comparison of intra- and interobserver variability with manual region of interest placement versus semiautomatic lesion segmentation and histogram analysis

    Radiology

    (2013)
  • T. Heye et al.

    Impact of precontrast T10 relaxation times on dynamic contrast-enhanced MRI pharmacokinetic parameters: T10 mapping versus a fixed T10 reference value

    J. Magn. Reson. Imaging

    (2014)
  • R.B. Ger et al.

    A multi-institutional comparison of dynamic contrast-enhanced magnetic resonance imaging parameter calculations

    Sci Rep.

    (2017)
  • Cited by (10)

    • Physics-informed neural networks for myocardial perfusion MRI quantification

      2022, Medical Image Analysis
      Citation Excerpt :

      While, in this work, different loss functions were tested and different inputs to the PINN were tried, there is also a huge potential to exploit the flexible framework further. For example, it would be possible to change the input from the image-derived concentration curves to the acquired k-t space and to learn both the image reconstruction and kinetic parameter estimation, similar to Dikaios (2020). However, this is also a potential limitation as there are a wide range of design choices and hyperparameters to be set, which have the potential to influence the results obtained.

    • Quantification of T1, T2 relaxation times from Magnetic Resonance Fingerprinting radially undersampled data using analytical transformations

      2021, Magnetic Resonance Imaging
      Citation Excerpt :

      Sophisticated reconstruction and fitting algorithms have been suggested in the literature to robustly perform the fitting in the reconstructed images. Model-based reconstructions [1–3] have been suggested that aim to recover magnet properties of tissue directly from k-space data. Such methods can account for noise and highly undersampled k-space data (to improve the temporal resolution).

    • Modeling dynamic radial contrast enhanced MRI with linear time invariant systems for motion correction in quantitative assessment of kidney function

      2021, Medical Image Analysis
      Citation Excerpt :

      Note that the algorithm used to fit the proposed LTI model based registration is essential to be robust to motion and streaking artifacts. Recent work proposed using TK models instead of an LTI model to perform joint parameter estimation and registration Dikaios (2020) using a data driven optimization scheme. Because a single TK model does not fit to all tissues types, they used a model selection step to determine which model is best suited for each different tissue type.

    View all citing articles on Scopus
    View full text