WO2015110909A1 - Motion correction to improve t10 estimation in dce-mri - Google Patents

Abstract

The present invention relates to the correction of the TlO-parameter, which is used in DCE-MRI-based pharmacokinetic analysis, for subject motion occurring between the acquisition of images with different flip angles or repetition times from which T10 is derived. Groupwise registration of these images is based on a new similarity metric that quantifies the T10 regression error.

Description

MOTION CORRECTION TO IMPROVE T10 ESTIMATION IN DCE-MRI

CROSS-REFERENCE TO RELATED APPLICATION

[0001] This application claims priority to co-pending U.S. Provisional Application entitled "MOTION CORRECTION TO IMPROVE T₁₀ ESTIMATION" having Serial No. 61/930,787, filed on January 23, 2014, which is incorporated herein by reference.

TECHNICAL FIELD

[0002] The present disclosure generally relates to magnetic resonance (MR) image sequences, in particular variable flip angle and variable pulse repetition time image sequences, and more particularly to motion correction between images in MR image sequence to provide improved T_w estimation.

BACKGROUND

[0003] Pharamacokinetic (PK) analysis of tumours is an effective tool to distinguish between malignant and benign neoplasms. It has shown potential to assess and predict anti-angiogenic treatment response. Microvasculature can be characterized using models for contrast agent perfusion, which allows for differentiation between "normal healthy vessels" and vasculature grown due angiogenesis [1 ]. Furthermore, the evolution of these vascular characteristics over time has been successfully used to monitor tumour growth and evaluate how tumours will respond to treatment. Treatment outcome prediction allows for the omission of stages that would not be effective, better regulation of dosages and adjustments for specific cancer characteristics (such as using radiation sensitizers) [2].

[0004] To perform pharmacokinetic analysis, imaging of contrast agent uptake within human tissue is performed. Dynamic Contrast Enhanced MRI (DCE-MRI) is considered one of the most suitable techniques for such task. DCE-MRI has many advantages over other dynamic imaging techniques: it is more readily available in clinics and hospitals (compared to dynamic PET), has high signal/contrast-to-noise-ratio and high spatial and temporal resolution, and more importantly, does not expose the patient to radiation (compared to CT perfusion) [3].

[0005] DCE-MRI consists of taking several sequential magnetic resonance images over a short period of time upon the injection of a contrast agent, which locally enhances the measured intensities according to its concentration. These can then be analyzed through pharmacokinetic models yielding quantitative information on the tissue microvascular properties. The most frequently used model, the Tofts model [4], provides relevant tissue information through two permeability parameters, K_trans and k_ep.

[0006] Correlating intensity changes on DCE-MRI to agent concentration is one of the obstacles found on this type of analysis. Unlike CT perfusion, agent concentration is not linearly related to the image intensity. Locally estimating the agent concentration also requires the estimation of the relaxation time (Γ₁₀), a tissue property, for the same region. Additional volumes of images need to be acquired to estimate this property for each voxel of the DCE-MRI volumes.

Several protocols have been developed for this task, such as inversion recovery

(IR) and saturation-recovery (SR) [5]. In this work, T₁₀ estimation is performed by the variable nutation angle method, which consists of acquiring a sequence of spoiled gradient echo (SPGR) volumes over several flip angles with fixed pulse repetition time [5]. The tissue relaxation time T_w at each voxel may then be estimated by analyzing the MRI intensities at this same voxel on all the different flip angle volumes. Alternatively, this estimation can be done by acquiring image sequences with variable pulse repetition time T_R and a fixed flip angle. The intensities over these sequences can also be fitted to estimate 7₁₀in a very similar approach.

[0007] Patient motion within these volumes, however, may cause errors in the T₁₀ estimation, which may then hinder the subsequent accuracy of the pharmacokinetic analysis. The relation between T₁₀ estimation errors on K_trans and k_ep was studied in [8], showing that inaccuracies in T₁₀ resulted in even greater errors in K_trans . Even though a lot of effort has been put into analyzing and correcting for motion within the dynamic images of a DCE-MRI acquisition [6] [7], very little work has been carried out in correcting motion within the variable flip angle or pulse repetition time sequence. Moreover, a nonlinear image registration method using B-Spline free-form deformations with mutual information was also applied in that work to correct for motion [8], but the actual effect and efficacy of this procedure were not evaluated.

[0008] In current clinical and research practice, no motion between flip angle images is assumed, which leads to an erroneous T₁₀ (model input) parameter estimation, an error which in turn propagates into the estimation of the clinical PK parameters of interest. [0009] Accordingly, there is a need to address the aforementioned deficiencies and inadequacies.

SUMMARY

[0010] The present systems and methods can be applied to correct for subject motion between images within a sequence of magnetic resonance (MR) images having different acquisition parameters. The motion correction can be applied to improve T₁₀ estimation. The image sequence can be a flip angle or pulse repetition time sequence involving the acquisition of three or more images (or volumes) with different flip angles or pulse repetition times. The sequence may be a spoiled gradient (also called a spoiled gradient recall) (SPGR), echo pulse sequence, or other pulse sequence such as a gradient echo pulse sequence or a spin echo pulse sequence. Our systems and methods correct for this motion, thus providing a more reliable model input parameter T₁₀, which in turn can increase the accuracy and reliability of, for example, clinical PK parameter estimation based upon T₁₀ estimation. In various aspects, formulation of a T₁₀ regression based error metric for motion correction of variable flip-angle or pulse repetition time MR images is provided, leading to improved accuracy of T₁₀ tissue maps.

[0011] Patient-specific T_w tissue maps are normally estimated using a so- called variable flip-angle or variable pulse repetition time technique, using linear regression about several flip-angle or pulse repetition time MR images. In various aspects we use this linear regression to correct for motion between the various flip-angle or pulse repetition time images (i.e., within the image sequence), for example by simultaneously minimizing the regression error and aligning these images. The alignment can be done in a groupwise manner. The output of this method can be a motion-corrected T₁₀ map, which can then form the input for standard PK analysis. The motion-corrected T₁₀ map can also be used for other purposes, such as tissue characterization, or image segmentation.

[0012] In an embodiment a method for motion correction for T₁₀ estimation in a sequence of magnetic resonance (MR) imaging is provided comprising the steps of: a) acquiring a sequence of MR images with varying acquisition parameters; b) applying a transformation model to at least one of the acquired images to produce a transformed image; c) using the transformed image to estimate T₁₀; d) generating a model of the transformed image using the estimated T₁₀; e) applying an error metric to determine the error between at least one acquired image and the modeled image; f) updating the transformation model based on the error determined; g) applying the updated transformation model to produce a new transformation of at least one acquired image; and h) using the new transformation of the at least one acquired image to provide a new T₁₀ estimation.

[0013] In an embodiment a system is provided, comprising: at least one computing device; and programming logic executable in the at least one computing device, the programming logic comprising logic that: a) acquires a sequence of MR images with varying acquisition parameters; b) applies a transformation model to at least one of the acquired images to produce a transformed image; c) uses the transformed image to estimate T₁₀; d) generates a model of the transformed image using the estimated T₁₀; e) applies an error metric to determine the error between at least one acquired image and the modeled image; f) updates the transformation model based on the error determined; g) applies the updated transformation model to produce a new transformation of at least one acquired image; and h) uses the new transformation of the at least one acquired image to provide a new T₁₀ estimation.

[0014] In any one or more of the aforementioned embodiments sub-parts (d)-(h) can be repeated until a desired minimized error has been determined. The logic in part (a) can acquire a variable flip angle or a pulse repetition time sequence of MR images. The transformation model can be selected from the group consisting of linear (rigid-body or affine) or non-linear (parametric or non- parametric) transformations, or combinations thereof. The error metric can be selected from the group consisting of sum of squared residuals, sum of mean residuals, maximum residual mean square or any other suitable error metric. The transformation model can include a set of parameters, and the transformation model can be updated by adjusting one or more of the parameters using an optimizer. The optimizer can be selected from the group consisting of a Levenberg-Marquardt, Gradient Descent, Hill Climbing, quasi- Newton, Conjugate Gradients, and Simulated Annealing optimizers. The acquisition parameter can be varied. The acquisition parameter can be either the flip angle or the pulse repetition time, and the sequence can be acquired using a spoiled gradient recall (SPGR), echo pulse sequence, a gradient echo pulse sequence or a spin echo pulse sequence. The new T₁₀ estimation can be used to provide a motion corrected T₁₀ map. The acquired MR images can be 2D or 3D images. The acquired images can be of tissue of a subject. The final T₁₀ estimate can be used for pharmacokinetic analysis, tissue characterization and/or image segmentation. The image sequence can include three or more MR images and the transformation model is applied groupwise to the images.

[0015] Other systems, methods, features, and advantages of the present disclosure for motion correction within a sequence of various flip angle MR images will be or become apparent to one with skill in the art upon examination of the following drawings and detailed description. It is intended that all such additional systems, methods, features, and advantages be included within this description, be within the scope of the present disclosure, and be protected by the accompanying claims.

BRIEF DESCRIPTION OF THE DRAWINGS

[0016] Many aspects of the disclosure can be better understood with reference to the following drawings. The components in the drawings are not necessarily to scale, emphasis instead being placed upon clearly illustrating the principles of the present disclosure. Moreover, in the drawings, like reference numerals designate corresponding parts throughout the several views.

[0017] FIG. 1 depicts a T₁₀ regression based on the intensity of a variable flip angle sequences in which the regression does not fit perfectly the data points.

[0018] FIG. 2 depicts a flow chart of an exemplary process using the T₁₀ similarity metric of the present disclosure.

