Nuclear magnetic resonance spectroscopy

Abstract

A method is provided for performing NMR spectroscopy. The method comprises positioning a sample in a homogeneous stationary magnetic field directed along an axis, preparing nuclei in at least a predetermined volume of the sample for resonant emission of an NMR signal and creating this NMR signal. This comprises irradiating the sample with at least one RF excitation pulse in accordance with an MRI sequence preparation and/or evolution module. The method also comprises sensing the NMR signal in the absence of frequency encoding magnetic field gradients such that analysis of the NMR signal yields a chemical shift spectrum from the nuclei. During this sensing, a plurality of intermittently blipped phase gradient pulses are applied to incrementally shift a position in k-space such that different time segments of the NMR signal, demarcated by the blipped phase gradient pulses, correspond to different predetermined locations in k-space.

Claims

1. A method for performing nuclear magnetic resonance spectroscopy, the method comprising the steps of: positioning a sample in a homogeneous stationary magnetic field directed along an axis; preparing nuclei in at least a predetermined volume of said sample for resonant emission of an NMR signal and creating said NMR signal, in which said preparing and creating comprises irradiating said sample with at least one radio-frequency excitation pulse in accordance with an MRI sequence preparation module and/or an MRI sequence evolution module; and sensing said NMR signal in the absence of frequency encoding magnetic field gradients such that analysis of said NMR signal yields a chemical shift spectrum from said nuclei; wherein, during said sensing, a plurality of intermittently blipped phase gradient pulses are applied to incrementally shift a position in k-space such that different time segments of said NMR signal, demarcated by said blipped phase gradient pulses, correspond to different predetermined locations in k-space.

2. The method of claim 1, wherein each of said plurality of intermittently blipped gradient pulses has a pulse width in the range of 1 s to 10 s.

3. The method of claim 1, wherein said plurality of intermittently blipped gradient pulses have a same pulse amplitude, wherein said plurality of intermittently blipped gradient pulses are spaced apart in time at regular time intervals, said time interval being smaller than a fourth of the duration of said NMR signal.

4. The method of claim 1, wherein said preparing of the at least one predetermined volume in the sample comprises irradiating said sample with at least one radio-frequency excitation pulse in accordance with an MRI sequence preparation module to prepare said at least one predetermined volume in the form of a predetermined imaging slice.

5. The method of claim 1, wherein said creating of said NMR signal comprises irradiating said sample with at least one radio-frequency refocusing pulse in accordance with an MRI sequence evolution module so as to invert spins of said nuclei, thereby creating said NMR signal in the form of at least part of an echo signal.

6. The method of claim 5, wherein said preparing and/or creating comprises irradiating the sample with at least three selective radiofrequency pulses while simultaneously applying magnetic field gradients, such as to select the predetermined volume at a predetermined location in the sample in accordance with a point resolved spectroscopy or stimulated echo acquisition mode MRI sequence.

7. The method of claim 1, further comprising applying a further blipped phase gradient pulse before creating said NMR signal to shift a k-space location corresponding to the first sensed time segment of said plurality of time segments to a peripheral region of k-space.

8. The method of claim 1, wherein said preparing and/or creating further comprises applying an MRI sequence module for water signal suppression and/or water signal elimination.

9. The method of claim 1, further comprising repeating at least the step of creating the NMR signal to create at least one further NMR signal and the step of sensing the at least one further NMR signal, said NMR signal and said or each at least one further NMR signal corresponding to a different predetermined volume in said sample.

10. The method of claim 1, further comprising, for each of said different time segments of said NMR signal, performing a signal processing step for extending said time segment in time over at least half of the time interval as said NMR signal was sensed.

11. The method of claim 10, further comprising a step of computing a complex conjugate of each extrapolated NMR signal to obtain a conjugated extrapolated NMR signal representative of a mirror location in k-space.

12. The method of claim 10, further comprising, for each of said different time segments of said NMR signal, extrapolating said time segment of the NMR signal in a forward sense and/or a backward sense of time such as to obtain an extrapolated NMR signal, representative of the location in k-space corresponding to said time segment, that is extrapolated over at least half of the time interval as said NMR signal was sensed, in which said extrapolating comprises applying an autoregressive model or generalized autoregressive model to extrapolate said time segment of the NMR segment in said forward sense of time and/or in said backward sense of time.

13. The method of claim 10, furthermore comprising acquiring a reference signal without applying said blipped phase gradient pulses during said acquisition of said reference signal.

14. The method of claim 13, wherein said signal processing step comprises estimating parameters of a predetermined plurality of spectral peaks, said parameters being substantially constant over the time segments, by fitting said reference signal to a signal model representative of said plurality of spectral peaks.

15. The method of claim 14, wherein said parameters comprise a frequency and a damping factor for each of said plurality of spectral peaks.

16. The method of claim 14, wherein said signal processing further comprises using said estimated parameters to generate a noiseless basis representing the predetermined spectral peaks, partitioning said basis into time segments of the basis corresponding to the time segments of said NMR signal and fitting each time segment of the NMR signal to a linear combination of the time segments of the basis to determine complex valued coefficients representative of phase modulated amplitudes of the predetermined plurality of spectral peaks for each time segment.

17. The method of claim 16, comprising extending each time segment in time by using the linear combination of said noiseless basis that corresponds to the complex valued coefficients determined for the time segment as the extended time segment and/or as an extrapolated time segment.

18. The method of claim 1, further comprising applying a Fourier transform to said time segments, or to extrapolated NMR signals, of both k-space and time domain dimensions, to provide chemical species spectra as function of spatial location in the sample.

19. The method of claim 1, further comprising obtaining an image of the sample representative of spatial variations in a chemical constitution of the sample and storing said image in a data storage memory.

20. A control unit comprising a controller configured for controlling a magnetic resonance scanner such as to perform the steps of a method according to claim 1.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1 shows an exemplary method in accordance with embodiments of the present invention.

(2) FIG. 2 illustrates an exemplary method in accordance with embodiments of the present invention, based on a PRESS sequence.

(3) FIG. 3 shows a flow diagram from encoding to reconstruction of data using a fast Fourier transform (FFT) algorithm, in accordance with embodiments of the present invention.

(4) FIG. 4 shows another flow diagram from encoding to reconstruction of data using a fast Fourier transform (FFT) algorithm, in accordance with embodiments of the present invention.

(5) FIG. 5 illustrates an encoding function for encoding time segments, relating to a method in accordance with embodiments of the present invention.

(6) FIG. 6 illustrates another encoding function for encoding time segments, relating to a method in accordance with embodiments of the present invention.

(7) FIG. 7 illustrates an exemplary method in accordance with embodiments of the present invention, based on a FID signal.

(8) FIG. 8 shows exemplary reconstruction results in an example relating to embodiments of the present invention.

(9) FIG. 9 shows a signal processing step in a method in accordance with embodiments of the present invention.

(10) FIG. 10 shows a forward and backward extrapolation in a method in accordance with embodiments of the present invention.

(11) FIG. 11 shows simulated results using spectra without added noise, in an example illustrating embodiments of the present invention.

(12) FIG. 12 shows simulated results using spectra with added noise, in an example illustrating embodiments of the present invention.

(13) FIG. 13 to FIG. 20 show results of creatine, choline and NAA peaks in a simulation example illustrating embodiments of the present invention. FIG. 13, FIG. 15, FIG. 17 and FIG. 19 correspond to noiseless spectra, while FIG. 14, FIG. 16, FIG. 18 and FIG. 20 correspond to noisy data. FIG. 13 and FIG. 14 correspond to data for a number of phase encoded time segments PE=4, FIG. 15 and FIG. 16 correspond to data for PE=8, FIG. 17 and FIG. 18 correspond to data for PE=16 and FIG. 19 and FIG. 20 correspond to data for PE=32.

(14) FIG. 21 to FIG. 24 show results of choline, lactate and NAA peaks in a phantom study illustrating embodiments of the present invention. Data are collected with both CSI and RAPID for comparison. FIG. 21 and FIG. 22 correspond to 1D and 2D datasets with PE=8 and R=8 along AP direction, while FIG. 23 and FIG. 24 correspond to 1D and 2D datasets with PE=16 and R=8 along AP direction. For comparison, CSI data is overlayed on top of the RAPID data in the zoomed voxels.

(15) The drawings are only schematic and are non-limiting. In the drawings, the size of some of the elements may be exaggerated and not drawn on scale for illustrative purposes.

(16) Any reference signs in the claims shall not be construed as limiting the scope.

(17) In the different drawings, the same reference signs refer to the same or analogous elements.

DETAILED DESCRIPTION OF ILLUSTRATIVE EMBODIMENTS

(18) The present invention will be described with respect to particular embodiments and with reference to certain drawings but the invention is not limited thereto but only by the claims. The drawings described are only schematic and are non-limiting. In the drawings, the size of some of the elements may be exaggerated and not drawn on scale for illustrative purposes. The dimensions and the relative dimensions do not correspond to actual reductions to practice of the invention.

(19) Furthermore, the terms first, second and the like in the description and in the claims, are used for distinguishing between similar elements and not necessarily for describing a sequence, either temporally, spatially, in ranking or in any other manner. It is to be understood that the terms so used are interchangeable under appropriate circumstances and that the embodiments of the invention described herein are capable of operation in other sequences than described or illustrated herein.

