Fast multidimensional NMR spectroscopy for sparse spectra

Fast multidimensional spectroscopy for sparse spectra Dany MERHEJ1,2,3, Hélène RATINEY1, Chaouki DIAB3, Mohamad KHALIL2, and Rémy PROST1

1

Université de Lyon, CNRS, Inserm, INSA-Lyon, CREATIS, UMR5220, U630, Villeurbanne, France

2

EDST, Azm research center, Lebanese University, Tripoli, Lebanon

3

ISAE – Cnam Liban, Beirut, Lebanon

Address correspondence to: Rémy PROST INSA-Lyon CREATIS CNRS UMR 5220 Address: Building Blaise Pascal, 7 Avenue Jean Capelle, 69621 Villeurbanne cedex, FRANCE. Tel:

(33) 4 72 43 80 72 (33) 4 72 43 82 27

Fax:

(33) 4 72 43 85 26

E-MAIL: [email protected]

This work was supported in part by a scholarship from the Agence Universitaire de la Francophonie (AUF), and the Région Rhône Alpes Cluster 2 ISLE, PP3, subproject SIMED: Simulation en Imagerie MEDicale pour le Diagnostic et la Thérapie. This work is in the scope of the scientific topics of the PRC-GDR ISIS research group of the French National Center for Scientific Research (CNRS).

Approximate word count:

194 (Abstract)

4429 (body)

Submitted to NMR in biomedicine as a Full Paper.

1

Abstract Multidimensional NMR spectroscopy, including magnetic resonance spectroscopic imaging, is widely used for studies of molecular and biomolecular structure. A major disadvantage of multidimensional NMR is the long acquisition time that, regardless of sensitivity considerations, may be needed for obtaining the final multidimensional frequency domain coefficients. In this paper, a method for undersampling multidimensional NMR of sparse spectra is presented. The approach uses a priori knowledge about the possible nonzero locations in the multidimensional frequency domain to remove the resulting under-determinacy when undersampled acquisitions are performed. In order to obtain this a priori knowledge at a low cost, we take a one dimensional NMR acquisition in the case of two dimensional NMR - which is usually always acquired-, and a single voxel acquisition in the case of multi-voxel spectroscopic imaging. If these low dimension acquisitions are intrinsically sparse, it is then possible to reconstruct the corresponding multidimensional acquisitions from far fewer observations than that imposed by the Nyquist criterion and subsequently reduce acquisition time. Simulated and experimental acquisitions demonstrate reliable and fast reconstructions with acceleration factor of order 3. The approach can be easily extended to higher dimensions. Key Words – Multidimensional NMR spectroscopy, undersampling, sparse spectra.

2

List of abbreviations Card

Cardinality

COSY

Correlation spectroscopy

CS

Compressed sensing

CSI

Chemical shift imaging

DFT

Discrete Fourier transform

FID

Free induction decay

FOV

Field of view

IS

Identified support

MR

Magnetic resonance

MRI

Magnetic resonance imaging

MRS

Magnetic resonance spectroscopy

MRSI

Magnetic resonance spectroscopic imaging

NMR

Nuclear magnetic resonance

SLS

Spatially localized spectra

SNR

Signal to noise ratio

SNRrec

Reconstruction SNR

TE

Echo time

TOCSY

Total correlation spectroscopy

TR

Repetition time

3

Material for the Graphical Abstract

Fast multidimensional spectroscopy for sparse spectra Dany MERHEJ, Hélène RATINEY, Chaouki DIAB, Mohamad KHALIL, and Rémy PROST* -0.2

-0.2 0

0

0.2

0.2

0.4

0.4

kF11 (ppm)

kF11 (ppm)

A method for undersampling multidimensional NMR of sparse spectra is presented. It uses a priori knowledge about the possible nonzero locations in the multidimensional frequency domain to remove the resulting under-determinacy from undersampled acquisitions. In order to obtain this a priori knowledge at a low cost, we take a one dimensional NMR acquisition in two dimensional NMR, and a single voxel acquisition in multi-voxel spectroscopic imaging. The approach can be easily extended to higher dimensions.

0.6 0.8

0.6 0.8

1

1

1.2

1.2 1.4

1.4 1.3

1.2

1.1

1

0.9

0.8

0.7

kF22 (ppm)

HH-COSY of octane from fully sampled measurments

1.3

1.2

1.1

1

0.9

0.8

0.7

kF22 (ppm)

HH-COSY of octane from 33% of the measurements

4