[0019] FIG. 3 depicts an example of the process of contrast agent leaking from the blood plasma into the extravascular extracellular space (EES) between the cells [8]. [0020] FIG. 4 depicts a two compartment model for the tracer kinetic behaviour in tumours in which the contrast agent enters the plasma compartment (vp) through the plasma flow (Fp) and a proportion of this tracer is extracted at a rate PS to the extravascular extracellular space (EES) volume (ve) and at the same rate the tracer flows back into vp.

[0021 ] FIG. 5 depicts sample agent concentration curves for the Weinmann AIF (top curve) and the Tofts Model using the Weinmann AIF (bottom curve).

[0022] FIG. 6 depicts exemplary MRI acquisitions from a spoiled gradient echo (SPGR) sequence with variable flip angles: (a) 30°; (b) 90°; and (c) 120°.

[0023] FIG. 7A is an exemplary flow chart depicting pharmacokinetic model parameters estimation using a signal enhancement curve from a DCE- MRI sequence in which an observed signal enhanced curve (S) is compared with an artificial (estimated) signal curve (S_est) generated using the Tofts model and an Arterial Input Function (AIF) model and the mean square error between the two curves is minimized as a function of the parameters being estimated (K_trans and k_ep ).

[0024] FIG. 7B depicts an example of use of the T₁₀ estimation depicted in FIG. 2 to estimate pharmacokinetic parameters.

[0025] FIG. 8A illustrates an example of a groupwise registration of a series of variable flip angle volumes by simultaneously registering the whole set of volumes.

[0026] FIG. 8B depicts a sample of a B-Spline free-form deformation transform: original image with a deformation free control mesh (left), transformed image with the respective control mesh (right). [0027] FIG. 9 is a flow chart depicting an overview of a method used to estimate T-io error due to non-linear motion in which ground truth Tw and Mw maps are used to create variable flip angle volume sequences with non-linear motion and these sequences are then used to estimate Ti₀ and M_w maps, which are then compared by the mean absolute error with the ground truth maps.

[0028] FIG. 10 is a flow chart depicting a procedure performed to evaluate similarity metrics for groupwise rigid registration of variable flip angle SPGR sequences in which ground truth T₁₀ and ₁₀ maps are used to create variable flip angle volume sequences with rigid motion and these sequences are used to correct for motion using a groupwise rigid registration framework and the post registration mean displacement error is then computed, indicating how well the registration performed.

[0029] FIG. 1 1 depicts estimated T₁₀ error due to motion within the SPGR sequence volumes.

[0030] FIG. 12 depicts K_trans and k_ep estimation error due to T₁₀ inaccuracies showing a magnification of the error occurs when using underestimated T_w values.

[0031] FIG. 13 depicts an example of incorrect pharmacokinetic model parameter estimation due to - 20% inaccuracy in T₁₀, showing that the computed signal enhancement (SE) curve (black) does not fit the observed signal enhancement curve (lower curve) and an overestimated K_trans parameter generates a matching signal enhancement curve.

[0032] FIG. 14 depicts K_trans estimation error due to the absolute T_w error in which the overall trend seen in this result was of a slight increase of the error from 7^~ ₀ being propagated to the estimated K_trans . [0033] FIG. 15 is a schematic block diagram of an apparatus in which embodiments for subject motion correction between images in a flip angle image sequence can be implemented.

DETAILED DESCRIPTION

[0034] Having summarized various aspects of the present disclosure, reference will now be made in detail to the description of the disclosure as illustrated in the drawings. While the disclosure will be described in connection with these drawings, there is no intent to limit it to the embodiment or embodiments disclosed herein. On the contrary, the intent is to cover all alternatives, modifications and equivalents included within the spirit and scope of the disclosure as defined by the appended claims.

[0035] Described in detail below are systems and methods for subject motion correction within a sequence of magnetic resonance (MR) images having different acquisition parameters, for example a sequence of flip angle or pulse repetition time magnetic resonance (MR) images, for T₁₀ estimation. In various aspects, the flip angle or pulse repetition time sequence involves the acquisition of three or more images (or image volumes) with different flip angles or pulse repetition time, for example a variable flip angle image sequence. The image sequence is preferably acquired using a spoiled gradient recall (or sometimes simply called a spoiled gradient) (SPGR) pulse sequence. Other pulse sequences may be used, however, such as gradient pulse sequence or a spin- echo pulse sequence. The motion correction within the image sequence can be used for T₁₀ estimation and subsequently for pharmacokinetic analysis, for example in Dynamic Contrast Enhanced-MRI (DCE-MRI), tissue characterization or image segmentation.

[0036] In various aspects we provide systems and methods for correcting for subject motion in variable flip angle or pulse repetition time sequences used for T₁₀ estimation. In one or more aspects, the image sequences can include a sequence of flip angle volumes or pulse repetition time taken of tissue of interest in the subject. Misalignment due to subject motion between two or more of the volumes of the image sequence can be identified. The amount of misalignment can be quantified. Systems and methods can then be applied to correct for the misalignment within the image sequence providing improved T₁₀ estimation which may subsequently be used, as noted above.

[0037] The goal in image registration is to align two or more images. In the present case, the images (which can be 2D or 3D images), are from the same subject, for example a patient, acquired during the same study, but have different acquisition parameters (for example with variable, or differing, flip angles or pulse repetition times). In various aspects, to find the best alignment it is preferable to evaluate and compare different possible solutions. A metric can be applied as a function that evaluates how good a solution is. For example, the metric can be an error metric and by minimizing it an optimal alignment may be found. Conversely, the metric can be a similarity metric and by maximizing the similarity an optimal alignment may be found. For simplicity herein we will refer to both types of metrics as an "error" metric. We present a new error metric suitable for MR image sequences described herein (explained below). Non- limiting examples of suitable error metrics include sum of squared residuals, sum of mean residuals, maximum residual, mean square or any other suitable error metric.

[0038] The transformation model defines how the images (or volumes) can be manipulated to be aligned. Various transformation models can be used, such as linear (rigid or affine), or non-linear/deformable models (parametric models such as B-spline Free-Form Deformations, or non-parametric/dense models such as optical flow or Demons variants), as well as combinations of global and local motion models. The most simple concept is translation, where the images are just moved along the x-, y- or z-axis in Cartesian coordinates. This, however, can be extended to rigid-body transformations, in which images not only translate but also rotate. More complex (non-linear) transformation models allow different regions of the image to move independently, deforming it locally.

[0039] A transformation consists of a transformation model and associated degrees of freedom. For example, a translation transformation would consist of a translation value for each of the 3 axes (x, y and z). The transformation model used in our examples below was rigid, but it is important to understand that the transformation model is quite independent to the error metric. Thus, our error metric could be used with any transformation model. They can be different, interchangeable components of a larger system. For example, a non-linear or deformable transformation model can be used.

[0040] Hence, in one or more aspects, given a set of images (/) being registered, a transformation model (T ) and an error metric, here a similarity metric (SF ), for each set of transformation parameters or degrees of freedom

(Θ), the quality of the registration may be evaluated as SF (T (I , Θ)), where T (I , P ) consists of applying the transformation T with parameters Θ to the images / . To obtain an optimal registration according to a given similarity metric, we then find the set of parameters Θ* which minimizes SF. The last component of a registration framework is the optimizer, which will iteratively test solutions P searching for the optimum similarity.

[0041] A suitable optimizer is a Levenberg-Marquardt optimizer. But as the transformation model this is another interchangeable part of the alignment or registration system. Other optimizers may be used, such as e.g. Gradient Descent, Hill Climbing, quasi-Newton, Conjugate Gradients, and Simulated Annealing.

[0042] As a non-limiting example, the T₁₀ regression error metric can be designed for the problem of estimating T₁₀ over sequences of MR images acquired using different acquisition parameters, for example variable flip angle sequences. These sequences consist of three to four (or more) MRI acquisitions from the same subject. In the case of variable flip angle sequences, they are taken with different flip angles (a), an MR acquisition parameter. While in variable repetition time sequences, the varying parameter is TR, another acquisition parameter. The intensities (S) at each of pixel (/^') of these images (volumes) within a sequence are related by the following equation:

where T_R is an known acquisition parameter. By fitting the intensity at a given pixel /^' for all the flip angles or pulse repetition times images in the above equation it is possible to estimate T_w and M₀. However, patient motion between the different MRI acquisitions may affect the Tw estimation as the correspondence between the pixels in the different flip angle and pulse repetition time volumes would be compromised. In this case, Equation 1 does not fit perfectly and will show a regression error. FIG. 1 shows an example of how the data points from the intensities from variable flip angle images might not be perfectly fitted using Equation 1 due to motion.

[0043] The concept behind the Tw regression similarity metric is that, as misalignment causes a misfit between 1 and the data points, a regression error can be used as a measure of misalignment. Thus, minimizing the regression error will simultaneously correct for motion and provide a better T_w estimation. This similarity metric is defined, for a variable flip angle sequence, as:

= I 8{i, (x(n))— 8(L ain), X½ , Mo ) ]

⁵ (2a) where N is the number of acquired flip angles, S is the MRI intensity at each of the variable flip angle volumes and S is the predicted intensity by the model fit.

[0044] For a variable pulse repetition time sequence, it is defined as,:

»===! (2b) where S now is the MRI intensity at each of the variable pulse repetition volumes and S is the predicted intensity by the model fit.

[0045] Referring to FIG. 2, an exemplary registration process is depicted, in which at each iteration (k) of the registration process 200, the following steps can be taken: 1. A sequence of MR images, for example aflip angle image sequence (such as a variable flip angle sequence) or pulse repetition time image sequence is acquired 203.

2. The transformation model T is applied 206 to the images (volumes) / with a set of parameters 9_K, generating / , a transformed set of volumes.

3. T- _O and _Q maps are estimated by fitting 209 Equation 1 to the data points from / .