(20) Moreover, the terms top, under and the like in the description and the claims are used for descriptive purposes and not necessarily for describing relative positions. It is to be understood that the terms so used are interchangeable under appropriate circumstances and that the embodiments of the invention described herein are capable of operation in other orientations than described or illustrated herein.

(21) It is to be noticed that the term comprising, used in the claims, should not be interpreted as being restricted to the means listed thereafter; it does not exclude other elements or steps. It is thus to be interpreted as specifying the presence of the stated features, integers, steps or components as referred to, but does not preclude the presence or addition of one or more other features, integers, steps or components, or groups thereof. Thus, the scope of the expression a device comprising means A and B should not be limited to devices consisting only of components A and B. It means that with respect to the present invention, the only relevant components of the device are A and B.

(22) Reference throughout this specification to one embodiment or an embodiment means that a particular feature, structure or characteristic described in connection with the embodiment is included in at least one embodiment of the present invention. Thus, appearances of the phrases in one embodiment or in an embodiment in various places throughout this specification are not necessarily all referring to the same embodiment, but may. Furthermore, the particular features, structures or characteristics may be combined in any suitable manner, as would be apparent to one of ordinary skill in the art from this disclosure, in one or more embodiments.

(23) Similarly it should be appreciated that in the description of exemplary embodiments of the invention, various features of the invention are sometimes grouped together in a single embodiment, figure, or description thereof for the purpose of streamlining the disclosure and aiding in the understanding of one or more of the various inventive aspects. This method of disclosure, however, is not to be interpreted as reflecting an intention that the claimed invention requires more features than are expressly recited in each claim. Rather, as the following claims reflect, inventive aspects lie in less than all features of a single foregoing disclosed embodiment. Thus, the claims following the detailed description are hereby expressly incorporated into this detailed description, with each claim standing on its own as a separate embodiment of this invention.

(24) Furthermore, while some embodiments described herein include some but not other features included in other embodiments, combinations of features of different embodiments are meant to be within the scope of the invention, and form different embodiments, as would be understood by those in the art. For example, in the following claims, any of the claimed embodiments can be used in any combination.

(25) In the description provided herein, numerous specific details are set forth. However, it is understood that embodiments of the invention may be practiced without these specific details. In other instances, well-known methods, structures and techniques have not been shown in detail in order not to obscure an understanding of this description.

(26) Where in embodiments of the present invention reference is made to nuclear magnetic resonance spectroscopy (NMR spectroscopy), reference is made to methods for determining physical and/or chemical properties of atoms and/or molecules using nuclear magnetic resonance. Intramolecular magnetic fields, which are determined by the molecular configuration, e.g. the electronic structure of a molecule, may affect the resonance frequency, such that, for example, data that relate to the structure, dynamics, reaction state, and/or chemical environment of molecules in a sample may be obtained by NMR spectroscopy.

(27) In a first aspect, the present invention relates to a method for performing nuclear magnetic resonance spectroscopy. This method comprises positioning a sample in a homogeneous stationary magnetic field directed along an axis. The method also comprises preparing nuclei in at least a predetermined volume of the sample for resonant emission of an nuclear magnetic resonance (NMR) signal and creating this NMR signal, in which this preparing and creating comprises irradiating the sample with at least one radio-frequency excitation pulse in accordance with an magnetic resonance imaging (MRI) sequence preparation module and/or an MRI sequence evolution module. The method also comprises, e.g. after this step of creating the NMR signal, sensing the NMR signal in the absence of frequency encoding magnetic field gradients such that analysis of the NMR signal yields a chemical shift spectrum from the nuclei. During this sensing, a plurality of intermittently blipped phase gradient pulses are applied to incrementally shift a position in k-space such that different time segments of the NMR signal demarcated by the blipped phase gradient pulses correspond to different predetermined locations in k-space.

(28) Where in present description reference is made to an MRI sequence, this refers to a plurality of radiofrequency pulses and magnetic gradients that are applied in a predetermined timed and coordinated manner such as to result in a detectable radiofrequency signal or plurality of radiofrequency signals that can be used, by suitable signal processing, to produce at least one image with particular appearance indicative of an internal structure of an object of interest. The MRI sequence may be represented, schematically, in the form of a pulse sequence timing diagram. Such MRI sequence is typically applied in the presence of a strong and substantially homogeneous magnetic field, e.g. a uniform field of at least 0.5 T, preferably at least 1.0 T, e.g. in the range of 1.5 T to 15 T. For example, the MRI sequence may be applied by operating gradient coils and/or radiofrequency antennae, e.g. by operating an MRI apparatus. When current is passed through such coils, in accordance with the MRI sequence, a secondary magnetic field may be created that slightly distorts the main magnetic field in a controlled manner, e.g. such as to shift resonance frequencies of protons in an object under study as a function of position.

(29) Where in present description reference is made to an MRI sequence preparation module, or simply to a preparation module, this refers to a part of an MRI sequence that prepares nuclei for resonant emission, for example by applying gradients to select a predetermined subvolume, e.g. a slice, of the object under study. The preparation module may comprise, for example, at least one 180 inversion pulse, spectral saturation pulses, driven equilibrium pulses and/or spatial saturation pulses. A preparation module may also comprise a timed delay before activating a next module of the MRI sequence, e.g. an evolution module and/or acquisition module. The preparation module may be as simple as applying only a single radiofrequency pulse for spin preparation, yet can also be substantially more complex for particular sequences known in the art.

(30) Where in present description reference is made to an MRI sequence evolution module, or simply to an evolution module, this refers to a part of an MRI sequence that manipulates the emissions of nuclei in the sample under study after a resonant emission of the NMR signal has been created, i.e. initiated by an excitation signal, e.g. in a preparation module of the sequence, e.g. by at least irradiating the sample with at least one radio-frequency excitation pulse in the preparation module. This manipulation may for example comprise applying refocusing radiofrequency pulses to create a train of echo signals. This may, in general, also comprise the application of field gradients, e.g. such as to manipulate the frequency and/or phase properties of the NMR signal(s) after the initial excitation, e.g. to spatially encode the echo signals in such echo train differently.

(31) Where in present description reference is made to the sensing of the NMR signal, this refers to the detection of a radiofrequency signal emitted by the object under study due to a nuclear magnetic resonance excitation, e.g. emitted by protons in the studied object in response to the application of a plurality of magnetic field gradients and/or radiofrequency pulses in accordance with an MRI sequence. This sensing may typically comprise receiving the NMR signal, e.g. using an antenna, and digitizing the received signal, e.g. generating a time series of digital data using an analog to digital converter. This sensing may also comprise other steps as known in the art of NMR signal acquisition, such as applying analog and/or digital filtering, and/or performing other data signal processing steps for, for example, gain control, noise reduction, suppression of other signal sources and/or combining signals acquired simultaneously by different antennae in a suitable manner.

(32) It is an advantage of embodiments of the present invention that short acquisition times can be achieved in NMR spectroscopy, e.g. in magnetic resonance spectroscopic imaging, such as in chemical shift imaging. It is a further advantage of embodiments of the present invention that motion artefacts can be reduced.

(33) Advantageously, a good signal to noise ratio (SNR) can be obtained while collecting the chemical shift information using spatial encoding. This may be achieved by performing a subdivision of an NMR signal, e.g. an MRS signal, acquired from a localized region into a certain number of sub-signals, e.g. the time segments, e.g. sub-echoes, using an inserted blipped phase encoding gradient in one, or several spatial directions. However, since shortening the time-length of the acquired sub-signals, e.g. the time segments of the NMR signal, e.g. sub-echoes, would normally affect both the SNR and the spectral resolution, particularly advantageous embodiments of the present invention may expand the encoded time segments to the full NMR signal length to overcome this problem, e.g. using signal processing to expand the (or each) time segment of the NMR signal in time, e.g. by extrapolation, e.g. using auto-regression methods.

(34) The time duration of the acquisition window of a suitable MRS sequence used in combination with a method in accordance with embodiments of the present invention may determine a maximum number of phase encoding steps that can be inserted, e.g. given a minimum allowed time duration of the time segments. Thus a link can be seen between the spatial and the spectral resolution, and consequently the SNR. Thus, for a single shot acquisition, increasing the spatial resolution would lower the spectral resolution, and vice versa. To ease the stress on the duration of the acquisition window, one can generate multiple NMR signals, e.g. multiple echoes, for the encoding. In such case, the method in accordance with embodiments of the present invention could be comparable to CSI insofar the number of inserted blips would be one per echo, and would reach its full speed when all the necessary blips can be played within a single NMR signal, e.g. a single echo (single shot).

(35) Referring to FIG. 1, an exemplary method 10 for performing nuclear magnetic resonance spectroscopy in accordance with embodiments of the present invention is shown.

(36) The method 10 comprises positioning 11 a sample in a homogeneous stationary magnetic field directed along an axis, e.g. an axis z. The method 10 may also comprise a step of shimming the magnetic field to correct for inhomogeneity of the magnetic field. This may comprise tuning different pulses in three orthogonal principal directions x, y, and z.

(37) The method further comprises preparing 12 nuclei in at least a predetermined volume of the sample for resonant emission of an NMR signal and creating 13 this NMR signal.

(38) This preparing 12 and creating 13 comprises irradiating the sample with at least one radio-frequency excitation pulse in accordance with an MRI sequence preparation module and/or an MRI sequence evolution module.