Introduction Nuclear magnetic resonance (NMR) spectroscopy is a rich source of molecular information. Almost every structure determination of any organic or biological molecule, as well as that of many inorganic molecules, often begins with nuclear magnetic resonance spectroscopy (MRS). Aside from X-ray crystallography, MRS is the most direct and general tool for identifying the structure of both pure compounds and mixtures as either solids or liquids. The process often involves performing several NMR experiments to deduce the molecular structure from the magnetic properties of the atomic nuclei and the surrounding electrons. The introduction of two-dimensional (2D) NMR spectroscopy during the early 1980s (1), revolutionized traditional one dimensional (1D) NMR spectroscopy and made it easier to study more complex molecular systems, such as proteins, nucleic acids, etc. MRS has been combined with magnetic resonance imaging (MRI) to provide spatially localized image of biochemistry in living tissue by using signals from organic molecules; and thus provides a powerful tool for quantitative non-invasive, in-vivo-detection, and imaging of biomarkers – also known as metabolites, for the evaluation of response to therapy, and the biochemical characterization of pathology, in the diseased brain, heart, prostate, and other parts of the human and animal body. MRS localization techniques include, monovoxel acquisition which detects signal from a single voxel, and magnetic resonance spectroscopic imaging (MRSI) also called Chemical Shift Imaging (CSI), allowing detection and localization of spectra from several spatially distributed voxels. High-dimensional NMR spectroscopy is however very time consuming. In conventional multi-dimensional NMR acquisition, the 1D-NMR is repeated as many times as required for the desired resolution in each additional dimension, resulting in very long acquisition time. These experiments often run for many hours, days and weeks (2), consuming spectrometer time and limiting the possible throughput. The same consideration can be made for multi-voxel MRSI where a major drawback is the long acquisition time required to gather necessary data to achieve satisfactory resolution at acceptable signal to noise ratio (SNR). In a conventional MRSI experiment, the k-space (Fourier domain) is sampled by incrementing the amplitude of a magnetic field gradient 5

pulse in a similar manner as phase encoding in MRI. At each application of the phase encoding gradient, a k-space location is selected and the spectroscopic signal encoding the chemical shift dispersion is acquired. This conventional acquisition technique results in long acquisition times, e.g. more than 10 min on the human brain and 40 min on the mouse or rat brain, and low spatial resolution, e.g. bigger than 1cc for the human brain, a few µL for the small animal brain, because of low SNR and low spatial encoding scheme. Therefore, reducing scan time is a primary goal in the development of multidimensional NMR technique. A scan time reduction, while preserving frequency resolution in the indirect dimension, would enable acquisitions in higher dimensions in multidimensional NMR, and better spectrometer throughput. A scan time reduction, while preserving sensitivity per units of time and volume in MRSI, would enable an improvement of patient comfort resulting from shorter exam duration. To this end, reducing scan time in multidimensional NMR has been the focus of many investigations. In 2D NMR, faster acquisition schemes have been proposed for example by limiting the number of indirect dimension increments and employing adequate subsequent reconstruction algorithm (Linear Prediction (3) or maximum entropy (4)) or irregularly sampling the indirect dimension such as in Hadamard spectroscopy (5), (6), or else by invoking radial sampling (7). In spectroscopic imaging, novel fast encoding methods are mainly based on high speed MRI such as Echo Planar Imaging, RARE (8), spiral encoding (9) and steady state free precession (10). In these cases, the reconstruction, regularisation, and data formatting stem from the designed acquisition schemes. In this paper, we propose a new approach based on modelling the whole process from standard nD data acquisition to spatially (CSI) or spectrally (2D/3D) resolved spectroscopy, in a matrix-vector notation. We will then use a priori knowledge, on the nonzero locations in the multidimensional spectra, to remove the under-determinacy

induced

from

data

undersampling.

This

will

allow

faster

multidimensional NMR acquisitions (11), including MRSI (12). As the proposed approach relies on matrix-vector notation, the following section will begin by modelling the acquisition and reconstruction procedures in that form.

6

Theory Two dimensional NMR data acquisition and matrix notation No matter what the nature of the interaction to be mapped, all 2D acquisition sequences have the same basic format, and can be subdivided into four well-defined units termed the preparation, evolution, mixing and detection periods as shown in fig. 1. In standard 2DNMR, the t2 time axis, known as direct dimension, is acquired directly and is entirely analogous to the detection period of any 1D-NMR experiment during which the spectrometer collects the free induction decay (FID) of the excited spins. The additional time axis, known as indirect dimension, is monitored as a discrete incrementation of the evolution time variable t1, throughout a series of repetitions (scans) of the basic experiment. Each scan will fill exactly one single line along the t2 axis in a 2D grid known as the 2D time domain. Suppose Nt scans are performed and Nt samples are acquired in each detection period, our time domain is a 2D Nt×Nt grid acquired line by line along the discrete time axis t2. Let

x  t1 , t2   , 0  t1 , t2  Nt  1 , represent the total acquired data set of N t2 samples in the time domain, the corresponding 2D-NMR spectra, i.e. the frequency domain, is computed by applying a 2D unitary discrete Fourier transform (DFT) as follows

X  k1 , k2  

1 Nt

Nt 1 Nt 1



11



t



2j (k t  k t )    x t , t  .exp  N t1  0 t2  0

1

2

2 2



[1]