4. Using ¾ and _Q a modeled set of volumes ¾ are generated (using Equation 1 ). ¾ represents the variable flip angle with a perfect fit 212 to the estimated T_W and M₀ parameter maps at the current iteration.

5. The error is determined. As an example, a similarity metric can be used and the similarity computed 215 by computing the mean square error between the transformed data volume (/£) and l_j . This is equivalent to computing the regression error during the model fit in step 2.

6. The transformation model is updated 218, for example by updating the model's parameters, according to an optimization algorithm, for the next iteration.

7. Once the optimizer has minimized the regression error, the 7i₀ estimation at that point is recorded 221 .

[0046] While the above example is applied to the set of the acquired images (volumes), it can be applied to one or more, or all, of the images (volumes) acquired.

[0047] We now provide a discussion of the background to our present systems and methods for motion correction within a sequence of flip angle or pulse repetition time magnetic resonance (MR) images, and a more detailed discussion of various non-limiting aspects of our present disclosure, including our present T₁₀ regression error similarity metric and uses of the T₁₀ estimation.

I. Background

A. Tofts Model

[0048] Applying pharmacokinetic models to describe data from dynamic imaging allows a biological understanding of tissue characteristics. In the case of tumour analysis, the functional dynamic images obtained (for example from DCE-MRI) can be translated by these models into relevant microvasculature information, which may then be used to better differentiate the type of tumour and the efficacy of the treatment being performed. The blood flow and permeability, which affect how the contrast agent penetrates and accumulates within the tissue, are particularly meaningful in oncology and are the target information in this type of analysis. FIG. 3 shows the process of contrast agent leaking (proportional to flow and permeability) from the blood plasma and accumulating in the space between cells. Pharmacokinetic modelling of contrast agent assumes that the agent concentration in the patient's tissue over time (ct (t)) relates to the injection of a contrast bolus according to:

where F_p is the plasma flow, R(t) the tissue characteristic time-course and c_a(t) the Arterial Input Function (AIF) (the concentration at the blood vessel feeding the tissue) [10]. [0049] The Tofts model is a two compartment model, namely the plasma volume (v_p) and the extravascular extracellular space (EES) volume (v_e) (FIG. 4) (contrast agent molecules are considered too large to penetrate the intracellular space) [10]. FIG. 4 depicts a model for the tracer kinetic behaviour in tumours [1 1 ] in which the contrast agent enters the plasma compartment (vp) through the plasma flow (Fp) and a proportion of this tracer is extracted at a rate PS to the extravascular extracellular space (EES) volume (ve) and at the same rate the tracer flows back into vp.

[0050] Thus, the contrast concentration in the tissue is given by the weighted sum of the concentrations in these two compartments (c_p and c_e):

' ^{■ '} - ^' (4)

[0051 ] The model also assumes that the contrast agent is exchanged between the arterial system and the plasma volume and also between the EES and the plasma volume. The contrast agent flows into the EES through leaky capillary walls and the permeability that drives this process is the information which may help characterize tumours.

[0052] This model further presupposes that the agent concentration accumulates in the EES and is constant and negligible in the plasma [1 1 ].

Equation 5 presents the observed concentration obtained under these assumptions:

Cf it) = Kt._r(,.r< * exp i—k^ ) *€„{t^'i

^{" ^ '} ' (5)

K_trans= FpE is the volume transfer coefficient, given by the product between the plasma flow and the proportion of contrast agent that extravasates the wall between the plasma and the EES. Thus, this parameter quantifies how leaky the micro- vessels are in this region and may be used to characterize tumours. K_ep= K_trans/v_e quantifies the rate that the contrast agent goes back from the EES to the plasma.

B. Weinmann Arterial Input Function

[0053] The Arterial Input Function (AIF) represents how the contrast agent is fed to the microvascular system being analysed in the pharmacokinetic model. It can be obtained by extracting blood from large supplying vessels over time and measuring the agent concentration, by processing images of the contrast uptake on the MRI volumes inside the supplying vessel or by approximation using known population based models. The Weinmann AIF is one of these models and is represented by a decaying bi-exponential function. The agent concentration profile in the Weinmann AIF is given by [8]:

= (At expf—m-jf} 4- A*t expi—m.it)) where and A2 are amplitude parameters, while m₁ and m₂ are rate constants, D is the bolus dosage. Applying it to the Tofts Model (Equation 5) yields:

C i i ) = (exp(— m\i j— expl— fcepl))

(7)

[0054] Other models can be used, however, such as the Orton AIF model.

C. Converting concentration to relation time (Ti )

[0055] The Tofts Model requires agent concentration uptake data to estimate _trans and k_ep. The presence of contrast agent causes a shortening of the 7^~i relaxation time which is reflected by an increase of the acquired MRI intensities. Equation 8 presents the relation between the agent concentration and the MRI intensity value [4].

(8)

where r_x is the relaxivity, a contrast agent property, and C the tissue concentration. From Equation 6 it is evident that the estimation of T₁₀, the pre- injection relaxation time, plays a vital role on obtaining physical information from DCE-MRI sequences.

[0056] An analogous dependency exists for r₂ ^* relaxation time, with another relaxivity constant r₂:

-L ₌ _L _{^ Γ}._ιΓ

- ½o ₍₉₎

[0057] FIG. 5 shows exemplary sample agent curves for the Weinmann AIF (upper curve) and the Tofts Model using the Weinmann AIF (lower curve) for the agent concentration uptake data.

D. T₁₀ estimation

[0058] The tissue relaxation time T₁₀ can be computed from a sequence of flip angle or pulse repetition time volumes, for example variable flip angle or pulse repetition time spoiled gradient echo (SPGR) volumes. Under the assumption of no motion, noise and changes on the magnetic properties of the tissue, the intensity of any voxel /^' in these volumes will change according to Equation 10 [5]:

where T_R and a are respectively the sequence pulse repetition time and the flip angle, both known and configurable parameters from the MRI scanner. ₀ is the equilibrium magnetization, which depends on the scanner configuration and material properties:

g is the MRI scanner gain, 7Έ the sequence echo time (both acquisition parameters) and the tissue proton density. In 7 DCE-MRI acquisitions (such as the ones used herein), the sequence echo time is short enough that the effects of concentration on 7₂ ^* are negligible to the observed MRI signal intensity [12].

[0059] Both tissue dependent parameters (7₁₀ and ₀) can be estimated with two or more volumes with variable flip angles a. This can be directly done by rearranging Equation 8 into a linear function (Equation 12):

Y = AX ÷ B

(12)

where X , Y , A and B are defined as:

(13)

(14)

[0060] At each voxel, A and B remain constant and can be estimated by a linear regression over the observed values on the variable flip angle volumes.

[0061] However, A and B can only be considered constant at a voxel i over all flip angles and pulse repetition times if T₁₀ and ₀ are the same at i for over all volumes. For this to be true, no motion may be present between the MRI acquisitions, as at a given position, there would be different parameters T₁₀ and ₀ related to different tissues on each MRI volume. As a consequence, the linear relation shown in Equation 12 does not hold true, a regression based on it does not provide a perfect fit and leads to inaccuracies in the T₁₀ estimation.

[0062] In a similar form, noise also adds errors to the dependency found in Equation 10 and hinders the estimation of the relaxation time.

[0063] FIG. 6 shows an example of three images from different flip angle acquisitions, at (a) 3 (b) 9 and (c) 12 degrees, from the same patient; despite all images showing the same tissues and morphology, the MRI intensities change according to the acquired flip angle.

E. Pharmacokinetic model parameters estimation

[0064] In Sections l-A to l-D above steps that can link acquired DCE-MRI data (both the dynamic and the variable flip angle SPGR volumes) to contrast agent concentration and pharmacokinetic modeling are presented. Combining these tools, these MRI volumes can be used to estimate the desired pharmacokinetic parameters from the Tofts model (K_trans and k_ep). [0065] The DCE-M RI sequence provides, at each point /^', intensity curves S, caused by the contrast agent uptake in the imaged region. It may be expressed as S(i, t(n)), as it is an intensity curve dependent on the position /^' and acquisition n, given at the time t(n).

[0066] Similarly, an estimated signal curve S(i, t(n), Q>, k_trans, k_ep) can be computed for any point /^', using the pharmacokinetic model, the AIF and the conversion between concentration to signal intensity (Equations 5, 6 and 8). The known parameters are expressed by Φ, and include acquisition settings, estimated T₁₀ and all of the model parameters (with exception of K_trans and k_ep).

[0067] Thus, a pair of tf_trans and k_ep can be evaluated as properly describing the tissue properties by analyzing how well does the estimated signal curve, S, generated by it approaches the observed curve, S. This can be accomplished by computing the mean square error of both of these signals: errof¾ = Y^*' (s(i,i{n)} - 8{i_t n. Φ,

(16)

[0068] The phamacokinetic (PK) parameters can then be estimated by finding the pair K_trans and k_ep which minimizes this error function.

[0069] This analysis is generally performed using the signal enhancement S_E{i, t{n)) = S(i, t( ))/S(i, t(0)) , where S(i, t(0)) is a pre-contrast acquisition. This transformation eliminates the need for M₀ estimation when comparing modelled signal curves to the ones observed in the DCE-M RI volumes. Combining Equations 7 and 9, the signal enhancement can be expressed explicitly as:

[0070] An overview of this procedure to estimate K_trans and k_ep is shown in FIG. 7A. This procedure 700 consists of finding the set of K_trans and k_ep values which show a better approximation to the DCE-MRI signal enhancement data. These parameters are estimated for each voxel, and this analysis may be performed independently for each voxel. The DCE-MRI data 701 is considered a signal enhancement time intensity curve S. The PK parameters which better describe this are estimated iteratively. At each iteration a pair of candidate parameters K_trans and k_ep are evaluated. This is done by generating 703 model concentration data, for example a model agent concentration curve C_est, from the candidate parameters by using a PK model 706 (such as the Tofts model) and an AIF function 709 (such as the Weinmann model). The model concentration data, for example C_est, is converted 712 to estimated signal enhancement data, for example to an estimated signal enhancement curve S_est, by using Equations 7 and 10, which requires the estimation of T₁₀. The estimated signal enhancement data, for example S_est, is then compared 715 to the actual signal enhancement data 718, for example to the actual signal enhancement curve S, by an error metric (such as the mean square error) to quantify how well the candidate PK parameters describe the real DCE-MRI data. This procedure is repeated iteratively and halts 721 when a minimum fit error is found. [0071 ] FIG. 7B depicts an example of use of the T₁₀ estimation depicted in FIG. 2, and described above, to estimate the pharmacokinetic parameters as depicted in FIG. 7A.

F. Rigid body transformation

[0072] Rigid body transformations were applied in this research as a means to both simulate and correct motion within the variable flip angle MRI sequence. Rigid body motion has the properties of preserving distances between points. In three dimensional volumes, such as the ones analysed here, a rigid body transformation consists of translations (trans_x, trans_y, and trans_z ) and rotations (rot_x, rot_y, and rot_z ) around the volume being transformed center on all the axes (x, y and z).

[0073] Any point p may be transformed according to:

^'i^' P) = I snt

(18)

[0074] The translation segment is expressed by:

trass.

(19)

While the rotation part of this transformation model is given by:

(20)

R_k is the rotational matrix transformation in the k axis. The R_k matrix is shown below, wherein (R_y and R_z are analogous): I 0 0

1L =

0 siairot_<;) eossCroi )

(21 )

G. Groupwise Normalized Mutual Information

[0075] The Groupwise Normalized Mutual Information (GNM I) was used as a gold standard similarity metric for the simultaneous registration of variable flip angle sequences [9]. It extends the Normalized Mutual Information (NM I) [22] for sets of volumes by summing the NMI of each volume I_n to a reference volume /. The reference volume is generated by taking the mean intensity at each voxel of all the volumes being registered.

[0076] Mutual information similarity metrics, in the context of image registration, assess how much information from image I₂ is gained with the knowledge of image I_t [13]. The mutual information between two volumes (I_t and l₂) is measured by [14], [15]:

Where H^) is the marginal entropy of I_t:

(23)

a is the set of intensities found in volume l .

[0077] In the context of two different SPGR volumes with different flip angles or pulse repetition times, the intensities at each voxel, under no motion or noise, would follow the relation shown in Equation 10 and there is a one to one mapping from the intensity on each voxel in image to the same position in in l₂. Thus, the knowledge of l₂ gives complete information of I_t and, consequently, a Mutual Information measure can be used to assess the misalignment between these two images.

[0078] This dependency between the intensities from two images can be evaluated with several other formulations with similar concepts [13]. Mutual Information measures show difficulties at evaluating alignment on images with large non-overlapping regions. In these cases there is a lower number of sample points (only the overlapping ones) which hinders the probability distribution analysis (PI^ and causes higher dependence measures to be found with a greater misalignment and lack of overlap. The Normalized Mutual Information (NM I) decreases this effect [16], and is given by:

Hi A) + H(B)

Hi 4 B)

K '' (24)

Where H(/i, /₂)_> the joint entropy between l and I₂, is given by Equation 25 [14]

( , b)

(25)

[0079] Groupwise Normalized Mutual Information is a similarity metric that extends the Normalized Mutual Information to sets of volumes being registered simultaneously [9]. FIG. 8A illustrates an example of a groupwise registration of a series of variable flip angle image volumes by simultaneously registering the whole set of volumes. It is performed by summing the Normalized Mutual Information between each of the volumes /„ and a reference volume /:

the reference volume is created as a the mean intensity of all the transformed volumes at each voxel.

H. B-Spline free-form deformations

[0080] B-Spline free-from deformations (FFDs) form a non-linear model. Unlike linear transformations which translate, rotate and scale all the image on the same manner, they affect the image locally, managing to achieve more complex deformations.

[0081] The basic concept of free-form deformations is to have a mesh of uniformly distributed control points around the volume [17]. By displacing a control point, a proportional deformation on the same direction is created around that region. Thus, the deformation of any point inside the volume is given by the composition between the deformation on surrounding controlling points.

[0082] B-Splines are one of the classes of functions used to govern this composition of displacements and it may be regarded as an interpolation function. It has a limited support, thus, deformations on each control point only affect a limited region around it and not the whole volume [17]. The greater the order of the B-Spline function the greater is the number of control points around the transformed point that are needed to define its deformation.

[0083] For any point [x, y, z] in the volume, its third order B-Spline transformation T (x. y. z) is defined by a linear combination of the 64 surrounding control points φ [19]:

3 ¾

where i =

v =

y/δ - Yy/S\ and w = ζ/δ - [ζ/δ\. δ is the uniform spacing of the original

(deformation free) control mesh.

[0084] This way, j, k indicate the index of a control point on the mesh and u, v, w shows how positioned is the point between the control points. The third order B-Spline function B is given by:

[0085] FIG. 8B shows an example of a B-Spline transformation, as well as the control point mesh before and after deformation. An original image with a deformation free mesh is shown on the left while a transformed image with the respective control mesh after deformation is shown on the right.

I. Interpolation

[0086] Image volumes are acquired discretely, as a collection of samples from the acquired signal. However, even though the volumes are restricted to uniformly spaced voxels, the image transformation is performed, for example in a continuous manner, shifting intensities found in a voxel to regions in between voxels. Hence, any transformation model requires an interpolation step to transform the discrete data found in volumes into a continuous function [20]. In this work, a tri-linear interpolation method was applied for its speed and simplicity. II. EXAMPLE

A. Original data

[0087] DCE-MRI scans from thirteen patients with advanced colorectal adenocarcinomas were acquired on a 1 .5T GE scanner using a T1 -weighted, gradient-echo, fat-suppressed sequence (LAVA). The actual dynamic sequence is not relevant in this work, but the SPGR variable flip angle acquired to estimate 7^"io were used. These sequences were acquired with T_R = 4.5ms and T_E = 2.2ms, each composed of three or four volumes with different flip angles (3°, 9°, 12° and/or 15°). Since some of the patients had two pairs of sets of acquisitions, before and after chemoradiotherapic treatment, there was a total of twenty sets of variable flip angle sequences.

[0088] These volumes had mean resolution (in millimetres per voxel) of 0.8828 x 0.8828 x 4.65 (between 0.7813 and 0.9375 in the x and y directions and 4 and 5 in the z direction). All images had 512 x 512 dimensions and the volumes between 20 to 28 images.

[0089] The colorectal region is the object of interest of these acquisitions. To reduce computation time, this region was cropped from the original volumes. Each image was reduced to an in-plane 80 x 80 voxels region around the colon and rectum. The upper and lower three slices were also removed from the volume.

B. Simulated data

[0090] To assess the effects of motion on T_w estimation, ground truth volumes with known relaxation time and motion are required. However, from acquired clinical exams this is not possible, as the motion, which is object of this study, will be present, but unknown. Thus, simulated ground truth variable flip angle volumes and Tw and Mo maps were generated.

[0091] The ground truth maps were created based on the original clinical data used: for each variable flip angle sequence, T_w and M₀ were estimated according to the procedure shown in Section l-D above and were considered as the ground truth maps. Creating ground truth data in this way ensures a similar distribution of values as those found in real data. The variable flip angle volumes were generated by applying these maps to Equation 10. For each ground truth generated sequence, the chosen flip angles were the same as the original sequence, thus, some of them were composed of four volumes while others only three.

[0092] 1) Applying rigid motion to simulated data: Ground truth data with randomly generated known rigid motion were applied in this work to compare different registration algorithms (Section IV-B below). The rigid transformations were generated using a random generator with a normal distribution and zero mean for each transform parameter. The rotation parameters (rot_x, rot_y and rot_z), however, have different scales of motion, according to the size of the volume in the corresponding direction. In other words, a translation of 1 mm will cause a mean displacement of the same magnitude, while a rotation of 15° will have a mean displacement proportional to the size of the volumes. For the rotation parameters a volume size dependent scale was applied to normalize the random motion. A normalizing parameter s_k was computed according to each volume size to compensate for this effect. Thus, each transform parameter 9{k) was generated by: (29)

σ is the distribution variance and controls the amount of displacement associated with this transformation.

[0093] 2) Applying non-linear motion to simulated data: Non-linear motion, in form of non-linear B-splines free-form deformations, was applied to the simulated variable flip angle sequences to evaluate the effects of motion on 7~io estimation. Similarly to the generation of rigid motion, these parameters were created by using a normal distribution with zero mean. To emulate the effect of motion in sequentially acquired volumes, we combined random motion at each sequential flip angle cumulatively. The deformations were generated by moving the mesh control points by a random value from a normal distribution with zero mean. To simulate the effect of motion in sequentially acquired volumes, we applied random motion at each sequential flip angle volume incrementally. Hence, for the m-th flip angle, the applied transformation was:

where T_m-_\ is the transform applied to the previous flip angle and T_g is a new randomly generated transform. This causes the effect of an incremental motion happening to the patient as each new volume is acquired.

[0094] The new random transform is then created by moving the original control mesh point k: III. METHODS

A. Groupwise rigid registration

[0095] A groupwise rigid registration framework was applied to evaluate different metrics for correcting motion within the SPGR volume sequence. The registration framework consists of the optimizer, similarity metric and transform. The transform defines which types of motions are allowed to be applied to the volumes while searching for the correct alignment. The similarity metric quantifies how well a given transformation corrects the motion between these volumes. Finally, the optimizer iteratively adjusts the transform parameters aiming to maximize (some minimizing) the similarity measure.

