Elsevier

Medical Image Analysis

Volume 49, October 2018, Pages 117-127
Medical Image Analysis

Unmixing dynamic PET images with variable specific binding kinetics

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

Highlights

  • An enhanced factor analysis of dynamic PET images is proposed, which explicitly models the variability of the specific binding factor time-activity-curve.

  • The decomposition relies on a new interpretation for the spatial heterogeneity of PET images.

  • Taking into account the spatial heterogeneity clearly improves the estimation of the specific binding factor and its spatial map.

Abstract

To analyze dynamic positron emission tomography (PET) images, various generic multivariate data analysis techniques have been considered in the literature, such as principal component analysis (PCA), independent component analysis (ICA), factor analysis and nonnegative matrix factorization (NMF). Nevertheless, these conventional approaches neglect any possible nonlinear variations in the time activity curves describing the kinetic behavior of tissues with specific binding, which limits their ability to recover a reliable, understandable and interpretable description of the data. This paper proposes an alternative analysis paradigm that accounts for spatial fluctuations in the exchange rate of the tracer between a free compartment and a specifically bound ligand compartment. The method relies on the concept of linear unmixing, usually applied on the hyperspectral domain, which combines NMF with a sum-to-one constraint that ensures an exhaustive description of the mixtures. The spatial variability of the signature corresponding to the specific binding tissue is explicitly modeled through a perturbed component. The performance of the method is assessed on both synthetic and real data and is shown to compete favorably when compared to other conventional analysis methods.

Introduction

Dynamic positron emission tomography (PET) is a non-invasive nuclear imaging technique that allows biological processes to be quantified and organ metabolic functions to be evaluated through the three-dimensional measure of the radiotracer concentration over time.

The analysis of dynamic PET images, in particular the quantification of the kinetic properties of the tracer, requires the extraction of tissue time-activity-curves (TACs) in order to estimate the parameters from compartmental modeling (Innis , 2007). Nevertheless, PET images are corrupted by a prominent statistical noise and elementary TACs are mixed up due to the partial volume effect. Therefore, inferring rigorous and reliable information from these images still remains a challenging issue.

Several generic methods have been applied to estimate elementary TACs and their corresponding proportions from dynamic PET images. These techniques have different denominations depending on the application context, but all aim at tackling blind source separation (BSS) problems. For instance, Barber (1980) and Cavailloles (1984) proposed matrix factorization-based PET analysis techniques, referred to as factor analysis of dynamic structures (FADS) (Chou , 2007). An improvement of FADS taking into account the nonnegative physiological characteristic of PET images was proposed subsequently by Wu (1995) and Sitek (2000), also appearing in other biomedical imaging domains (Martel , 2001). Nonnegative matrix factorization (NMF) techniques pursue the same objective, under nonnegativity constraints on the latent factors to be recovered, and have been intensively applied in dynamic PET studies (Lee, Lee, Choi, Park, Lee, 2001, PadillaP.andLpez, Grriz, Ramrez, Salas-Gonzlez, lvarez, 2012, Schulz, 2013). The works of Sitek (2002) and El Fakhri (2005) improved nonnegative FADS with a penalization that promoted non-overlapping regions in each voxel.

The approach proposed in this paper follows the same line as NMF or nonnegative FADS. It aims at decomposing each PET voxel TAC into a weighted combination of pure physiological factors, representing the elementary TAC associated with the different tissues present within the voxel. This factor modeling is enriched with a sum-to-one constraint to the factor proportions, so that they can be interpreted as tissue percentages within each voxel. In particular, this additional constraint explicitly solves the scaling ambiguity inherent to any NMF models, which has proven to increase robustness as well as interpretability. This BSS technique, referred to as unmixing or spectral mixture analysis, originates from the geoscience and remote sensing literature (Bioucas-Dias , 2012) and has proven its interest in other applicative contexts, such as microscopy (Huang , 2011) and genetics (Dobigeon and Brun, 2012).

The linearity assumption for voxel decomposition underlined by the above-mentioned PET analysis methods can be envisioned in the light of compartment modeling, a tool widely employed to describe the kinetic behavior occurring within the voxel. The selection of a given model should be based on the radiotracer under study (Gunn, Lammertsma, Hume, Cunningham, 1997, Innis, Cunningham, Delforge, Fujita, Gjedde, Gunn, Holden, Houle, Huang, Ichise, et al., 2007) but most of the models consider that the measured signal in a given voxel is the sum of the comprising compartments.