with k1 , k2   , 0  k1 , k2  Nt  1 , and j  1 . In order to present our proposed approach in its simplest form, we would like to have a matrix-vector notation relating the unknown frequency coefficients and time domain samples. Let xs, Xs  n , with n  Nt2 , represent the stacked column ordered vectors of the time domain samples, and unknown frequency domain coefficients, respectively (see fig. 2). Let F1, F2, represent the one dimensional unitary Fourier transform matrices, of size Nt×Nt, related to the indirect and direct dimensions respectively, and let F be the inverse of their Kronecker product of size n×n 7

F   F2  F1 

1

[2]

where  denotes the kronecker product. The mapping from the unknown 2D-NMR frequency domain to the 2D time domain measurements can be expressed as xs  FX s

[3]

Fig. 2 gives a synoptic of this mapping. Equation [3] is of the form:

y  Fx

[4]

where the vector yCn represents our acquired measurements x s , and xCn the unknown 2D frequency domain coefficients X s .

MRSI data acquisition and matrix notation MRSI combines both spectroscopy and imaging methods to produce spatially localized spectra (SLS) from within a sample or a patient. The MRSI spatial resolution is however much lower than that of an MR image, limited by the required acquisition SNR. Standard 2D-MRSI data acquisition and processing is performed as follows. Data are acquired, on a 2D kx×ky grid known as the k-space, i.e. the discrete 2D Fourier domain, along a third discrete time axis t, t  . The procedure is repeated for each (kx, ky) location in the k-space grid. Suppose for example a 2D grid of size Nk×Nk, and  samples acquired along the t time axis in each (kx, ky) location, then the total acquired data set is N k2 × samples (Fig. 3). MRSI computation goes then through an intermediate image domain, calculated individually for each acquired k-space. Let X  k x , k y , t   represent a unique 2D grid along time axis t with



Nk

2

 kx , k y 

Nk

2

 1 and

0  t    1 , then the corresponding image,

x  k x , k y , t  , is computed by using the 2D unitary inverse DFT

1 x  nx , ny , t   Nf

Nk 1 2

Nk 1 2

 

 Nk  Nk ky  kx  2 2

 2j (nx k x  n y k y )  X  k x , k y , t  .exp   Nk  

[5]

8

with nx , ny   ,

 Nk

2

 nx , ny 

Nk

2

 1 , and j  1 as previously defined. Final SLS

are obtained by performing a 1D unitary DFT at each (nx,ny) location in the image plane, along the time axis t, as follows X SLS  nx , ny , f  

 2jft   

 1

 x  n , n , t  .exp  

1

x

t 0

y

[6]

with f   and 0  f    1 . Here  represents both, the number of time domain samples in the t direction, and the number of chemical shift coefficients in each SLS. The calculated chemical shift coefficients are located at each

f

 frequency interval. An optional

apodization window w(t) can also be used. Fig. 3 gives a synoptic of a conventional 2D– MRSI processing from acquired raw data to SLS. Also here, we would like to have a matrix-vector notation relating the unknown SLS coefficients and K-space samples. Let Xs(t)  Nk , with 0  t    1 , represent the stacked 2

column ordered vectors of the 2D k-space domain samples at index t, and let Xss n , with n  Nk2  , represent the stacking of all Xs(t) vectors. Let XSLSs n represent the stacking

of all corresponding unknown SLS coefficients and let F1,F2 and Ft, represent the 1D unitary DFT matrices, of appropriate sizes related to dimensions kx, ky and t respectively. The mapping from the unknown SLS coefficients to the k-space measurements can then be expressed as: X ss  FX SLSs

[7]

F=(Ft)-1F2F1

[8]

with and F nn . Fig. 4 gives a synoptic of this mapping, with an example for Nk=3. Equation [7] is of the form of equ. [4], i.e. y=Fx, where here the vector yCn with n  Nk2  , represents our acquired k-space measurements, and x n the unknown SLS coefficients.

9

A priori information and system undersampling We aim to undersample the acquisition domain, in order to perform faster experiments in MRSI and multidimensional NMR, with acceptable reconstruction SNR. However undersampling would result in under-determinacy in the previously proposed acquisition model, i.e. the y=Fx model, and subsequently the impossibility to identify the original signal, using linear reconstruction methods. A newly introduced signal processing scheme, known as Compressed sensing (CS)(13-18), has succeeded in making a breakthrough in this problem, when the signal x is sparse or compressible, using nonlinear L1-minimization reconstructions (18,19). The vector x does not have to be intrinsically sparse, but could be sparse in any other known transform domain, incoherent with the acquisition domain (20). Our proposal relies however on the intrinsic sparsity of the acquired spectra. More precisely, in the case of multidimensional NMR and MRSI, if the performed lower dimensional acquisitions (as will be described in the two following sections) are intrinsically sparse, we will be able to efficiently identify a set T{1,2,..n}, of indices of possible nonzero coefficients in the final frequency domain x. We will refer to this set T by the identified support (IS). Eventually we should have T0  T , and Card(T0) ≤ Card(T), where T0 corresponds to the true and unknown set of indices, of nonzero coefficients of x. The more our one dimensional acquisitions are sparse, the more we will be able to reduce the set T such that Card(T)