[0096] As a side note, some metrics present decreasing results when there is a greater similarity (mean intensity difference between two images), while others, such as GNMI, have a high degree of similarity reflected with higher output values. Thus, depending on the metric, the optimizer will act by maximizing or minimizing this function.

[0097] One of the N flip angle volumes was considered static, while the other N - 1 were subjected to rigid body transformations (see Section l-F above). Thus, the set of transform parameters Θ being optimized in the search for the registration between these volumes is given as a composition of all these N - 1 rigid body transformations, totaling 6(N - 1 ) parameters:

θ !

trans, i-n [0098] In this work, the optimizer utilized was the Levenberg-Marquandt, which is a damped version of the Gauss-Newton Method, from the alglib library [21 ], through other optimizers may be used, such as e.g. Gradient Descent, Hill Climbing, quasi-Newton, Conjugate Gradients, and Simulated Annealing.

[0099] Two similarity metrics were used and compared: the GNM I ll-G and a new metric based on the T_w regression error, which is further explained in Section lll-B below.

[0100] Moreover, a pre-processing filtering stage was performed before registration to reduce the effect of noise over the image. This step consisted of applying a Gaussian blurring filter with a standard deviation of 1 mm.

B. 7~io regression error similarity metric

[0101] Estimating T_w consists of fitting a line over two or more points from the SPGR volume sequence at each voxel. In the presence of more than two points, this regression will yield a residual error which may be caused by motion between the volumes, noise and/or (to a much lesser extent) changes in the magnetic properties of the tissue. Hence, the lack of fit found in this regression may, among other causes, be a consequence of misalignment between these volumes. Conversely, the misfit found could be used to assess the amount of patient motion between the MRI acquisitions.

[0102] Taking advantage of a modelled relationship between the intensities over a series of images has previously been proposed as a similarity metric to determine results in image registration. In [6], the error between the data and the PK model was used as a non-linear motion correction similarity metric on DCE-MRI volumes. Moreover, the Tofts Model parameter was fitted simultaneously to the motion correction by assuming that better alignment between the images would lead to a better model fit, and vice versa.

[0103] In this work, a similar approach was taken to develop a similarity metric for the variable flip angle volume sequence. The T_w regression should have a better fit when there is a better alignment. A non-limiting example of a suitable similarity metric is the mean square error between the observed values Y {n) (Equation 12) and the fitted values (Y(n) = AX(n) + B): sin½_¾ = -^- I — F{n) ^'

(33)

IV. EXPERIMENTS

A. Evaluating the effects of motion on the estimation of pharmacokinetic model parameters

[0104] Two experiments were performed to assess the effects of motion within the SPGR volume sequence of a DCE-MRI sequence used to estimate T_w and subsequent pharmacokinetic (PK) model parameter estimation. First, the error in T_w caused by motion within these volumes was estimated by applying varying levels of random levels of B-splines free-form deformations to simulated MRI variable flip angle sequences with known T_w values (see Section l l-B). After applying the deformations, 7~i₀ was estimated and compared to the ground truth, providing a measure of the estimation error as a function of disalignment.

[0105] The second experiment consisted of evaluating how deviated 7~io values affected PK model parameter estimation. Simulated DCE-MRI signals were generated with known K_trans and k_ep values and T_w = 1 .0s. Erroneous T_w values were used to estimate these PK parameters and the obtained results were compared to the original ones, providing a measure of how 7io error propagates to PK analysis. In this experiment K_trans ranged from 0.01 to 2.00 and k_ep from 0.01 to 5, providing an average estimation error over different curves. Other parameters of this experiment were: r_l=4.5(mM.s)^'1 , r₂=5.5(mM.s)^'1 , T_R=4.5ms, T_E=0.22ms and a 15° flip angle. The synthetic DCE- MRI curves were generated with samples at every 12 seconds for 5.2 minutes.

[0106] 1) Evaluation of motion on the estimated T_w: We assessed the effect of motion on T_w estimation by applying non-linear transformations to uncorrupted SPGR volume sequences and measuring this error for different degrees of displacement. The procedure to perform this test is presented below, with an overview of the procedure 900 presented in FIG. 9.

1 ) A volume 903 of ground truth parameter maps ( ₁₀ and M₀) was created, as described in Section ll-B.

2) Non-linear B-Splines were applied 906 sequentially to these maps to simulate patient motion during the scans (as described in Section II-B2).

3) Each of the motion corrupted maps were then used to create 909, or generate, MRI volumes, in particular SPGR volume sequences, with variable flip angles 912.

4) The relaxation time, f_w, was then estimated 915 from the motion corrupted MRI volumes using the procedure explained in Section ll-D.

5) The error between f₁₀ and r₁₀ was then computed 918 as a function of the displacement caused by the non-linear motion applied. [0107] The B-Splines control point mesh grid used was 5 x 5 x 5, with increasing standard deviations for the random motion generator, yielding mean incremental displacements from 0mm to 2mm. Ground truth T-|₀ and M-|₀ maps are used to create variable flip angle volume sequences with non-linear motion and these sequences are then used to estimate Ti₀ and Mi₀ maps, which are then compared by the mean absolute error with the ground truth maps.

[0108] 2) Evaluation of pharmacokinetic model parameters errors due to inaccurate 7~io: The following experiment was created to quantify the effect of inaccurate T_w in the estimation of the Tofts Model parameters K_trans and k_ep.

1 ) A set of ground truth parameters were chosen: K_trans, K_ep and f₁₀.

2) The two model parameters were used in conjunction with the Tofts model and the Arterial Input Function to create the true concentration (c_t) curve using Equation 7.

3) The accurate relaxation time ( ₁₀) was used to create, from the concentration curve, a DCE-MRI signal curve S (Equation 8).

4) An inaccurate relaxation time was created: T₁₀ = f₁₀ + dev.

5) The model parameters were estimated (K_trans and k_ep) as shown in Section l-E above.

6) The estimated parameters were compared to the ground truth ones.

[0109] The ground truth pharmacokinetic parameters varied from 0.01 to

2.00 (K_trans) and 0.01 and 5 (k_ep), giving an average estimation error over different curves. Table I shows the other parameters used. TABLE f

PA RA ETERS USED O THE EVALUATION OF THE EFFECTS OP JHACCTI A ES »Ϊ<^" ¾ O TO ARAM ACOKI ETFC MODEL P AM TER

ESTIMATION,

[0110] The experiment was performed with three values of T_w: 0.5s, 1 .0s and 1.5s, which are within the most common range of relaxation times found on colorectal tumour tissue [8]. The curves were generated with samples at every 12 seconds for 5.2 minutes, which is similar to DCE-MRI sequences acquisition rates.

B. Evaluation of similarity metrics for rigid body motion correction

[011 1] The 7^"io regression and GNMI similarity metrics were tested using the rigid registration framework (shown in Section ll l-A above) with ground truth sequences of variable flip angle MRI volumes. The procedure to perform this test is presented below, with an overview of the procedure 1000 presented in FIG. 10.:

1 ) Volumes 1003 of ground truth parameter maps ( ₁₀ and M₀) were created from the original clinical data, as explained in Section ll-B.

2) For each flip angle, a random rigid transform was created 1006 and applied to the T_w and ₀ volumes as described in Section I I-B1 . 3) The transformed T_w and M₀ volumes were used to create 1009, or generate, motion corrupted MRI volumes, in particular SPGR volume sequences, with variable flip angle 1012 using Equation 10.

4) The volumes were registered 1015 (Section lll-A) using both similarity metrics (Sections lll-B and l-G).

5) The mean displacement was computed 1018 before and after registration as a metric to compare the similarity measures.

6) The 7~io estimation error was also computed before and after registration for both the similarity measure results.

[0112] This procedure was performed over twenty simulated ground truth parameter maps with rigid simulated mean displacement varying from 0.42mm to 2.1 1 mm. This procedure provided a method for evaluating how well the groupwise rigid body motion registration correction procedure performed.

C. Estimation of the rigid motion found on SPGR variable flip angle volume sequences

[0113] A third experiment was performed to estimate the amount of rigid motion found within SPGR variable flip angle sequences (Sections IV-A1 , IV-A2 and IV-B), assuming there is some motion within the SPGR variable flip angle sequence of MRI volumes and then evaluated the effects and methods to correct it.

[0114] Rigid motion was quantified by applying the groupwise rigid registration method (Section lll-A) to several SPGR sequences. This was performed using both similarity metrics (T_w regression error and Groupwise

Normalized Mutual Information); the original motion was considered as the displacement caused by the obtained registration transforms. [0115] Since this procedure assumes that the registration framework is indeed providing a better alignment to the volumes (which may not always be true), a third measure of this motion was performed by considering, for each sequence, the minimum displacement result between transforms obtained with the different similarity measures. This avoids results with larger motions, which may be found when the registration process fails, thus providing a lower bound to the natural misalignment of these sequences.

V. RESULTS AND DISCUSSION

A. Evaluating the effects of motion on the estimation of pharmacokinetic model parameters

[0116] 1) Evaluation of motion on the estimated T_w: The T_w estimation error was determined as a function of motion. FIG. 1 1 presents this experiment results for each tested motion. FIG. 1 1 shows an estimated T-|₀ error due to motion between the SPGR sequence volumes for each tested motion.