(39) In a method in accordance with embodiments of the present invention, the step of preparing 12 of the at least one predetermined volume in the sample may comprise irradiating the sample with at least one radio-frequency excitation pulse in accordance with an MRI sequence preparation module to prepare a predetermined imaging slice, e.g. to prepare the at least one predetermined volume in the form of the predetermined imaging slice.

(40) A method in accordance with embodiments of the present invention may be a single-shot method, where only one NMR signal, e.g. a FID or echo, is sufficient to apply all required phase encoding steps, in the form of the intermittently blipped phase encoding gradient, to acquire an entire chemical shift image, e.g. a 2D image, or an entire chemical shift profile, e.g. a 1D image.

(41) However, the method 10 may also comprise repeating 18 at least the step of creating the NMR signal to create at least one further NMR signal and the step of sensing the at least one further NMR signal.

(42) The NMR signal and the or each at least one further NMR signal may correspond to a different predetermined volume in the sample. For example, the repeating 18 may comprise selecting different imaging slices in different repetitions. It is an advantage of embodiments of the present invention that efficient spectroscopic 2D imaging in a slice, or 3D imaging by iterating over a plurality of slices, can be obtained by isolating nuclear spins in a planar slice of the sample.

(43) However, this repeating may also comprise selecting the same volume, e.g. imaging slice, in multiple repetitions, in which a different trajectory in k-space is selected in different repetitions by varying blip amplitudes and/or phase encoding gradients used for generating the intermittently blipped phase gradient pukes, as described further hereinbelow. Furthermore, both different selected slices and different k-space trajectories may be combined in embodiments comprising such repetitions, as shall be clear to the skilled person.

(44) The step of creating 13 of the NMR signal may comprise irradiating the sample with at least one radio-frequency refocusing pulse in accordance with an MRI sequence evolution module so as to invert spins of the nuclei, thereby creating 16 the NMR signal in the form of at least part of an echo signal, e.g. a spin echo signal or a decaying half of an echo signal, e.g. creating the NMR signal in the form of at least part of a spin echo signal.

(45) It is an advantage of embodiments that multiple points in k-space, and hence multiple voxel locations, can be acquired from a single spin echo signal. For example, a larger number of encoding steps can be achieved than the number of echo signals that are generated, e.g. the number of encoding steps may be a multiple of the number of generated echo signals. Furthermore, the number of voxels may be doubled for the same number of data points collected in k-space using the principle of complex conjugate of the Fourier Transform, as discussed further hereinbelow.

(46) For example, the steps of preparing 12 and creating 13 may comprise irradiating the sample with at least three selective radiofrequency pulses while simultaneously applying magnetic field gradients, such as to select the predetermined volume at a predetermined location in the sample in accordance with a point resolved spectroscopy (PRESS) or stimulated echo acquisition mode MRI (STEAM) sequence.

(47) In PRESS, a spin echo signal may be collected to obtain chemical shift information from a box-shaped volume, e.g. a voxel, located at the intersection of three orthogonal planes and selected by using the three RF pulses. The echo signal thus originates from only the box-shaped volume submitted to the three RF pulses. Phase encoding gradients can be activated in between the three RF pulses to enable spatial encoding of the spins during their evolution according to the Fourier Transform principle. For example, in a prior-art PRESS sequence, one needs as many phase encoding gradient steps as the number of voxels that need to be sampled, for each direction to be sampled. The prior-art approach may thus require the same number of echoes for a complete and correct spatial encoding. Therefore, the running time of the prior-art sequence may be equal to the time duration of the base sequence times the number of voxels required in each direction.

(48) However, advantageously, a method in accordance with embodiments of the present invention may reduce the acquisition time of a PRESS sequence, or equivalently of another suitable sequence such as STEAM, by a factor that is at least equal to the number of time segments generated for each NMR signal, e.g. by sub-dividing the acquired echo signal into a number of sub-signals, e.g. sub-echoes, e.g. the time segments of the NMR signal, and inserting the phase encoding gradient in between these sub-signals in the form of a short blip.

(49) The steps of preparing 12 and creating 13 may further comprise applying an MRI sequence module for water signal suppression and/or water signal elimination, such as a CHESS sequence module. For example, because water molecules contain hydrogen and relative concentrations of water to an arbitrary metabolite can be about 10000:1 or even higher in biological samples, the water signal may be advantageously suppressed to allow the metabolite peaks to become discernible in the generated spectra. However, in other embodiments of the present invention, MRS may be achieved without necessarily performing such water suppression.

(50) The method 10 further comprises, e.g. after creating 13 the NMR signal, sensing 14 the NMR signal in the absence of frequency encoding magnetic field gradients such that analysis of the NMR signal yields a chemical shift spectrum from the nuclei. It is known in the art that frequency encoding will adversely affect the chemical shift spectral information in most cases, e.g. unless the MRI sequence is specifically designed to allow the preservation of this spectral information. For example, in a proton echo-planar spectroscopic imaging (PEPSI) sequence as known in the art, a succession of zig-zag readout gradients can be applied to preserve chemical shift information. However, it is an advantage of embodiments of the present invention that a fast and efficient NMR spectroscopy sequence is provided that does not require the application of frequency encoding gradients during sensing that could potentially deteriorate the chemical shift information. It is another advantage that a simple NMR sequence is provided, e.g. that does not require careful optimization of frequency encoding gradients to allow spectroscopic data collection.

(51) During the sensing 14, a plurality of intermittently blipped phase gradient pulses are applied 15 to incrementally shift a position in k-space such that different time segments of the sensed NMR signal correspond to different predetermined locations in k-space. These time segments are demarcated by the blipped phase gradient pulses.

(52) Thus the NMR signal, such as a single MRI echo signal, the decaying part of a single MRI echo signal or a single FID signal, is sub-divided into the time segments. An echo signal may, for example, describe the spin response of the protons of metabolites in a tissue to a radio frequency perturbation. The number of sub-signals, e.g. sub-echoes, e.g. the time segments of an echo signal, may be equal to the number of encoding steps, and thus equal to the number of voxels that are acquired from a single echo signal. Advantageously, this allows performing spatial encoding during the acquisition of a single echo signal, or of half of a single echo signal, e.g. the decaying half of a single echo signal, or likewise, of another NMR signal, such as a single FID signal.

(53) For example, assuming that the applied blips are substantially identical in length, amplitude and encoding gradient direction, a first blip may cause a jump from an initial position in k-space, e.g. from the centre of k-space, in a direction determined by phase encoding gradient used for generating the blip. A second blip may then cause another jump over a same distance as the first jump to arrive at a k-space position twice as far removed from the initial position, e.g. the centre of k space, than the resulting position of the first blip, and so on.

(54) For example, each of the plurality of intermittently blipped gradient pulses may have a pulse width in the range of 1 s to 10 s, such as in the range of 2 s to 6 s, e.g. 3 s to 5 s, e.g. a pulse width of 4 s. For example, the pulse width, i.e. the duration of each blip, may be preferably as short as can be achieved, or in the shorter end of the range that can be achieved, in view of gradient switching performance of the NMR system used for generating the plurality of blipped gradient pulses.

(55) In embodiments according to the present invention, the plurality of intermittently blipped gradient pulses may have a same pulse width.

(56) The plurality of intermittently blipped gradient pulses may have a same pulse amplitude.

(57) The plurality of intermittently blipped gradient pulses may be applied at regular time intervals, e.g. at a constant time offset with respect to a previous blip. Thus, the plurality of intermittently blipped gradient pulses may be spaced apart in time at regular time intervals.

(58) For example, the time interval between consecutive blips, e.g. being a constant time interval in accordance with some embodiments or a variable time interval in other embodiments, may be preferably smaller than a fourth of the duration of the NMR signal. This time interval is preferably shorter, to allow more encoded segments, and thus more sampled voxels, yet long enough to obtain sufficient time samples points per segment.

(59) For example, blips of the same amplitude, and the same pulse width, i.e. the same length in time, may be activated at regular intervals, to obtain signal segments of equal length and corresponding to a uniform sampling of k-space points (e.g. uniformly sampled along a single direction in k-space corresponding to the phase encoding gradient used for applying the blips). Signal segments of equal length may have the advantage of providing uniformity in the extrapolation result, e.g. to avoid phase shifts in reconstructed spectra. Furthermore, a uniform sampling of k-space points, such as a one-dimensional Cartesian sampling or a two-dimensional Cartesian sampling, may provide data that can be easily processed, e.g. displayed as an image or image volume.

(60) However, the amplitude, pulse width and the time interval between blips may be non-uniform and/or non-constant, e.g. variable, for example such as to achieve a custom, non-uniform sampling in k-space, e.g. a non-cartesian sampling, and/or to collect more data points for some segments as opposed to other segments. For example, a polar sampling can be achieved, in which the steps of a method in accordance with embodiments of the present invention are repeated for different phase encoding gradients, e.g. to sample along different lines through the centre of k-space, e.g. spaced at regular angular intervals. Furthermore, for each sampled line in k-space, the distance in k-space between locations corresponding to consecutively applied blips may be varied, e.g. by selecting appropriate blip amplitudes, such as to, for example, sample more densely near the centre in k-space.

(61) A method in accordance with embodiments of the present invention may also comprise applying 17 a further blipped phase gradient pulse before creating said NMR signal to shift a k-space location corresponding to the first sensed time segment of the plurality of time segments to a peripheral region of k-space, e.g. to an edge of k-space.