However, factor TACs to be recovered cannot always be assumed to have constant kinetic patterns, as implicitly considered in conventional methods. Fig. 1 depicts an example of a 2-tissue compartmental model, where the radioligand is assumed to move between three compartments: Cp represents the radioligand concentration in arterial plasma, CF+NS represents the free plus non-specific compartment, and CS represents the specifically bound compartment. The exchange between compartments are subject to rate constants kj (j=1,,4) (Häggström , 2016). Considering both the 2-tissue and reference compartment models, the assumption of constant kinetic patterns seems appropriate for the blood compartment as well as non-specific binding tissues, since they present some homogeneity besides some perfusion difference (e.g. white matter versus gray matter). Therefore their contribution to the voxel TAC should be fairly proportional to the fraction of this type of tissue in the voxel. However, things get different regarding the specific binding class, as the TAC associated with this tissue is nonlinearly dependent on both the perfusion and the concentration of the radiotracer target. The spatial variation on target concentration is in part governed by differences in the k3 and k4 kinetic parameters, which nonlinearly modify the shape of the TAC characterizing this particular class. Muzi (2005) discussed the accuracy of parameter estimates for tumor regions and underlined high errors for the parameters related to specific binding, namely 26% for k3 and 49% for k4. These results were further confirmed by Schiepers (2007). More specifically, they studied the kinetics of lesioned regions that were tumor and treatment change predominant, showing that variations on k3 and k4 may allow for differentiation. Bai (2013) further discussed nonuniformity in intratumoral uptake and its impact on predicting treatment response and tumor aggressiveness. Nonetheless, this fluctuation phenomenon has not been taken into account by the decomposition models from the literature.

The main motivation of this paper is to propose a more accurate description of the tissues and kinetics composing the voxels in dynamic PET images, in particular for those affected by specific binding. To this end, this work proposes to explicitly model the nonlinear variability inherent to the TAC corresponding to specific binding, by allowing the corresponding factor to vary spatially. This variation is approximated by a linear expansion over the atoms of a dictionary, which have been learned beforehand by conducting a principal component analysis on a learning dataset.

The sequel of this paper is organized as follows. The proposed mixing-based analysis model is described in Section 2. Section 3 presents the corresponding unmixing algorithm able to recover the factors, their corresponding proportions in each voxel and the variability maps. Section 4 details synthetic data generation and real data acquisition. Simulation results obtained with synthetic data and experimental results on real data are reported in Section 5. Section 7 concludes the paper.

Section snippets

Specific binding linear mixing model (SLMM)

Consider N voxels of a 3D dynamic PET image acquired at L successive time-frames. First, we omit the spatial blurring induced by the point spread function (PSF) of the instrument and any measurement noise. The TAC in the nth voxel (n{1,,N}) over the L time-frames is denoted xn=[x1,n,,xL,n]T. Akin to various BSS techniques and following the linear mixing model (LMM) for instance advocated in the PET literature by Barber (1980), each TAC xn is assumed to be a linear combination of K elementary

Algorithm implementation

Given the nature of the optimization problem (9), which is genuinely nonconvex and nonsmooth, the adopted minimization strategy relies on the proximal alternating linearized minimization (PALM) scheme (Bolte , 2013). PALM is an iterative, gradient-based algorithm which generalizes the Gauss-Seidel method. It consists in iterative proximal gradient steps with respect to A, M and B and ensures convergence to a local critical point A*, M* and B*. The principle of PALM is briefly recalled by

Synthetic data generation

To illustrate the accuracy of our algorithm, experiments are first conducted on synthetic data for which the ground truth of the main parameters of interest (i.e., factor TACs and factor proportion maps) is known. In the clinical PET framework, ground truth concerning the tracer kinetics and uptake is never completely known. Meanwhile, simulations benefit from an entire knowledge of the patient properties and kinetics, and their degree of complexity and details can be selected according to the

Evaluation on synthetic data

The factor proportion maps recovered by the compared algorithms are shown in Fig. 3. Each column corresponds to a specific factor: SBF, white matter, non-specific gray matter, blood (from left to right, respectively). The six rows contain the factor proportion maps of the ground truth, and those estimated by K-means, NMF, VCA, LMM and the proposed SLMM (from top to bottom, respectively). A visual comparison suggests that the factor proportion maps obtained with LMM and SLMM are more consistent