[0117] The motion simulated in this experiment was applied incrementally to each volume of the flip angle sequence. Thus the last angle volume is created by adding the deformation from the previous volume to a new randomly generated incremental deformation. The estimated T_w error presented in FIG. 11 shows that there is a linear increase of error for an increasing amount of motion applied increase of error for an increasing amount of motion applied to the images. For the highest motion value tested, 2.02mm, the mean absolute 7~io error surpassed 30%. For a mean displacement of 0.8mm, close to the in- plane voxel length, the error found was around 17% and, within half standard deviation, reaching 25%. Comparing the overestimated with the underestimated 7^~io values, both had similar mean error and standard deviation, yet, in 59% of the time the value was overestimated.

[0118] 2) Evaluation of pharmacokinetic model parameters errors due to inaccurate T_w: FIG. 12 exposes how inaccuracies in T_w affect the estimation of ^Ktrans and /c_ep. FIG.12 depicts K_trans and k_ep error due to T-io inaccuracies. To develop the data in FIG. 12 the Orton AIF was used instead of the Weinmann AIF. The most notable result is that under estimations of T_w have a significantly greater effect on the pharmacokinetic (PK) model parameter estimation. This discrepancy is so high that may even motivate the development and application of biased T_w estimation techniques. Moreover, there is an almost linear dependency between both errors.

[0119] For K_trans, under estimated values of 7^~ ₀ show a higher propagation of error when estimating the model parameter, a behaviour which is not repeated on k_ep estimation.

[0120] An important aspect of this experiment is that different sets of conditions for pharmacokinetic analysis (such as other relaxivity constants and dosages) result in very different model parameters error curves. If the dosage changed to 1 mM/kg, K_trans would go to over 100% with a -30% 7^~io inaccuracy. This indicates that different exam conditions may propagate errors at different scales and that there may be better settings which could make pharmacokinetic analysis more robust.

[0121] An example of incorrect estimation of K_trans and k_ep due to errors in 7^~io is shown in FIG. 13. A - 20% inaccuracy is shown in T-|₀. As seen, the computed signal enhancement (SE) curve (black) does not fit the observed signal enhancement curve (lower curve) and an overestimated K_trans parameter generates a matching signal enhancement curve. The signal enhancement curve generated with the correct pharmacokinetic parameters fails to fit the observed ground truth curve, while incorrect parameters show a good match.

[0122] These results can be further combined with the ones regarding motion and its consequences on T_w (Section V-A1 ) to give an idea of the errors seen on the final clinical parameter estimation due to motion. Considering the estimated relaxations that show a negative bias (which happens 41 % of the time), with a mean displacement of 0.8mm), a K_trans error of over 20% is caused by the 16% inaccuracy on Tw estimation. Doubling the motion to 1 .6mm, the K_trans ^error exceeds 40%. As for the estimated Tw .with values higher than the real ones, the error does not increase when estimating the model parameters, but motion close to voxel size will still cause an average of 15% error in K_trans.

[0123] Considering that in 59% of the cases there is an overestimation of Tw (as shown in Section V-A1 ), an average K_trans error was estimated as a function of the absolute T_w error. This is illustrated in FIG. 14 and shows an overall trend of the error increasing slightly from Tw to K_trans estimation.

B. Evaluation of similarity metrics for rigid body motion correction

[0124] The similarity metrics were evaluated by the mean displacement found between the original (pre-motion) and the registered volumes. The amount of motion that was simulated, the target registration and Tw estimation errors for each similarity metric before and after registration are shown in Table II. TAfii.fi I

fti.i>iwi.s wrm SWUIATEB tama m a».

7.4¼ ^" Γ ^'<;.«*^■ SS..¾ ± ftiSsrssi «. ¾

S5.SS ■2.S« ¾S2 :·: i- !&tt S2.8¾

2&J% ^■i il i

2i,;¾ it .14^■·^■ a i ij.7¾ ^■,07 I" S8.!¾

2.i¾ 2«.7¾ i« 4'S i 72 :+ ¾S¾sn _iiS«

2.42 i!.53siis 2S.2S ■7.2S 1.5:: 25 ¾

[0125] In all cases, the present Tw regression metric showed a lower registration error than the GNMI metric, resulting on average 18% better alignment. The difference between the registration outcomes was almost negligible for low levels of motion, but became very significant as the motion increased. A displacement of the length of an in-plane pixel was almost reached on the most extreme case tested.

[0126] The new metric presented herein also showed a lower standard deviation, indicating greater consistency. The standard deviation of the mean displacement error from the registrations using T_w regression metric was, on average, lower than when the GNM I was utilized.

[0127] The T_w error was also evaluated for different levels of motion and registration (Table II). Its estimation almost always improved as the registration error decreased. It almost always showed decreasing values as the motion was corrected. The exception to this trend was found for 0.40mm of simulated motion, where the registration did correct for motion but did not improve T_w estimation. This was due to the fact that the interpolation step of the transformation does introduce some errors to the transformed volume. Thus, even if the transformation found does actually correct for motion, the corrected volume intensities may not be accurate. At the same time, this means that a better interpolation method could decrease even more the T_w error after motion correction.

[0128] In terms of pharmacokinetic (PK) model parameter estimation, when subjected to 1 .29mm mean simulated motion, after performing rigid registration with the proposed similarity metric, the T_w estimation error dropped by 6.4%, falling down to 14.1 %. This would result in even greater decrease on

K_trans 6ΓΓ0Γ.

C. Assessing the rigid motion found on SPGR variable flip angle volume sequences.

[0129] Table III presents the estimated motion found by applying rigid registration to the original MRI volumes. A significant discrepancy was found between the motion estimated using the T_w regression and the GNM I similarity metrics. The former found a mean displacement almost twice as large as the latter. This result not only shows that there may be considerable motion found within these volumes, but also that the choice of the metric when correcting for motion may have a relevant effect on the final registration. The results in Section V-B show a lower registration error when the T_w regression error approach was applied.

TABLE El

SlA! i S3^;L .^':i:i E :T H ^'i.S K A APO¾ «KOUi' ¾E KJittSD

OF SEAL DATA,

[0130] Another approach taken to estimate the motion found in these volumes was to evaluate the lowest of the registration displacement found using each metric for each volume. This lower bound for the existing patient motion was of 0.38mm, which is close to the motion using the GNM I.

[0131] Considering the lowest estimated motion using the GNMI metric, an average error of around 10% is expected for T_w leading to about 15% K_trans. error. Alternatively, using the motion estimated by the newly developed 7~io regression metric, 0.78mm of displacement is expected to cause an estimation error of more than 20% on the pharmacokinetic model parameter

Ktrans■

[0132] This analysis only considered rigid motion, thus, a higher average displacement could be present in colorectal MRI due to non-rigid movement.

VI. CONCLUSION

[0133] This work shows how patient motion within a flip angle sequence involving the acquisition of three or more images (or volumes) with different flip angles or pulse repetition times, for example a variable flip angle or pulse repetition time sequence can affect DCE-MRI analysis of tumours (particularly colorectal). A clear relation between motion and K_trans error was shown, quantifying how much misalignment may be present during relaxation time estimation without compromising the pharmacokinetic analysis. In contrast, k_ep estimation was very robust to T_w error. Conversely, it also indicates how much these volumes must be corrected for motion.

[0134] Our experiments also showed that there exists a relevant amount of motion within variable flip angle MRI sequences, which would cause, according to our estimations, an error between 15% and 25% on K_trans , motivating the application of the present motion correction techniques during PK analysis of DCE-MRI exams for example our new similarity metric for motion correction in SPGR variable flip angle sequences based on the Tw regression error. Our novel motion correction technique consistently showed lower registration error on sets of simulated volumes with simulated motion when compared to the GN Ml metric.

[0135] We now describe an embodiment of a system providing for motion correction to improve T₁₀ estimation, such as those described above. In an embodiment our present method and system may be carried out by programming logic executed in a computing environment. With reference to FIG. 15, shown is a schematic block diagram of a system including a computing environment 101 according to an embodiment of the present disclosure. The computing environment 101 may include one or more computing devices 401 . Each computing device 401 may include at least one processor circuit, for example, having a processor 402, a memory 404, a user interface 406, and an input/output (I/O) device 408, which may be coupled to a local interface 407. To this end, each computing device 401 may comprise, for example, at least one server computer or like device. The local interface 407 may comprise, for example, a data bus with an accompanying address/control bus or other bus structure as can be appreciated.

[0136] The computing environment 101 may comprise, for example, a handheld device, a portable device, a computer, server, dedicated processing system, or other system, as can be appreciated. The hand-held device can be, for example, a smart mobile phone, a tablet, or other mobile smart device. The computing environment 101 of such device may include various input devices such as a keyboard, microphone, mouse, touch screen, or other device, as can be appreciated. By way of example, the system can comprise a stand-alone device or part of a network, such as a local area network (LAN), GPRS cellular network or wide area network (WAN).

[0137] The processor 402 may include a central processing unit (CPU) or a semiconductor-based microprocessor in the form of a microchip. In addition, the processor 402 may represent multiple processors and the memory 404 may represent multiple memories that operate in parallel. In such a case, the local interface 407 may be an appropriate network that facilitates communication between any two of the multiple processors, between any processor and any one of the memories, or between any two of the memories, etc. The processor 402 may be of electrical or optical construction, or of some other construction as can be appreciated by those with ordinary skill in the art.

[0138] Stored in the memory 404 are both data and several components that are executable by the processor 402. In particular, stored in the memory 404 and executable by the processor 402 may be programming logic 1 10, and potentially other applications. Also stored in the memory 404 may be a data store 1 1 1 storing for example acquired MR image data, modeled data, and other data. In addition, an operating system may be stored in the memory 404 and executable by the processor 402.

[0139] It is understood that there may be other applications that are stored in the memory 404 and are executable by the processor 402 as can be appreciated such as an operating system. Where any component discussed herein is implemented in the form of software, any one or more of a number of programming languages may be employed such as, for example, C, C++, C#, Objective C, Java^®, JavaScript^®, Perl, PHP, Visual Basic^®, Python^®, Ruby, Flash^®, or other programming languages.