(62) It is an advantage of embodiments that multiple points in k-space, and hence multiple voxel locations, can be acquired from a single NMR signal, e.g. created by a single preparation-evolution module of excitation and/or refocusing prior to the signal sensing. For example, the NMR signal may decay over a relatively long time, such that several segments can be collected. For example, the number of different voxel measurements that may be collected may thus be substantially proportional to the decay time.

(63) However, the NMR signal, e.g. a free induction decay (FID) or echo signal, or part thereof, may have a frequency content that changes over time. Therefore, a first time segment may contain substantially all the resonance frequencies generated in the sample (e.g. in a selected volume of the sample), while a second time segment that is acquired later than the first time segment may contain less information due to the decay of signal components that have a short decay time. For example, near the end of the NMR signal, signal components, such as components corresponding to lipids, may disappear. Thus, lipid signal peaks may be reduced in the last segment or segments. However, this can be advantageous in applications where such lipid peaks may be a source of problems, e.g. causing baseline distortion and artefacts in data quantification. Nevertheless, most signal components of interest may have a decay time that is sufficiently long to allow recovery of such signal components from each segment and thus for each sampled k-space location. Furthermore, a deconvolution technique may be applied in the k-space domain to allow signal peak recovery in all sampled k-space locations, including those signals peaks characterized by a rapid decay.

(64) FIG. 2 illustrates an exemplary method in accordance with embodiments of the present invention. In a top part 201 of FIG. 2, a simulated NMR signal 20 is shown that corresponds to a spin echo generated from a PRESS sequence 21, while, similarly, in a bottom part 202 of FIG. 2, a simulated NMR signal 20 is shown that corresponds to only the decaying half of a spin echo generated from a PRESS sequence 21. The PRESS sequence allows the generation of such an echo from a voxel volume with variable size and location. The encoding scheme in accordance with embodiments of the present invention may thus be inserted into a PRESS sequence. The PRESS sequence provides an echo from a voxel with variable size and location.

(65) The intermittently blipped phase gradient pulses 22, e.g. the gradient blips, may be generated at equidistant time locations. These blips divide the collected NMR signal 20 into a number of phase encoded (PE) time segments 23, e.g. sub-echoes, with equal time duration. A further blip gradient 24, e.g. which may typically have a larger amplitude than a gradient blip 22, may be played before the second refocusing RF pulse 25 to move the cursor to the edge of the k-space.

(66) Each time segment is thus phase encoded by the blip gradient and represents one k-space data point, as illustrated in FIG. 3. The collected PE k-space data points may thus fill one line in a given spatial direction, determined by the activation axis of the blipped gradient, e.g. the x axis corresponding to the gradient Gx in present example. FIG. 3 further illustrates a flow diagram from encoding to reconstruction of data using a fast Fourier transform (FFT) algorithm. The results 31 after autoregression show both the reconstructed data and the original data on which the simulation was performed for comparison, in which the original data is shown as the upper plot line of each pair, and the reconstructed data is shown as the lower plot line of each pair. This example relates to the exemplary simulated NMR signal corresponding to a spin echo signal as illustrated in the top part 201 of FIG. 2.

(67) FIG. 4 shows another exemplary flow diagram from encoding to reconstruction of data using a fast Fourier transform (FFT) algorithm, in accordance with embodiments of the present invention. In this example, the results 31 after applying a signal processing as discussed further hereinbelow show both the reconstructed data and the original data on which the simulation was performed for comparison, in which the original data is shown as the upper plot line of each pair, and the reconstructed data is shown as the lower plot line of each pair. This example relates to the exemplary simulated NMR signal corresponding to the decaying half of a spin echo signal as illustrated in the bottom part 202 of FIG. 2.

(68) FIG. 5 illustrates the encoding function for encoding the time segments, relating to the first example illustrated in the top part 201 of FIG. 2 and the exemplary flow diagram in FIG. 3, showing both imaginary 402 and real 401 components of the function.

(69) Likewise, FIG. 6 illustrates the encoding function for encoding the time segments, relating to the second example illustrated in the bottom part 202 of FIG. 2 and the exemplary flow diagram in FIG. 4.

(70) The blip gradients may be activated along the x axis and may be played at constant time intervals during the acquisition of the echo signal or, respectively, the decaying half of the echo signal to encode the localized spins based on their spatial locations according to the harmonics shown in respectively FIG. 5 and FIG. 6. The NMR signals of these spatially encoded spins are summed, as shown in plots 32 in respectively FIG. 3 and FIG. 4, and sequentially collected in form of the sub-signals 33, thus filling a k-space x line. A Fourier transform may be applied in both the spatial and spectral domain to produce the localized spectra 34.

(71) In another exemplary method in accordance with embodiments of the present invention, the NMR signal 20 may be generated using a radiofrequency pulse along with a slice selection gradient and a rephasing gradient lobe 51, as a Free Induction Decay (FID) signal from a localized region, which is divided into time segments, e.g. sub-FIDs, in a similar manner as the previous examples using an echo signal or part thereof. A flowchart illustrating the signal generation for this example is shown in FIG. 7.

(72) After extrapolating the time segments, e.g. the sub FIDs, in a similar manner as in previous example, a Fourier transform can be performed to reconstruct the original signals. An exemplary result of such reconstruction is shown in FIG. 8. The original reference data is shown as the lower plot line of each pair, and the reconstruction using FFT and auto-regression is shown as the upper plot line of each pair. It can be concluded that in this exemplary simulation the metabolites are localized accurately in both spectral and spatial domain.

(73) Referring back to FIG. 1, a method 10 in accordance with embodiments of the present invention may also comprise, for each of the different time segments of the NMR signal, a signal processing step 19 for processing the time segment of the NMR signal to extend it in time, e.g. to a full time length, e.g. such as to obtain an extrapolated NMR signal representative of the location in k-space corresponding to the time segment that is being extended.

(74) Such processing may be performed to expand the sub-signal, e.g. sub-echo, e.g. in embodiments where the NMR signal is an echo signal or a part of an echo signal, e.g. the decaying part of the echo signal. In such embodiments, each time segment, e.g. each sub-signal, may be expanded to the time-length of a reference signal. For example, an FID signal may be collected with no gradient activation, e.g. without applying blips. This FIB may be used as reference data for the signal processing. Alternatively, another echo or half-echo signal may be used as the reference signal, e.g. an (half-)echo signal without activation of the blip gradients to shift the k-space coordinate, e.g. at a k=0 origin.

(75) For example, this signal processing 19 may comprise extrapolating the time segment of the NMR signal in a forward and/or a backward sense of time such as to obtain an extrapolated NMR signal, representative of the location in k-space corresponding to the time segment that was extrapolated. This extrapolated NMR signal may be extrapolated over at least half of the time interval as the NMR signal was sensed, e.g. over the half of the time covered by the NMR signal, e.g. a half echo, or over the time covered by the NMR signal in its entirety. Thus, each segment of the NMR signal may be extrapolated over at least half of the time frame of the NMR signal.

(76) An extrapolation, such as an auto-regression method, may be used to extrapolate and expand sub-echoes, e.g. in embodiments where the NMR signal is an echo. In such embodiments, each time segment, e.g. each sub-echo, may be expanded to the time-length of a half echo. For example, one half echo may be collected with no gradient activation, e.g. without applying blips. This half echo may be used as reference data for the extrapolation.

(77) For example, the blip gradients may be activated to collect PE/2 encoded sub-echoes during the increasing part of the echo. Thereafter, the blip gradients may be switched off to acquire the decaying half echo to be used as a reference in the extrapolation process. However, it shall be clear to the skilled person that other embodiments may relate to the use the increasing part of the echo as reference and the decreasing part for acquiring the time segments.

(78) The principle of the complex conjugate of the Fourier Transform, referred to further hereinbelow, may then be used on the collected PE/2 encoded sub-echoes to recover the conjugated other PE/2 ones. This is similar to the HASTE MRI method, known in the art, where half of phase encoding k-space data are collected, then complex conjugated to recover the opposite half. In this manner, the acquisition time may be kept short while the sub-echo data are extrapolated to the full time duration of a half echo. As a consequence, problems of peak broadening may be addressed.

(79) A method in accordance with embodiments of the present invention may further comprise a step of computing a complex conjugate of each of the extrapolated NMR signals to obtain a conjugated extrapolated NMR signal representative of a mirror location in k-space. For example, the extrapolation 19 may be performed on each time segment of the NMR signal to extend the segment data to cover the same length as the full NMR signal, and, in a following step, each extrapolated segment may be complex conjugated to generate another half of k-space data. For example, a half-Fourier approach may be applied, taking the relation S(k)=S*(k) into account, in which S refers to the extrapolated NMR signal, S* to its complex conjugate and k to k-space coordinates.

(80) In a step of signal processing 19 in a method in accordance with embodiments of the present invention, the time segment of the NMR signal may be expanded, e.g. extended, e.g. extrapolated. As described hereinabove, this extrapolation may comprise, for example, an auto-regression or similar regression algorithm. Alternatively, the signal processing may comprise a model-based time-expansion method as described hereinbelow and illustrated by the exemplary flow diagram in FIG. 9.

(81) Each time segment of the NMR signal, e.g. each acquired sub-signal, can be described by a combination of two sets of functions, in which the first functions are representative of the amplitude of a predetermined, e.g. selected, plurality of spectral peaks, for example related to a predetermined plurality of metabolites, while the second functions are a time-shifted functions representative of other spectral parameters than the amplitude of the predetermined spectral peaks, e.g. representative of the remaining chemical shift parameters, such as frequency and damping factor, of the predetermined spectral peaks.