Performance of the method

This work proposes a novel unmixing method, called SLMM, that takes into account the spatial variability of high-uptake tissues by modeling the SBF with an additional degree of freedom at each pixel. It also introduces a simplified version of this model with no variability, that consists of a regularized and constrained unmixing algorithm, herein named LMM.

For the cases studied in this paper, SLMM and LMM always provide physically interpretable results, though not always the best, on the

Conclusion and future works

This paper introduced a new model to conduct factor analysis of dynamic PET images. It relied on the unmixing concept accounting for specific binding TACs variation. The method was based on the hypothesis that the variations within the SBF can be described by a small number of basis elements and their corresponding proportions per voxel. The resulting optimization problem is extremely non-convex with highly correlated factors and variability basis elements, which leads to a high number of

References (47)

  • J. Bolte et al.

    Proximal alternating linearized minimization for nonconvex and nonsmooth problems

    Math. Program.

    (2013)
  • E. Bullmore et al.

    Colored noise and computational inference in neurophysiological (fMRI) time series analysis: resampling methods in time and wavelet domains

    Hum. Brain Mapp.

    (2001)
  • F. Cavailloles et al.

    Factor analysis in gated cardiac studies

    J. Nuclear Med.

    (1984)
  • Y.C. Cavalcanti et al.

    Unmixing dynamic PET images with variable specific binding kinetics – Supporting materials

    Technical Report

    (2017)
  • L. Condat

    Fast projection onto the simplex and the l1-ball

    Math. Program.

    (2015)
  • N. Dobigeon et al.

    Spectral mixture analysis of EELS spectrum-images

    Ultramicroscopy

    (2012)
  • J.A. Fessler

    Penalized weighted least-squares image reconstruction for positron emission tomography

    IEEE Trans. Med. Imag.

    (1994)
  • R.N. Gunn et al.

    Parametric imaging of ligand-receptor binding in PET using a simplified reference region model

    Neuroimage

    (1997)
  • A. Halimi et al.

    Unsupervised unmixing of hyperspectral images accounting for endmember variability

    IEEE Trans. Image Process.

    (2015)
  • S. Henrot et al.

    Does deblurring improve geometrical hyperspectral unmixing?

    IEEE Trans. Image Process.

    (2014)
  • Y. Huang et al.

    Temporal dynamics of host molecular responses differentiate symptomatic and asymptomatic influenza a infection

    PLoS Genet.

    (2011)
  • R.B. Innis et al.

    Consensus nomenclature for in vivo imaging of reversibly binding radioligands

    J. Cereb. Blood Flow Metab

    (2007)
  • N. Keshava

    A survey of spectral unmixing algorithms

    Lincoln Lab. J.

    (2003)
  • Cited by (8)

    • Compartment model-based nonlinear unmixing for kinetic analysis of dynamic PET images

      2023, Medical Image Analysis
      Citation Excerpt :

      However, standard factor analysis techniques rely on the assumption that the elementary response of each tissue to tracer distribution is spatially homogeneous. To overcome this limitation, Cavalcanti et al. (2018) introduced a factor analysis model that handled fluctuations in specific binding (SB) kinetics with a spatially indexed variability. The proposed model tried to extract a factor for the blood input function, a factor for each non-specific binding (nSB) tissue of the region under study and assigned a spatially varying factor for high-uptake tissues.

    • Matrix cofactorization for joint representation learning and supervised classification – Application to hyperspectral image analysis

      2020, Neurocomputing
      Citation Excerpt :

      More precisely, the observed measurements can be approximated by mixtures of dictionary elements able to simultaneously capture the variability and redundancy in the dataset. Representation learning can be tackled from different perspectives, in particular known as dictionary learning [13], source separation [14], compressive sensing [15], factor analysis [16], matrix factorization [17] or subspace learning [18]. Various models have been proposed to learn a dedicated representation relevant to the field of interest, differing by specific assumptions and/or constraints.

    • Unmixing Dynamic Pet Images: Combining Spatial Heterogeneity and Non-gaussian Noise

      2019, ICASSP, IEEE International Conference on Acoustics, Speech and Signal Processing - Proceedings
    View all citing articles on Scopus

    Part of this work has been supported by CAPES.

    View full text