[0140] A number of software components may be included in the programming logic 1 10 and may be stored in the memory 404 and executable by the processor 402. In this respect, the term "executable" means a program file that is in a form that can ultimately be run by the processor 402. Examples of executable programs may be, for example, a compiled program that can be translated into machine code in a format that can be loaded into a random access portion of the memory 404 and run by the processor 402, source code that may be expressed in proper format such as object code that is capable of being loaded into a random access portion of the memory 404 and executed by the processor 402, or source code that may be interpreted by another executable program to generate instructions in a random access portion of the memory 404 to be executed by the processor 402, etc. An executable program may be stored in any portion or component of the memory 404 including, for example, random access memory (RAM), read-only memory (ROM), hard drive, solid-state drive, USB flash drive, memory card, optical disc such as compact disc (CD) or digital versatile disc (DVD), floppy disk, magnetic tape, or other memory components.

[0141] The operating system may control the execution of these programs as well as other programs and provides scheduling, input-output control, file and data management, memory management, and communication control and related services. The operating system may be executed to control the allocation and usage of hardware resources such as the memory, processing time and peripheral devices in the computing environment 101 . In this manner, the operating system may serve as the foundation on which applications depend as is generally known by those with ordinary skill in the art.

[0142] The memory 404 is defined herein as including both volatile and nonvolatile memory and data storage components. Volatile components are those that do not retain data values upon loss of power. Nonvolatile components are those that retain data upon a loss of power. Thus, the memory 404 may comprise, for example, random access memory (RAM), read-only memory (ROM), hard disk drives, solid-state drives, USB flash drives, memory cards accessed via a memory card reader, floppy disks accessed via an associated floppy disk drive, optical discs accessed via an optical disc drive, magnetic tapes accessed via an appropriate tape drive, and/or other memory components, or a combination of any two or more of these memory components. In addition, the RAM may comprise, for example, static random access memory (SRAM), dynamic random access memory (DRAM), or magnetic random access memory (MRAM) and other such devices. The ROM may comprise, for example, a programmable read-only memory (PROM), an erasable programmable read-only memory (EPROM), an electrically erasable programmable read-only memory (EEPROM), or other like memory device.

[0143] The memory 404 typically comprises a native operating system, one or more native applications, emulation systems, or emulated applications for any of a variety of operating systems and/or emulated hardware platforms, emulated operating systems, etc. For example, the applications may include application specific software or code which may be configured to perform some or all of the systems and methods for providing motion sensitized and motion suppressed quantitative imaging for the medical and neuroscience applications described herein. In accordance with such embodiments, the application specific software is stored in memory 404 and executed by the processor 402. One of ordinary skill in the art will appreciate that the memory 404 can, and typically will, comprise other components which have been omitted for purposes of brevity.

[0144] Also, the processor 402 may represent multiple processors 402 and/or multiple processor cores and the memory 404 may represent multiple memories 404 that operate in parallel processing circuits, respectively. In such a case, the local interface 407 may be an appropriate network that facilitates communication between any two of the multiple processors 402, between any processor 402 and any of the memories 404, or between any two of the memories 404, etc. The local interface 407 may comprise additional systems designed to coordinate this communication, including, for example, performing load balancing. The processor 402 may be of electrical or of some other available construction.

[0145] The user interface 406 may comprise the components with which a user interacts with the processing circuit and therefore may comprise, for example, a keyboard, mouse, and a display, such as a liquid crystal display

(LCD) monitor. The user interface 406 can also comprise, for example, a touch screen that serves both input and output functions. The one or more I/O devices

408 are adapted to facilitate communications with other devices or systems and may include one or more communication components such as a modulator/demodulator (e.g., modem), wireless (e.g., radio frequency (RF)) transceiver, network card, etc. For example, where the system comprises a personal computer, these components may interface with one or more user input devices. A display may be provided that may comprise a computer monitor, a plasma screen for a PC, a liquid crystal display (LCD) on a hand held device, or other display device.

[0146] The system can be designed to communicate with a computing environment 101 or processor circuit, providing the measured data to the computing environment 101 . The computing environment can include programming logic 1 10 for receiving or acquiring the measured data. The system can be included in a mobile smart device such as a smart phone or a tablet. A smart device, such as a smart phone or tablet, may include a processor circuit on which the programming logic 1 10 can be executed. For example, the programming logic may take the form of an app that can be downloaded onto the smart device and executed by the processing system.

[0147] As described herein, the processor circuit, in particular the software provided on the processor circuit, may be configured to receive acquired data.

[0148] With further reference to FIG. 15, the system may further include a network interface device that may comprise various components that can be used to transmit and/or receive data over a network environment to other components in an MRI system. For example, the network interface may include a device that can communicate with both inputs and outputs, for instance, a modulator/demodulator (e.g., a modem), wireless (e.g., radio frequency (RF)) transceiver, a telephonic interface, a bridge, a router, network card, etc.). The system may communicate with one or more computing devices (not shown) via the network interface over a network. The system may further comprise a peripheral interface that may support various interfaces including, but not limited to IEEE-1394 High Performance Serial Bus (Firewire), USB, a serial connection, and a parallel connection.

[0149] The system shown in FIG. 15 may be embodied, for example, as a magnetic resonance apparatus that includes an image acquisition system, and which includes a processing module or logic for performing conditional data processing, and may be implemented either off-line or directly in a magnetic resonance (MR) apparatus. For such embodiments, the system may be implemented as a multi-channel, multi-coil system with advanced parallel image processing capabilities, and direct implementation makes it possible to generate immediate images available for viewing immediately after image acquisition, thereby allowing re-acquisition on-the-spot if necessary. Examples of apparatus in which the present systems and methods may be implemented are described in U.S. Patent No. 5,993,398 and No. 6,245,027 and U.S. Publication No. 201 1/0181285, which are incorporated by reference as if fully set forth herein expressly in their entireties.

[0150] Although the programming logic 1 10, and other various systems described herein may be embodied in software or code executed by general purpose hardware as discussed above, as an alternative the same may also be embodied in dedicated hardware or a combination of software/general purpose hardware and dedicated hardware. If embodied in dedicated hardware, each can be implemented as a circuit or state machine that employs any one of or a combination of a number of technologies. These technologies may include, but are not limited to, discrete logic circuits having logic gates for implementing various logic functions upon an application of one or more data signals, application specific integrated circuits (ASICs) having appropriate logic gates, field-programmable gate arrays (FPGAs), or other components, etc. Such technologies are generally well known by those skilled in the art and, consequently, are not described in detail herein.

[0151] The flowcharts herein, for example, FIGs. 2 and 7A-B, show the functionality and operation of an implementation of portions of the programming logic. If embodied in software, each block may represent a module, segment, or portion of code that comprises program instructions to implement the specified logical function(s). The program instructions may be embodied in the form of source code that comprises human-readable statements written in a programming language or machine code that comprises numerical instructions recognizable by a suitable execution system such as a processor 402 in a computer system or other system. The machine code may be converted from the source code, etc. If embodied in hardware, each block may represent a circuit or a number of interconnected circuits to implement the specified logical function(s).

[0152] Although the flowcharts of FIGs. 2 and 7A-B show a specific order of execution, it is understood that the order of execution may differ from that which is depicted. For example, the order of execution of two or more blocks may be scrambled relative to the order shown. Also, two or more blocks shown in succession in FIGs. 2 and 7A-B may be executed concurrently or with partial concurrence. Further, in some embodiments, one or more of the blocks shown in FIGs. 2 and 7A-B may be skipped or omitted. In addition, any number of counters, state variables, warning semaphores, or messages might be added to the logical flow described herein, for purposes of enhanced utility, accounting, performance measurement, or providing troubleshooting aids, etc. It is understood that all such variations are within the scope of the present disclosure.

[0153] Also, any logic or application described herein, including the programming logic 1 10, that comprises software or code can be embodied in any non-transitory computer-readable medium for use by or in connection with an instruction execution system such as, for example, a processor 402 in a computer system or other system. In this sense, the logic may comprise, for example, statements including instructions and declarations that can be fetched from the computer-readable medium and executed by the instruction execution system. In the context of the present disclosure, a "computer-readable medium" can be any medium that can contain, store, or maintain the logic or application described herein for use by or in connection with the instruction execution system.

[0154] The computer-readable medium can comprise any one of many physical media such as, for example, magnetic, optical, or semiconductor media. More specific examples of a suitable computer-readable medium would include, but are not limited to, magnetic tapes, magnetic floppy diskettes, magnetic hard drives, memory cards, solid-state drives, USB flash drives, or optical discs. Also, the computer-readable medium may be a random access memory (RAM) including, for example, static random access memory (SRAM) and dynamic random access memory (DRAM), or magnetic random access memory (MRAM). In addition, the computer-readable medium may be a read-only memory (ROM), a programmable read-only memory (PROM), an erasable programmable readonly memory (EPROM), an electrically erasable programmable read-only memory (EEPROM), or other type of memory device. [0155] It should be emphasized that the above-described embodiments of the present disclosure are merely possible examples of implementations set forth for a clear understanding of the principles of the disclosure. Many variations and modifications may be made to the above-described embodiment(s) without departing substantially from the spirit and principles of the disclosure. All such modifications and variations are intended to be included herein within the scope of this disclosure and protected by the following claims.

[0156] References

[1 ] T. E. Yankeelov and J. C. Gore, "Dynamic contrast enhanced magnetic resonance imaging in oncology: theory, data acquisition, analysis, and examples," Current medical imaging reviews, vol. 3, no. 2, p. 91 , 2009.