(82) The spectral peaks in the first functions may be phase modulated due to the activation of the blipped gradients (e.g. thus being dependent on the k-space coordinate). This first function therefore can differ from one sub-signal to another.

(83) The second functions may be substantially identical across the acquired sub-signals. Therefore, the constant parameters of the spectral peaks, e.g. frequency and damping factor, can be estimated by fitting a full time length reference signal acquired at the center of k-space (k=0) to an appropriate signal model representing the selected spectral peaks.

(84) For example, a time segment, representative of a sample of k-space data along a direction x in which the blipped phase encoding gradients are applied, may be expressed by:
s(k,t.sub.j)=s(x,t.sub.j).Math.e.sup.(i.2..k.x).Math.dx
with
s(x,t.sub.j)=.sub.m=1.sup.MA.sub.m.sub.x.Math.e.sup.t.sup.j.sup.(d.sup.m.sup.+2if.sup.m.sup.)+(t.sub.j),

(85) in which t.sub.j is the constant time duration of the j.sup.th acquired sub-signal, proportional to the number of data points of the sub-signal N.sub.sub (N.sub.sub=N.sub.full/R), which may run from 1+(j1).N.sub.sub to j.N.sub.sub. Furthermore, N.sub.full refers to the total number of points of the acquired NMR signal, and R may be referred to as an acceleration factor, in which R is equal to the number of the sub-signals into which the NMR signal is partitioned. The variable j sets the time occurrence of each sub-signal, i.e. indexes the time segments in chronological order, and assumes integer values from 1 to R. M is the number of predetermined peaks, e.g. indicative of predetermined signal metabolites, with A, d and f describing amplitude, damping factor and frequency for each peak m, e.g. for each metabolite, and (t) denotes a noise term. By combining the two equations hereinabove, the encoded sub-signal s, for a data point k.sub.j in k-space, can be written as:

(86) $\begin{array}{c}s\left(k_j,t_j\right)=\underset{m=1}{\overset{M}{.Math.}}{A}_{m_x}.Math.e^{-t_j\left(d_m+2{\mathrm{if}}_m\right)}.Math.e^{\left(i.Math.2.Math.\mathrm{kj}.Math.x\right)}.Math.\mathrm{dx}+\mathrm{.Math.}\left(t_j\right)\\ \underset{m=1}{\overset{M}{.Math.}}{A}_{m_{x_j}}^{k_j}.Math.{}_{\mathrm{sub}}\left(m,t_j\right)+\mathrm{.Math.}\left(t_j\right)\end{array}$

(87) where

(88) $A_{m_{k_j}}^{{}_{k_j}}={\mathrm{Am}}_x.Math.e^{1.2k_{j}x}\mathrm{dx}$
are the first functions specific to each sub-signal due to their dependence on the k-space variable

(89) $k_j(k_j=\left(\frac{}{2.Math.}\right){}_0^{{}_G}{G}_x\left(t\right).Math.\mathrm{dt},$
with .sub.G representative of the gradient blip duration) and the second functions .sub.sub(m,t.sub.j)=e.sup.t.sup.j.sup.(d.sup.m.sup.+2if.sup.m.sup.) are dependent on the parameters d.sub.m and f.sub.m of the spectral components during the time t.sub.j. and independent of the k-space coordinate.

(90) The constant parameters of the predetermined spectral peaks, e.g. frequency and damping factor, can be estimated by fitting a full time length reference signal acquired at the center of k-space (k=0) to the signal model representing the predetermined spectral peaks.

(91) Therefore, the step of signal processing may further comprise estimating parameters of said spectral peaks that are substantially constant over the time segments, e.g. the other spectral parameters than amplitude, e.g. parameters such as frequency and damping factor, by fitting a reference signal, e.g. being defined over a longer time interval than each of the time segments cover, e.g. being defined over a full time scale of the NMR signal to a signal model representative of the predetermined spectral peaks. The reference signal, e.g. a FID signal, may be acquired 91, e.g. a pre-acquired reference signal may be stored and used in the step of signal processing 19, at the center of k-space (e.g. where k=0).

(92) For example, To estimate d.sub.m and f.sub.m, a full time length reference signal may be acquired at k.sub.j=0 and fitted to the corresponding full time length signal model given by:

(93) $s\left(k_j=0,N_{\mathrm{full}}\right)=\underset{m=1}{\overset{M}{.Math.}}{A}_{m_{k_j=0}}^{{}_{k_j}=0}\left(m,N_{\mathrm{full}}\right)$

(94) where t.sub.j is replaced by the time length N.sub.full for expanding the time segments to, and k.sub.j is set to 0. The normalized basis .sub.full(m,N.sub.full) over the time length N.sub.full may be described by an N.sub.fullM matrix:

(95) ${}_{\mathrm{full}}\left(m,N_{\mathrm{full}}\right)=\left(\begin{array}{ccc}{e}^{-t_1\left(d_1+2{\mathrm{if}}_1\right)}& \mathrm{.Math.}& e^{-t_1\left(d_M+2{\mathrm{if}}_M\right)}\\ \mathrm{.Math.}& & \mathrm{.Math.}\\ e^{-t_{N_{\mathrm{full}}}\left(d_1+2{\mathrm{if}}_1\right)}& \mathrm{.Math.}& e^{-t_{N_{\mathrm{full}}}\left(d_M+2{\mathrm{if}}_M\right)}\end{array}\right)$

(96) and for each t.sub.j, a corresponding sub-basis function .sub.sub(m,t.sub.j) reduced to N.sub.subM may be extracted:

(97) ${}_{\mathrm{sub}}\left(m,t_j\right)=\left(\begin{array}{ccc}{e}^{-t_{\left(1+\left(j-1\right)*N_{\mathrm{sub}}\right)}\left(d_1+2{\mathrm{if}}_1\right)}& \mathrm{.Math.}& e^{-t_{\left(1+\left(j-1\right)*N_{\mathrm{sub}}\right)}\left(d_M+2{\mathrm{if}}_M\right)}\\ \mathrm{.Math.}& & \mathrm{.Math.}\\ e^{-t_{\left(j*N_{\mathrm{sub}}\right)}\left(d_1+2{\mathrm{if}}_1\right)}& \mathrm{.Math.}& e^{-t_{\left(j*N_{\mathrm{sub}}\right)}\left(d_M+2{\mathrm{if}}_M\right)}\end{array}\right)$

(98) Therefore, an encoded sub-signal s(k.sub.j,t.sub.j) can be modeled by a noiseless matrix form as:

(99) $s\left(k_j,t_j\right)={}_{\mathrm{sub}}\left(m,t_j\right).Math.\left(\begin{array}{c}{A}_{1_{k_j}}^{{}_{k_j}}\\ \mathrm{.Math.}\\ \mathrm{.Math.}\\ A_{M_{k_j}}^{{}_{k_j}}\end{array}\right)$

(100) Or, by substitution by the equation for .sub.sub(m,t.sub.j) provided hereinabove, s(k.sub.j,t.sub.j) can be expressed as:

(101) $\left(\begin{array}{c}{s}_{(1+\left(\left(j-1\right)*N_{\mathrm{sub}}\right)}\\ \mathrm{.Math.}\\ s_{\left(j*N_{\mathrm{sub}}\right)}\end{array}\right)=\left(\begin{array}{ccc}{e}^{-t_{\left(1+\left(j-1\right)*N_{\mathrm{sub}}\right)}\left(d_1+2{\mathrm{if}}_1\right)}& \mathrm{.Math.}& e^{-t_{\left(1+\left(j-1\right)*N_{\mathrm{sub}}\right)}\left(d_M+2{\mathrm{if}}_M\right)}\\ \mathrm{.Math.}& & \mathrm{.Math.}\\ e^{-t_{\left(j*N_{\mathrm{sub}}\right)}\left(d_1+2{\mathrm{if}}_1\right)}& \mathrm{.Math.}& e^{-t_{\left(j*N_{\mathrm{sub}}\right)}\left(d_M+2{\mathrm{if}}_M\right)}\end{array}\right).Math.\left(\begin{array}{c}{A}_{1_{k_j}}^{{}_{k_j}}\\ \mathrm{.Math.}\\ \mathrm{.Math.}\\ A_{M_{k_j}}^{{}_{k_j}}\end{array}\right)$

(102) where the complex amplitudes

(103) $\left(\begin{array}{c}{A}_{1_{k_j}}^{{}_{k_j}}\\ \mathrm{.Math.}\\ \mathrm{.Math.}\\ A_{M_{k_j}}^{{}_{k_j}}\end{array}\right)$
can be estimated using a matrix inversion or approximate matrix inversion algorithm, e.g. a pseudo-inverse matrix operation, and then injected into the equation

(104) 0$s\left(k_j,N_{\mathrm{full}}\right)=\underset{m=1}{\overset{M}{.Math.}}{A}_{m_{k_j}}^{{}_{k_j}}\left(m,N_{\mathrm{full}}\right)$
to construct full time noiseless encoded MRS signals. A Fourier Transform can be applied on both spatial and spectral domain to obtain the spatial distribution of the metabolite peaks.

(105) The processing method can be comparable to CSI when the number of inserted blips is one per acquired signal (R=1), and reaches a high speed gain when all the necessary blips are played within a single acquisition (R=PE), with PE equal to the total number of phase encoding steps. Therefore, the acquisition time reduction is proportional to R.

(106) The efficiency of a method in accordance with embodiments of the present invention may be linked to the accuracy of the estimation of the sub-signal parameters. If frequency and damping factors may be assumed to be generally accurately estimated from the reference signal, as described hereinabove, the estimation of the phase modulated amplitudes values, from which precise metabolite maps may be derived, might still be challenging due to the decreasing SNR of the sequentially acquired sub-signals. The SNR per voxel of a 1D CSI is given by:

(107) ${\mathrm{SNR}}_{\mathrm{CSI}}=\frac{A.Math.e^{-{}_G/T_2^{*}}.Math.T_2^{*}.Math.f.Math.\left(1-e^{-N_{\mathrm{full}}/T_2^{*}}\right)}{a.Math.\sqrt{N_x}.Math.\sqrt{N_{\mathrm{full}}f}}$

(108) where A denotes the peak amplitude of a metabolite. N.sub.x is the number of voxels along x direction, f is the bandwidth of the MRS spectrum, and a being a constant term contributing to the noise standard deviation. The term (1e.sup.N.sup.full.sup./T*.sup.2) is related to the duration of the acquisition window. Given the acceleration is performed along x direction, one can derive the SNR per voxel of this exemplary method in accordance with embodiments of the present invention:

(109) ${\mathrm{SNR}}_{\mathrm{RAPID}}=\frac{\begin{array}{c}A.Math.e^{-{}_G/T_2^{.Math.}}.Math.(1+e^{-N_{\mathrm{full}}/R.Math.T_2^{*}}+e^{-2.Math.N_{\mathrm{full}}/R.Math.T_2^{*}}+\mathrm{.Math.}+\\ e^{-\left(R-1\right)N_{\mathrm{full}}/R.Math.T_2^{*})}.Math.T_2^{*}.Math.f.Math.\left(1-e^{-N_{\mathrm{full}}/R.Math.T_2^{*}}\right)\end{array}}{a.Math.R.Math.\sqrt{N_x}.Math.\sqrt{N_{\mathrm{sub}}f}}.$

(110) Which can be rewritten as:

(111) ${\mathrm{SNR}}_{\mathrm{RAPID}}=\frac{{\mathrm{Ae}}^{-{}_G/T_2^{*}}.Math.T_2^{*}f.Math.\left(1-e^{-N_{\mathrm{full}}/T_2^{*}}\right)}{a.Math.R.Math.\sqrt{N_x}.Math.\sqrt{N_{\mathrm{full}}f/R}}$

(112) Therefore:

(113) ${\mathrm{SNR}}_{\mathrm{RAPID}}=\frac{{\mathrm{SNR}}_{\mathrm{CSI}}}{\sqrt{R}}$

(114) The sensitivity for RAPID is known to approach .sub.RAPID=.sub.CSI. The SNR can be increased by shortening N.sub.full, which could be set to the shortest T.sub.2* value of the collected metabolite signals (e.g: N.sub.fullT*.sub.2,shortest). However, reducing the number of points N.sub.full, lowers the acceleration factor R since a minimum data length of N.sub.sub is required to allow for an accurate estimation of the phase modulated amplitudes

(115) $A_{m_{k_j}}^{{}_{k_j}}.$
Therefore, to allow for higher acceleration and minimise the SNR loss, N.sub.full may generally be set equal to the longest T*.sub.2 value of the collected metabolite signals (N.sub.fullT*.sub.2,longest).

(116) In a processing method in accordance with embodiments of the present invention, a reference signal and the sub-signals may be used to estimate the parameters and to produce noiseless reconstructed encoded full time signals, which can be used to create metabolite maps. As a consequence, both SNR and spectral resolution may be increased.

(117) For an accurate estimation of the signal parameters, both N.sub.full, and R are preferably selected appropriately, as described hereinabove. A robust post-processing technique, e.g. a HLSVD-PRO method, may preferably be applied to achieve both data denoising and signal parameter estimation.

(118) The step of signal processing may further comprise using the estimated parameters to generate a noiseless basis representing the predetermined spectral peaks over a longer time interval than each of the time segments cover individually, e.g. over the full time scale of the NMR signal, e.g. in which the amplitudes of the predetermined spectral peaks are normalized.

(119) The step of signal processing may further comprise partitioning the basis into time segments of the basis corresponding to the time segments of the NMR signal.

(120) The step of signal processing may further comprise fitting each time segment of the NMR signal to a linear combination of the time segments of the basis to determine complex valued coefficients representative of phase modulated amplitudes of the predetermined spectral peaks for each time segment.

(121) For example, the normalized basis may be time-partitioned into sub-bases similar to the encoded sub-signals, which may be represented in a matrix form, and may be fitted to the corresponding collected sub-signals. Thus, the corresponding phase modulated amplitudes of the peaks can be estimated.

(122) The step of signal processing may further comprise extending each time segment in time, e.g. extrapolating each time segment in time, e.g. in a forward and/or backward sense of time, by using the linear combination of the basis corresponding to the complex valued coefficients for the time segment as the extended and/or extrapolated time segment.

(123) For example, once the complex amplitudes of the peaks are estimated for each time segment, they can be inserted into the full basis, e.g. the normalized basis without time-partitioning applied thereto, to generate a full time-length encoded MRS signal for each time segment.

(124) Alternatively, the step of extrapolating 19 may comprise applying an autoregressive model or a generalized autoregressive model to extrapolate the time segment of the NMR segment in the forward sense of time and/or in the backward sense of time.

(125) In a method in accordance with embodiments, the NMR signal may be an echo signal, and the step of extrapolating 19 may comprise applying a first autoregressive model to extrapolate the time segment in a time sense away from a time of occurrence of the maximum of the echo signal envelope, thus producing a first extrapolated signal, and applying a second autoregressive model to extrapolate the first extrapolated signal in a time sense toward the time of occurrence of the maximum of the echo signal, thus producing the extrapolated NMR signal. The same auto regressive model to extend a time segment in a backward sense towards the maximum of the echo and in a forward sense towards the end of the echo may be used.

(126) The order of the autoregressive model may, for example, be selected as the order that minimizes an error variance, an Akaike information criterion, a Bayesian information criterion, another information-theoretic information criterion or a similar figure-of-merit criterion. For example, the model order may be selected by fitting a plurality of autoregressive models of different order to a reference signal. The model order which results in lowest error variance, e.g. as used to determine an Akaike information criteria, may be selected as the model order. The model may thus generate autoregressive coefficients, for the selected model order, that may be used to extrapolate the time segments.

(127) For example, Burg's method may be applied to determine filter coefficients in an auto-regression model. This method has been applied in signal processing, as described by Uchiyama et al., in In Vivo potassium-39 NMR spectra by the burg maximum-entropy method, Journal of Magnetic Resonance (1969), 1991. 94(3): p. 449-454.

(128) An auto-regression model can be formulated as:
x(n)=a.sub.1x(n1)a.sub.2x(n2) . . . a.sub.Mx(nM)+(n), . . . n=1,2, . . . , N, . . .
where x(n) is a reference signal, for example one half of an echo signal, M is the order of the prediction filter polynomial a=[a.sub.1,a.sub.2,a.sub.3, . . . ,a.sub.M], and v(n) is white Gaussian noise. The filter coefficients can be calculated using Yule Walker equations:

(129) $\left[\begin{array}{cc}r\left(0\right)& r\left(-1\right)\mathrm{.Math.}\mathrm{.Math.}\mathrm{.Math.}r\left(-M+1\right)\\ r\left(1\right)& r\left(0\right)\mathrm{.Math.}\mathrm{.Math.}\mathrm{.Math.}r\left(-M+2\right)\\ & \mathrm{.Math.}\\ r\left(N-1\right)& r\left(N-2\right)\mathrm{.Math.}\mathrm{.Math.}\mathrm{.Math.}r\left(0\right)\end{array}\right]\left[\begin{array}{c}{a}_1\\ a_2\\ a_3\\ \mathrm{.Math.}\\ a_M\end{array}\right]=\left[\begin{array}{c}r\left(1\right)\\ r\left(2\right)\\ r\left(3\right)\\ \mathrm{.Math.}\\ r\left(M\right)\end{array}\right]$
where r=[r(0),r(1),r(2), . . . ,r(M)], are the autocorrelation estimate of x(n).

(130) These equations may, for example, be solved using the Levinson-Durbin recursion method, e.g. as implemented in MATLAB. The Akaike Information Criterion (AIC) may be used to determine the model order. The model order is used in the AR model described hereinabove and the filter coefficients may be calculated using the Burg's method for better accuracy in the determination of the coefficients. The model along with the filter coefficients can thus be used to extrapolate the acquired time segments of the NMR signal in both forward and backward directions as demonstrated in FIG. 10. A time segment 23 may be first extrapolated 71 in a forward direction and then extrapolated 72 in a backward direction. For comparison, the reference signal 73, e.g. a half echo, is shown.

(131) However, embodiments of the present invention are not necessarily limited to extrapolation using auto-regressive methods. For example, other methods, which are known in the art for the purpose of extrapolation, may be applied, such as singular value decomposition (SVD) based methods for solving the filter coefficients as described in Barkhuijsen et al., Retrieval of frequencies, amplitudes, damping factors, and phases from time-domain signals using a linear least-squares procedure, Journal of Magnetic Resonance (1969), 1985. 61(3): p. 465-481 and in Lin et al., A Novel Detection-Estimation Scheme for Noisy NMR Signals: Applications to Delayed Acquisition Data, Journal of Magnetic Resonance, 1997. 128(1): p. 30-41.

(132) Furthermore, ARMA based modelling may also be applied, e.g. as has been used for spectroscopy signals as mentioned in Tomczak et al., A subband ARMA modeling approach to high-resolution NMR spectroscopy, Journal of Magnetic Resonance, 2002. 158(1-2): p. 86-98.

(133) A method in accordance with embodiments of the present invention may further comprise applying a Fourier transform to the time segments, or to the extrapolated NMR signals, in both the k-space dimensions and the time domain dimension, such as to provide chemical species spectra as function of spatial location in the sample.

(134) A method in accordance with embodiments of the present invention may further comprise obtaining an image of the sample representative of spatial variations in a chemical constitution of the sample and storing the image in a data storage memory. For example, the chemical species spectra as function of spatial location may determine, or may be used to determine, such image.

(135) In a second aspect, the present invention also relates to a control unit that comprises a controller configured for controlling a magnetic resonance scanner such as to perform the steps of a method 10 in accordance with embodiments of the present invention.

(136) In a third aspect, the present invention also relates to a magnetic resonance scanner system that comprises such control unit in accordance with embodiments of the present invention.

(137) Hereinbelow, theoretical considerations are described relating to embodiments of the present invention. These theoretical considerations are provided to aid the skilled person in understanding the invention and reducing the invention to practice, and should not be construed as limiting the invention in any way.

(138) Methods in accordance with embodiments of the present invention may be based on a sub-division of an NMR signal, e.g. a MRS echo signal, that may be collected from a defined excited region using any suitable MRS sequence, e.g. a PRESS or STEAM sequence. This signal is divided into a predetermined number of time segments, e.g. sub-signals, such as sub-echoes, by inserting blipped phase encoding gradients, e.g. applied by a gradient corresponding to a predetermined spatial direction, for example in a direction x. The inserted blipped gradients may be applied at a constant time interval, and may have a constant amplitude, such that the phase encoding principle is satisfied, thereby linearly spanning the k-space:

(139) $k_x\left(t\right)=\frac{}{2}{}_0^{T_G}{G}_x\left(t\right).Math.\mathrm{dt}$

(140) where G.sub.x, is the blipped gradient amplitude applied along the x axis during the gradient duration time T.sub.G. Therefore, each acquired sub-signal may be a simultaneous collection of one k-space data point and the chemical shift information, e.g. a collection of spatial and spectral information at the same time, which is described by:
s(k.sub.x,t.sub.j)=s(x,t.sub.j).Math.e.sup.(i.2..k.sup.x.sup..x).Math.dx
where t.sub.j set to [(j1),j] determines the time occurrence of the sub-acquisition window (=T/PE) for the j.sup.th sub-signal, with PE the number of phase encodings (and hence the number of voxels) and T the time duration of the whole acquisition window. The variable j runs from 1 to (PE-1) (1 to PE). This equation can be approximated as a discrete sum:

(141) $s\left(k_x,t_j\right)=\underset{n=0}{\overset{\mathrm{PE}-1}{.Math.}}\stackrel{\_}{s}\left(n.Math.x,t_j\right).Math.e^{(i.Math.2.Math..Math.k_x.Math.\left(n.Math.x\right)}.Math.\mathrm{dx}$

(142) The sub-signal s may be described by:

(143) $\stackrel{\_}{s}\left(x,t_j\right)=\underset{m=1}{\overset{M}{.Math.}}{A}_{m.x}.Math.e^{-\left(d_{m.Math.x.Math.}{t}_j+2{\mathrm{if}}_{m.Math.x}{t}_{m.Math.x}\right)}$
where M is the number of metabolites, with A, d, and f describing amplitude, damping factor and frequency for each metabolite m.

(144) The k-space line can thus be filled by collecting PE equidistant data points s(k.sub.x,t.sub.j). This is achieved by linearly moving the k.sub.x variable in the encoded direction by increasing the accumulated blipped gradient amplitude, whereas the chemical shift information can be stored during the acquisition of the sub-signals. A Fourier Transform may then be performed on both the spatial and spectral domain to obtain the spatial distribution of the metabolite peaks.

(145) Spatial and spectral resolution may be limited due to a small number of data points that can be collected in both spatial and spectral domain. The spatial resolution may be set by the number of phase encodings PE, whereas the spectral resolution may be given by the sub-signal time-length set to N/PE, with N the total number of points of the collected MRS signal. In order to increase the spectral resolution, an auto-regressive approach may be used, in accordance with embodiments of the present invention, e.g. in a forward and/or a backward time sense, on each collected sub-signal to extrapolate it, e.g. to the total number N. As a consequence, a good SNR and spectral resolution can be achieved. As for the spatial direction a simple zero filling and a filtering procedure may be applied to reduce pixel contamination.

(146) However, k-space coverage is not necessarily restricted to a linear span. Other trajectories in k-space may be applied, e.g. by varying the blip amplitudes, the blip durations, varying the time interval in between blips, and/or varying the gradient direction in which the blips are applied. For example, in order to preserve the SNR, more data points may be collected around the centre of k-space and less points near the edges of k-space.

(147) To show the performance of a method in accordance with embodiments of the present invention, corresponding to the PRESS-based sequence illustrated in FIG. 3, a number of PE voxels containing each a noiseless proton brain spectrum containing the three main metabolites, namely, Creatine, Choline and NAA, were simulated in the example provided hereinbelow.

(148) The amplitude of the peaks of interests are varied in each voxel to duplicate real case spectra of healthy and brain tumor tissues. In a second set, noise with random distribution as well as macromolecules peaks describing baseline distortion are added into these spectra. The method was simulated in a Matlab code and tested against the reference spectra referred to hereinabove to recover the correct spatial distribution along with the estimation of the metabolite peaks.

(149) In this example, the echo time length N was fixed, e.g. to 512 points, then the number of voxels PE was varied, respectively for values 1, 4, 8 and 16, thus setting the time-length of each sub echo to N/PE. At each PE value the encoded sub-echoes were extrapolated to N/2 using a linear regression method and zero filled to the full length N Fourier transform is then used in both space and frequency domain to recover the spatial distribution of the spectra. The metabolite peaks of the spectra are fitted to a Lorentzian model to estimate their amplitudes.

(150) FIG. 11 and FIG. 12 show four recovered spectra without and with noise, corresponding to PE=4. FIG. 11 shows the noiseless original simulated spectra, e.g. the gold standard spectra, as the solid line of each pair, and the corresponding reconstructed spectra as the dashed line of each pair. FIG. 12 shows noisy original simulated spectra, with a SNR=20 dB, as the solid line of each pair, and the corresponding reconstructed spectra as the dashed line of each pair.

(151) As can be observed, the method in accordance with embodiments can accurately recover the correct spatial position of each spectrum, regardless of whether it belongs to a healthy or tumoral simulated tissue. The amplitudes of the peaks namely, creatine, choline and NAA are calculated using a spectral fitting method and compared to the expected ones. The results are shown in FIG. 13 to FIG. 20.

(152) Expected and estimated peak amplitudes of the three main metabolites, creatine, choline and NAA, are shown in these figures. FIG. 13, FIG. 15, FIG. 17 and FIG. 19 correspond to noiseless spectra, while FIG. 14, FIG. 16, FIG. 18 and FIG. 20 correspond to noisy data. FIG. 13 and FIG. 14 correspond to data for PE=4, FIG. 15 and FIG. 16 correspond to data for PE=8, FIG. 17 and FIG. 18 correspond to data for PE=16 and FIG. 19 and FIG. 20 correspond to data for PE=32

(153) The echo acquisition parameters for the simulations performed in this example were: an ADC bandwidth of 2500 Hz, an acquisition time of 204 ms, 512 sample points collected per signal, and one echo acquired per sequence.

(154) These exemplary results demonstrate that the method is able to preserves not only the correct spatial position of the spectra but also their amplitudes values for PE values ranging from 4 to 32. The means and standard deviations of the difference in metabolite peak amplitudes in percentage between the estimated and expected ones calculated as follows:

(155) 0$\mathrm{Diff}=\left(\frac{A_{\mathrm{ex}}-A_{\mathrm{es}}}{A_{\mathrm{ex}}}\right)\%$

(156) Where A.sub.ex and A.sub.es are the expected and estimated metabolite peak amplitude values respectively. These values are shown in the tables hereinbelow. For 4 encoding steps:

(157) TABLE-US-00001 Without noise With noise (20 dB) Peak Mean Standard deviation Mean Standard deviation NAA 0.0004 0.0001 0.08 0.09 Choline 0.005 0.003 0.10 0.04 Creatine 0.018 0.003 0.11 0.07 For 8 encoding steps:

(158) TABLE-US-00002 Without noise With noise (20 dB) Peak Mean Standard deviation Mean Standard deviation NAA 0.004 0.002 0.10 0.088 Choline 0.006 0.003 0.106 0.047 Creatine 0.02 0.007 0.11 0.07 For 16 encoding steps:

(159) TABLE-US-00003 Without noise With noise (20 dB) Peak Mean Standard deviation Mean Standard deviation NAA 0.004 0.004 0.129 0.10 Choline 0.017 0.006 0.10 0.04 Creatine 0.045 0.008 0.11 0.10 For 32 encoding steps:

(160) TABLE-US-00004 Without noise With noise (20 dB) Peak Mean Standard deviation Mean Standard deviation NAA 0.02 0.015 0.17 0.19 Choline 0.033 0.012 0.166 0.14 Creatine 0.05 0.013 0.30 0.15

(161) It can be observed that the data accuracy may slightly drop at high phase encoding steps, which may shorten the sub-echo time duration, and consequently broaden the shapes of the peaks in the frequency. The auto-regression methods described hereinabove may be used to extrapolate the sub-echoes and expand them, e.g. to the time-length of a half echo.

(162) To show the performance of another method in accordance with embodiments of the present invention, corresponding to the sequence illustrated in FIG. 4 and the processing illustrated in FIG. 9, the sequence, for convenience referred to as RAPID-SI, configured for an acceleration in the anterior-posterior (AP) direction, was implemented on a Siemens Tim-Trio 3T MRI scanner.

(163) The acquired NMR data is processed offline using an in-house Matlab software, in which a HLSVD technique is applied on a full MRS reference signal, e.g. without application of the blipped gradients, e.g. at the centre of k-space (k.sub.j=0). The reference signal is acquired prior to the encoding using blipped phase gradient pulses in accordance with embodiments of the present invention. The signal parameters d.sub.m, and f.sub.m are estimated from the reference signal, and the sub-bases .sub.sub(m,t.sub.j) are determined, e.g. as illustrated by the workflow in FIG. 9. These sub-bases are then used to estimate the complex amplitudes

(164) $A_{m_{k_j}}^{{}_{k_j}}$
from the collected sub-signals. A Fourier Transform is performed on both spatial and time domains to obtain the spatial distribution of the spectra. A Hamming filter is applied in the spatial domain to reduce voxel contamination.

(165) To accurately estimate signal parameters using the HLSVD approach, the duration of the acquisition window is set to the average T*.sub.2 values of the main brain metabolites (NAA, Creatine, Choline, Myo-inositol, and lactates) obtained at 3 Tesla, which was found to be approximately 250 ms, and given the receiver bandwidth set to 1600 Hz, the adequate total number of points N.sub.full for the highest R possible, was found to be 512. This is confirmed by 1D RAPID-SI tests performed on phantoms using different N.sub.full values described further hereinbelow.

(166) Tubes were filled with brain metabolites solutions (Cho, NAA and lactate) having known concentrations and are placed inside a cylindrical water solution phantom, used to acquire 1D and 2D RAPID-SI and CSI data sets using the following acquisition parameters: TE=35 ms, TR=2000 ms, BW=1600 Hz, N=512, N.sub.avg=2. The FOV dimensions are AP=100 mm, HF=10 mm, whereas RL is set to 55 mm for 1D and to 100 mm for 2D acquisitions. The spatial resolution is varied from 14, 18 and 116 for the 1D data sets and from 44, 88, to 1616, for the 2D data sets. The acceleration factor ranges from 4 to 8 depending the data size.

(167) To check the accuracy on the RAPID-SI data, absolute quantification results from the voxels covering the tubes of the metabolites are calculated using the equation below:

(168) $C_{\mathrm{metabolite}}=C_{\mathrm{water}}.\left(\frac{2}{\mathrm{np}}\right).\left(\frac{A_{\mathrm{metabolite}}}{A_{\mathrm{water}}}\right).\frac{\exp \left(\frac{\mathrm{TE}}{T{2}_{\mathrm{metabolite}}}\right)}{\exp \left(\frac{\mathrm{TE}}{T{2}_{\mathrm{water}}}\right)}.\frac{\left\{1-\exp \left(\frac{-\mathrm{TR}}{T{1}_{\mathrm{water}}}\right)\right\}}{\left\{1-\exp \left(\frac{-\mathrm{TR}}{T{1}_{\mathrm{metabolite}}}\right)\right\}}$

(169) and compared to those obtained from the full reference spectrum which contains metabolite peaks signals from all the tubes of the phantom. The variable C.sub.metabolite is the metabolite concentration (mM), C.sub.water is the water concentration used as internal reference (110 M), A.sub.metabolite is the area under the metabolite peak calculated using the HLSVD approach, A.sub.water is the area under the water peak calculated from the unsupressed water data, np is the number of protons in a given metabolite, and T1.sub.metabolite, T2.sub.metabolite, T1.sub.water, T2.sub.water are the relaxation times of metabolite and water respectively. Metabolite maps for both RAPID-SI and CSI methods are calculated, and voxel bleed are calculated and compared.

(170) From the obtained 1D RAPID-SI spectra, acquired with different N.sub.full values at different spatial resolutions, N.sub.full set to 512 is found to be a good compromise between data acceleration, SNR, and accurate signal parameter estimation. With this setting, the highest acceleration factor R reached the value of 8, while high spectral resolution, e.g. better separation between the lactate doublet, and SNR values are obtained.

(171) Compared to CSI, the RAPID-SI reduces the acquisition time by R, which is ranging from 4 to 8, given the spatial resolution varying from 4 to 16. A 1616 data acquired in around 17 min with CSI, given the setting above, is reduced to about 2 min with the RAPID-SI. The table hereinbelow summarizes acquisition times for both techniques at different spatial resolutions.

(172) TABLE-US-00005 Acquisition (1 4) (1 8) (1 16) (4 4) (8 8) (16 16) time R = 4 R = 8 R = 8 R = 4 R = 8 R = 8 CSI 16 s 32 s 64 s 64 s 256 s 1024 s RAPID 4 s 4 s 4 s 16 s 32 s 128 s

(173) In terms of data quantification, RAPID-SI provides accurate results. The table hereinbelow shows the estimated concentration values vs. the expected ones from the following phantom data sets.

(174) TABLE-US-00006 Lac Concentration (contamination per (mM) NAA voxel %) Cho Expected 25 17 7 CSI (4 4) 26 17.88 (0.2%) 6.53 RAPID (4 4) 26.2 18 (0.4%) 7.68 CSI (8 8) 23 17 (0.02%) 6 RAPID (8 8) 27 17.14 (0.04%) 7.20 CSI (16 16) 22 16 (0.2%) 7.8 RAPID (16 16) 22 16.70 (1.2%) 8

(175) A qualitative comparison of the 1D and 2D RAPID-SI results displayed in FIGS. 21 to 24 shows a very good match to the corresponding CSI results. However, from the figures above, one can notice that voxel contamination in RAPID-SI is slightly higher than CSI, especially for the lactate doublet. This contamination is 0.2% per voxel in the 1616 results for the CSI and increases to 1.2% per voxel for the RAPID-SI. This is due to the overestimation of the

(176) $A_{m_{k_j}}^{{}_{k_j}}$
values of the lactate doublet, which appear as singlet in the collected sub-signals caused by the decrease in the spectral resolution. Therefore, a more reliable post-processing may be preferable to better estimate these complex amplitudes to reduce voxel contamination. Voxel contamination is calculated by dividing the amplitude of any metabolite detected in any voxel by its amplitude in the expected voxel.

(177) As expected, the overall SNR is lower in RAPID-SI than CSI by approximately the factor R in the 44, 88, and 1616, results when processing reduced size sub-signals only, e.g. see the table hereinbelow. These values, however, are increased in the fully reconstructed data as shown in FIG. 22 and FIG. 24. The SNR is calculated as the ratio of the NAA peak intensity and the standard deviation of noise. It is to be noted that the encoding efficiency of the metabolite signals is subject to their T*.sub.2 values. Given that the whole k space is filled within the duration of the acquisition window, the number of k space data points used to spatially encode each metabolite signal is proportional to its corresponding T*.sub.2 value. Signal components with long T*.sub.2 are sampled by several data points in k space, whereas those with short T*.sub.2 have less points to characterize them. This may cause a spatial blurring of the of metabolite maps with short T*.sub.2. To overcome this problem, a deconvolution between sinc-like functions with width at half maximum (WHM) inversely proportional to the T*.sub.2 values of the blurred metabolite signals and their corresponding maps is performed to recover the original spatial distribution of the metabolites.

(178) TABLE-US-00007 SNR (4 4) (8 8) (16 16) R = 4 R = 8 R = 8 CSI 32 35.6 36 RAPID 15 12.8 12.6

(179) While present examples relate to tests performed with Cartesian k-space trajectories, embodiments of the present invention are not limited thereto, e.g. other k-space paths may be implemented in a method in accordance with embodiments of the present invention. A method in accordance with embodiments can advantageously be integrated into many localized MRS sequences as known in the art that generate an FID or echo signal, e.g. PRESS or STEAM. Furthermore, a method in accordance with embodiments of the present invention can be combined with other acceleration methods, such as parallel imaging, compressed sensing and/or the complex conjugate principle.

(180) In conclusion, a method in accordance with embodiments of the present invention may use a partition of the acquired MRS signal into encoded sub-signals by inserting blipped gradients in a predetermined, e.g. selected, direction. A method in accordance with embodiments may provide MRSI results with a good balance between SNR, spectral resolution and acquisition time. For example, compared to CSI, a RAPID-SI approach may achieve a significant time reduction in acquiring accurate metabolite maps. Even though an SNR decrease may be implied, this can be compensated by a signal processing in accordance with embodiments of the present invention to reduce the noise contribution and to estimate the signal parameters.