[2] M. a. Zahra, K. G. Hollingsworth, E. Sala, D. J. Lomas, and L. T. Tan, "Dynamic contrast-enhanced MRI as a predictor of tumour response to radiotherapy." The lancet oncology, vol. 8, no. 1 , pp. 63-74, Jan. 2007.

[3] A. R. Padhani, "Dynamic contrast-enhanced MRI in clinical oncology: Current status and future directions," Journal of Magnetic Resonance Imaging, vol. 16, no. 4, pp. 407-422, 2002.

[4] P. S. Tofts, "Modeling Tracer Kinetics in Dynamic GD-DTPA MR Imaging," Journal of magnetic resonance imaging, vol. 7, no. 1 , pp. 91-101 , 1997.

[5] S. C. L. Deoni, B. K. Rutt, and T. M. Peters, "Rapid combined T1 and T2 mapping using gradient recalled acquisition in the steady state." Magnetic resonance in medicine, vol. 49, no. 3, pp. 515-26, Mar. 2003.

[6] M. Bhushan, J. Schnabel, L. Risser, M. Heinrich, J. Brady, and M. Jenkinson, "Motion correction and parameter estimation in dceMRI sequences: application to colorectal cancer," Medical Image Computing and Computer-Assisted Intervention, MICCAI 2011, pp. 476-483, 201 1 .

[7] G. a. Buonaccorsi, C. Roberts, S. Cheung, Y. Watson, J. P. B. O'Connor, K. Davies, A. Jackson, G. C. Jayson, and G. J. M. Parker, "Comparison of the performance of tracer kinetic model-driven registration for dynamic contrast enhanced MRI using different models of contrast enhancement." Academic radiology, vol. 13, no. 9, pp. 1 1 12-23, Sept. 2006.

[8] L. N. Tanner, "Functional Imaging Markers for Tumour Characterisation," Ph.D. dissertation, University of Oxford, 2010.

[9] K. K. Bhatia, J. Hajnal, A. Hammers, and D. Rueckert, "Similarity metrics for groupwise non-rigid registration." Medical Image Computing and Computer-Assisted Intervention, MICCAI 2007, vol. 10, no. Pt 2, pp. 544-52, Jan. 2007.

[10] M. Ingrisch and S. Sourbron, "Tracer-kinetic modeling of dynamic contrast-enhanced MRI and CT: a primer." Journal of pharmacokinetics and pharmacodynamics, vol. 40, no. 3, pp. 281-300, June 2013.

[1 1 ] "On the scope and interpretation of the Tofts models for DCE-MRI." Magnetic resonance in medicine, no. 3, pp. 735-45, Sept. 201 1 .

[12] A. Garpebring, P. Brynolfsson, J. Yu, R. Wirestam, A. Johansson, T. Asklund, and M. Karlsson, "Uncertainty estimation in dynamic contrast-enhanced MRI." Magnetic resonance in medicine, vol. 69, no. 4, pp. 992-1002, Apr. 2013.

[13] J. P. W. Pluim, J. B. A. Maintz, and M. a. Viergever, "Mutual- information-based registration of medical images: a survey." IEEE transactions on medical imaging, vol. 22, no. 8, pp. 986-1004, Aug. 2003.

[14] F. Maes, a. Collignon, D. Vandermeulen, G. Marchal, and P. Suetens, "Multimodality image registration by maximization of mutual information." IEEE transactions on medical imaging, vol. 16, no. 2, pp. 187-98, Apr. 1997. [15] P. Viola and W. M. Wells, III, "Alignment by maximization of mutual information," Int. J. Comput. Vision, vol. 24, no. 2, pp. 137-154, Sept. 1997.

[16] C. Studholme, D. Hill, and D. Hawkes, "An overlap invariant entropy measure of 3D medical image alignment," Pattern Recognition, vol. 32, no. 1 , pp. 71-86, Jan. 1999.

[17] D. Rueckert, L. I. Sonoda, C. Hayes, D. L. Hill, M. O. Leach, and D. J. Hawkes, "Nonrigid registration using free-form deformations: application to breast MR images." IEEE transactions on medical imaging, vol. 18, no. 8, pp. 712-21 , Aug. 1999.

[18] Z. Xie and G. E. Farin, "Image registration using hierarchical B- splines." IEEE transactions on visualization and computer graphics, vol. 10, no. 1 , pp. 85-94, 2004.

[19] F. Ino, K. Ooyama, and K. Hagihara, "A data distributed parallel algorithm for nonrigid image registration," Parallel Computing, vol. 31 , no. 1 , pp. 19-43, 2005.

[20] A. Amanatiadis and I. Andreadis, "A survey on evaluation methods for image interpolation," Measurement Science and Technology, vol. 20, no. 10, p. 104015, Oct. 2009.

[21 ] S. Bochkanov, "Alglib," 2010, http ://^'m loss . orq/softwa re/^'vi ew/231 /.

[22] C. Studholme, D.L.G. Hill, and D.J. Hawkes. "An overlap invariant entropy measure of 3D medical image alignment. Pattern Recognition, 32(1 ):71 - 86, January 1999.

Claims

What is claimed:

A method for motion correction for T-|₀ estimation in a sequence of magnetic resonance (MR) imaging comprising the steps of: a) Acquiring a sequence of MR images with varying acquisition

parameters;

b) Applying a transformation model to at least one of the acquired images to produce a transformed image;

c) Using the transformed image to estimate T₁₀;

d) Generating a model of the transformed image using the estimated

e) Applying an error metric to determine the error between at least one acquired image and the modeled image;

f) Updating the transformation model based on the error determined; g) Applying the updated transformation model to produce a new

transformation of at least one acquired image; and

h) Using the new transformation of the at least one acquired image to provide a new Tw estimation.

The method of claim 1 , wherein steps (d) - (h) are repeated until a desired minimized error has been determined.

The method of claims 1 or 2, wherein step (a) comprises acquiring a variable flip angle or a pulse repetition time sequence of MR images. The method of any of claims 1 -3, wherein the transformation model is selected from the group consisting of linear (rigid-body or affine) or non-linear (parametric or non-parametric) transformations, or combinations thereof.

5. The method of any of claims 1 - 4, wherein the error metric is selected from the group consisting of sum of squared residuals, sum of mean residuals, maximum residual, mean square or any other suitable error metric.

6. The method of any of claims 1 - 5, wherein the transformation model includes a set of parameters and the transformation model is updated by adjusting one or more of the parameters using an optimizer.

7. The method of claim 6, wherein the optimizer is selected from the group consisting of a Levenberg-Marquardt, Gradient Descent, Hill Climbing, quasi-Newton, Conjugate Gradients, and Simulated

Annealing optimizers.

8. The method of any of the foregoing claims, wherein the acquisition parameter is varied, the acquisition parameter being either the flip angle or the pulse repetition time, and the sequence is acquired using a spoiled gradient recall (SPGR) echo pulse sequence, a gradient echo pulse sequence or a spin echo pulse sequence.

9. The method of any of the foregoing claims, wherein the acquired MR images are 2D or 3D images and the new T_w estimation is used to provide a motion corrected T₁₀ map.

10. The method of any of the foregoing claims, wherein the acquired

images are of tissue of a subject and the final T₁₀ estimate is used for pharmacokinetic analysis, tissue characterization and/or image segmentation.

1 . The method of any of the foregoing claims, wherein the image sequence includes three or more MR images and the transformation model is applied groupwise to the images.

12. A system, comprising:

at least one computing device; and

programming logic executable in the at least one computing device, the programming logic comprising logic that:

a) Acquires a sequence of MR images with varying acquisition parameters;

b) Applies a transformation model to at least one of the acquired images to produce a transformed image;

c) Uses the transformed image to estimate T_w ;

d) Generates a model of the transformed image using the estimated 7^~ ₀ ;

e) Applies an error metric to determine the error between at least one acquired image and the modeled image;

f) Updates the transformation model based on the error determined;

g) Applies the updated transformation model to produce a new transformation of at least one acquired image; and

h) Uses the new transformation of the at least one acquired image to provide a new T₁₀ estimation.

13. The system of claim 12, wherein (d) - (h) are repeated until a desired minimized error has been determined.

14. The system of claims 12 or 13, wherein the logic in (a) acquires a variable flip angle or a pulse repetition time sequence of MR images.

15. The system of any of claims 12 - 14, wherein the transformation

model is selected from the group consisting of linear (rigid-body or affine) or non- linear (parametric or non-parametric) transformations, or combinations thereof.

16. The system of any of claims 12 - 15, wherein the error metric is

selected from the group consisting of sum of squared residuals, sum of mean residuals, maximum residual mean square or any other suitable error metric.

17. The system of any of claims 12 - 16, wherein the transformation model includes a set of parameters and the transformation model is updated by adjusting one or more of the parameters using an optimizer.

18. The system of claim 17, wherein the optimizer is selected from the group Consisting of a Levenberg-Marquardt, Gradient Descent, Hill Climbing, quasi-Newton, Conjugate Gradients, and Simulated Annealing optimizers.

19. The system of any of the foregoing claims, wherein the acquisition parameter is varied, the acquisition parameter being either the flip angle or the pulse repetition time, and the sequence is acquired using a spoiled gradient recall (SPGR) echo pulse sequence, a gradient echo pulse sequence or a spin echo pulse sequence.

20. The system of any of the foregoing claims, wherein the acquired MR images are 2D or 3D images the new T₁₀ estimation is used to provide a motion corrected T₁₀ map.

21 . The system of any of the foregoing claims, wherein the acquired

images are of tissue of a subject and the final T_w estimate is used for pharmacokinetic analysis, tissue characterization and/or image segmentation.

22. The system of any of the foregoing claims, wherein the image

sequence Includes three or more MR images and the transformation model is applied groupwise to the images.