System and method for model-based reconstruction of quantitative images

Abstract

A system and method for estimating a physiological parameter from data acquired with a medical imaging system includes acquiring data with the medical imaging system. A physiological parameter is estimated from the acquired data using an iterative estimation in which a model of the medical imaging system is decoupled from a physics-based model of the acquired data.

Claims

1. A system for generating medical images comprising: a communications connection configured to receive image data acquired from a subject using a medical imaging system; a processor configured to receive the image data and to: a) estimate a physiological parameter from the image data by inputting the image data to an estimation model comprising a model of the medical imaging system that is decoupled from a physics-based model of the image data, generating output as the estimate of the physiological parameter, wherein the physics-based model comprises a nonlinear signal model functional that models the image data such that the image data have a nonlinear dependence on at least one parameter of the physics-based signal model, and the model of the medical imaging system is actively decoupled from the physics-based model by applying a multi-level splitting routine to the model of the medical imaging system and the physics-based model; b) iteratively minimize a cost function that includes the model of the medical imaging system and the physics-based model of the acquired data to quantify the physiological parameter of the subject from the image data; and c) reconstruct a set of medical images of the subject from the image data at least including the physiological parameter of the subject.

2. The system of claim 1 wherein the model of the medical imaging system is a large-scale system operator that describes data acquisition using the medical imaging system.

3. The system of claim 2 wherein the medical imaging system is a magnetic resonance imaging (MRI) system and the model of the medical imaging system includes a Fourier transform.

4. The system of claim 1 wherein to estimate the physiological parameter the processor is further caused to apply the multi-level splitting routine in order to prospectively enable the use of a dimensionality reduction strategy to reduce overall complexity of the estimation due to non-uniqueness of data.

5. The system of claim 1 wherein to estimate the physiological parameter the processor is further caused to form separable nonlinear least squares problems.

6. The system of claim 1 wherein to iteratively minimize the cost function, the processor is further caused to control an influence of ancillary parameters, f, that are not associated with the physiological parameter.

7. The system of claim 6 wherein to control an influence of ancillary parameters, f, the processor is further caused to marginalize f out of estimations by constructing a specific alternating direction method of multipliers (ADMM) paradigm.

8. The system of claim 7 wherein to iteratively minimize the cost function, the processor is further caused to solve a linear system of: ${}_{\stackrel{\_}{u}}{J}_{\mathrm{AL}}=A^{*}\left(A\hat{u}-g\right)+{}_u\left(\hat{u}-v-{}_u\right)=0\hat{u}={\left[A^{*}A+{}_{u}I\right]}^{-1}\left[A^{*}g+{}_u\left(v+{}_u\right)\right]$ where J.sub.AL is quadratic with respect to , A is the system model, g is the image data, u and v are dummy variables that extended cost functions, .sub.u,.sub.v>0, and wherein .sub.u and .sub.v are Lagrange multiplier vectors.

9. The system of claim 2 wherein the medical imaging system is a magnetic resonance imaging (MRI) system and the physics-based model includes a magnetic resonance signal model.

10. A method for estimating a physiological parameter from data acquired with a medical imaging system, the steps of the method comprising: a) acquiring data with the medical imaging system; b) estimating a physiological parameter from the acquired data by inputting the acquired data to an estimation model comprising a model of the medical imaging system that is decoupled from a physics-based model of the acquired data, generating output as the estimate of the physiological parameter, wherein the physics-based model comprises a nonlinear signal model functional that models the acquired data such that the acquired data have a nonlinear dependence on at least one parameter of the physics-based signal model, and the model of the medical imaging system is actively decoupled from the physics-based model by applying a multi-level splitting routine to the model of the medical imaging system and the physics-based model; and c) reconstructing a set of medical images from the data, the medical images at least including the physiological parameter.

11. The method as recited in claim 10 wherein step b) includes iteratively minimizing a cost function that includes the model of the medical imaging system and the physics-based model of the acquired data.

12. The method as recited in claim 11 wherein step b) includes marginalizing unwanted variables out of the cost functional.

13. The method as recited in claim 12 wherein step b) in which one level of the multi-level splitting routine includes decoupling the model of the medical imaging system from the physics-based model of the acquired data.

14. The method as recited in claim 10 wherein step b) includes using an alternating direction method of multipliers routine to decouple the model of the medical imaging system from the physics-based model of the acquired data.

15. The method as recited in claim 10 wherein the medical imaging system is a magnetic resonance imaging (MRI) system, the model of the medical imaging system includes a Fourier transform, and the physics-based model of the acquired data is a magnetic resonance signal model.

16. A magnetic resonance imaging (MRI) system, comprising: a magnet system configured to generate a polarizing magnetic field about at least a portion of a subject arranged in the MRI system; a magnetic gradient system including a plurality of magnetic gradient coils configured to apply at least one magnetic gradient field to the polarizing magnetic field; a radio frequency (RF) system configured to apply an RF field to the subject and to receive magnetic resonance signals from the subject using a coil array; a computer system processor programmed to: control the magnetic gradient system and the RF system to acquire image data; estimate a physiological parameter from the image data by inputting the image data to an estimation model comprising a model of the MRI system that is decoupled from a physics-based model of the image data, generating output as the estimate of the physiological parameter, wherein the physics-based model comprises a nonlinear signal model functional that models the image data such that the image data have a nonlinear dependence on at least one parameter of the physics-based signal model, and the model of the medical imaging system is actively decoupled from the physics-based model by applying a multi-level splitting routine to the model of the medical imaging system and the physics-based model; and generate a quantitative image of the subject with the physiological parameter.

17. The MRI system of claim 16 wherein to estimate the physiological parameter the processor is further caused to apply the multi-level splitting routine in order to prospectively enable the use of a dimensionality reduction strategy to reduce overall complexity of the estimation due to non-uniqueness of data.

18. The MRI system of claim 16 wherein to estimate the physiological parameter the processor is further caused to form separable nonlinear least squares problems.

19. The MRI system of claim 16, further comprising iteratively minimizing a cost function, the processor programmed to control an influence of ancillary parameters, f, that are not associated with the physiological parameter.

20. The MRI system of claim 19 wherein to control an influence of ancillary parameters, f, the processor is further caused to marginalize f out of estimations by constructing a specific alternating direction method of multipliers (ADMM) paradigm.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1 is a block diagram of an example of a magnetic resonance imaging (MRI) system.

(2) FIG. 2 is a flow chart setting forth examples of steps that may be performed in accordance with the present disclosure.

DETAILED DESCRIPTION

(3) Described here are systems and methods for that provide a numerical strategy for reconstructing images from limited amounts of data, such as reconstructing images from magnetic resonance imaging (MRI) data acquired using an accelerated data acquisition that undersamples k-space. Although some aspect of the present disclosure will discuss an application for quantitative MRI, the systems and methods described herein can be applied to other technological areas both within and outside of medical imaging.

(4) MRI is widely used to diagnostically investigate the human body. However, unlike other medical imaging modalities (e.g., x-ray computed tomography), MRI data is not explicitly quantitative. That is, pixel values do not directly represent physiological information. There is strong clinical interest in transforming MRI into a quantitative modality, and to measure quantities such as flow, diffusion, perfusion, T1 relaxation times, and fat/water concentrations. For these applications, series of images are collected (e.g., at different acquisition settings), the acquired data is fit to physics-based models, and relevant physiological parameters are extracted. During this process, however, MRI data typically passes through a complex processing pipeline. Although errors in raw MRI data are well-understood, how they propagate through such a pipeline is not, and current results on these applications are typically biased and error-prone.

(5) The systems and methods of the present disclosure overcome these drawbacks by providing a generally applicable framework that directly can be used to estimate physiological parameters from the acquired imaging data, even if this data is incomplete (i.e., undersampled), thereby enabling statistically optimal error handling and increased accuracy. Given the strong statistical motivation behind this framework, as well as robust handling of challenging mathematical issues like signal bias and non-uniqueness, it is contemplated that the systems and methods described here will facilitate bringing quantitative MRI techniques into routine clinical use. This framework, and the computational methods that practically apply it, is a significant advance towards making MRI a clinically-useful quantitative diagnostic tool.

(6) As will be described, the systems and methods of the present disclosure may be used with a variety of different imaging modalities. As one non-limiting example, a magnetic resonance imaging (MRI) system will be described. Referring particularly now to FIG. 1, an example of an MRI system 100 is illustrated. The MRI system 100 includes an operator workstation 102, which will typically include a display 104; one or more input devices 106, such as a keyboard and mouse; and a processor 108. The processor 108 may include a commercially available programmable machine running a commercially available operating system. The operator workstation 102 provides the operator interface that enables scan prescriptions to be entered into the MRI system 100. In general, the operator workstation 102 may be coupled to four servers: a pulse sequence server 110; a data acquisition server 112; a data processing server 114; and a data store server 116. The operator workstation 102 and each server 110, 112, 114, and 116 are connected to communicate with each other. To this end, the operator workstation 102 (or, as will be described, other workstations 142) are connected to receive data acquired by the MRI system 100. For example, the servers 110, 112, 114, and 116 may be connected via a communication system 140, which may include any suitable network connection, whether wired, wireless, or a combination of both. As an example, the communication system 140 may include both proprietary or dedicated networks, as well as open networks, such as the internet.

(7) The pulse sequence server 110 functions in response to instructions downloaded from the operator workstation 102 to operate a gradient system 118 and a radiofrequency (RF) system 120. Gradient waveforms necessary to perform the prescribed scan are produced and applied to the gradient system 118, which excites gradient coils in an assembly 122 to produce the magnetic field gradients G.sub.x, G.sub.y, and G.sub.z used for position encoding magnetic resonance signals. The gradient coil assembly 122 forms part of a magnet assembly 124 that includes a polarizing magnet 126 and a whole-body RF coil 128.

(8) RF waveforms are applied by the RF system 120 to the RF coil 128, or a separate local coil (not shown in FIG. 1), in order to perform the prescribed magnetic resonance pulse sequence. Responsive magnetic resonance signals detected by the RF coil 128, or a separate local coil (not shown in FIG. 1), are received by the RF system 120, where they are amplified, demodulated, filtered, and digitized under direction of commands produced by the pulse sequence server 110. The RF system 120 includes an RF transmitter for producing a wide variety of RF pulses used in MRI pulse sequences. The RF transmitter is responsive to the scan prescription and direction from the pulse sequence server 110 to produce RF pulses of the desired frequency, phase, and pulse amplitude waveform. The generated RF pulses may be applied to the whole-body RF coil 128 or to one or more local coils or coil arrays (not shown in FIG. 1).

(9) The RF system 120 also includes one or more RF receiver channels. Each RF receiver channel includes an RF preamplifier that amplifies the magnetic resonance signal received by the coil 128 to which it is connected, and a detector that detects and digitizes the I and Q quadrature components of the received magnetic resonance signal. The magnitude of the received magnetic resonance signal may, therefore, be determined at any sampled point by the square root of the sum of the squares of the I and Q components:
M={square root over (I.sup.2+Q.sup.2)}(1);

(10) and the phase of the received magnetic resonance signal may also be determined according to the following relationship:

(11) $\begin{array}{cc}={\tan }^{-1}\left(\frac{Q}{I}\right).& \left(2\right)\end{array}$

(12) The pulse sequence server 110 also optionally receives patient data from a physiological acquisition controller 130. By way of example, the physiological acquisition controller 130 may receive signals from a number of different sensors connected to the patient, such as electrocardiograph (ECG) signals from electrodes, or respiratory signals from a respiratory bellows or other respiratory monitoring device. Such signals are typically used by the pulse sequence server 110 to synchronize, or gate, the performance of the scan with the subject's heart beat or respiration.

(13) The pulse sequence server 110 also connects to a scan room interface circuit 132 that receives signals from various sensors associated with the condition of the patient and the magnet system. It is also through the scan room interface circuit 132 that a patient positioning system 134 receives commands to move the patient to desired positions during the scan.

(14) The digitized magnetic resonance signal samples produced by the RF system 120 are received by the data acquisition server 112. The data acquisition server 112 operates in response to instructions downloaded from the operator workstation 102 to receive the real-time magnetic resonance data and provide buffer storage, such that no data is lost by data overrun. In some scans, the data acquisition server 112 does little more than pass the acquired magnetic resonance data to the data processor server 114. However, in scans that require information derived from acquired magnetic resonance data to control the further performance of the scan, the data acquisition server 112 is programmed to produce such information and convey it to the pulse sequence server 110. For example, during prescans, magnetic resonance data is acquired and used to calibrate the pulse sequence performed by the pulse sequence server 110. As another example, navigator signals may be acquired and used to adjust the operating parameters of the RF system 120 or the gradient system 118, or to control the view order in which k-space is sampled. In still another example, the data acquisition server 112 may also be employed to process magnetic resonance signals used to detect the arrival of a contrast agent in a magnetic resonance angiography (MRA) scan. By way of example, the data acquisition server 112 acquires magnetic resonance data and processes it in real-time to produce information that is used to control the scan.

(15) The data processing server 114 receives magnetic resonance data from the data acquisition server 112 and processes it in accordance with instructions downloaded from the operator workstation 102. Such processing may, for example, include one or more of the following: reconstructing two-dimensional or threedimensional images by performing a Fourier transformation of raw k-space data; performing other image reconstruction algorithms, such as iterative or backprojection reconstruction algorithms; applying filters to raw k-space data or to reconstructed images; generating functional magnetic resonance images; calculating motion or flow images; and so on.

(16) Images reconstructed by the data processing server 114 are conveyed back to the operator workstation 102 where they are stored. Real-time images are stored in a data base memory cache (not shown in FIG. 1), from which they may be output to operator display 112 or a display 136 that is located near the magnet assembly 124 for use by attending physicians. Batch mode images or selected real time images are stored in a host database on disc storage 138. When such images have been reconstructed and transferred to storage, the data processing server 114 notifies the data store server 116 on the operator workstation 102. The operator workstation 102 may be used by an operator to archive the images, produce films, or send the images via a network to other facilities.

(17) The MRI system 100 may also include one or more networked workstations 142. By way of example, a networked workstation 142 may include a display 144; one or more input devices 146, such as a keyboard and mouse; and a processor 148. The networked workstation 142 may be located within the same facility as the operator workstation 102, or in a different facility, such as a different healthcare institution or clinic.

(18) The networked workstation 142, whether within the same facility or in a different facility as the operator workstation 102, may gain remote access to the data processing server 114 or data store server 116 via the communication system 140. Accordingly, multiple networked workstations 142 may have access to the data processing server 114 and the data store server 116. In this manner, magnetic resonance data, reconstructed images, or other data may exchanged between the data processing server 114 or the data store server 116 and the networked workstations 142, such that the data or images may be remotely processed by a networked workstation 142. This data may be exchanged in any suitable format, such as in accordance with the transmission control protocol (TCP), the internet protocol (IP), or other known or suitable protocols.

(19) Referring to FIG. 2, a process 200 in accordance with the present disclosure begins at process block 202 by acquiring medical imaging data. Thereafter, at process block 204, a model for the data is selected. For example, using a medical imaging system, such as the above-described MRI system 100, an MRI signal model is selected. Suppose T imaging studies are performed, with each study corresponding to a different time point or, in the case of an MRI acquisition, a pulse sequence parameterization. During each study or acquisition, the C receiver channels each collect K data samples. This KTC data set can be modeled as:

(20) $\begin{array}{cc}g\left[k,t,c\right]=\underset{n=0}{\overset{N-1}{.Math.}}A\left[k,t,c\right]f\left[n,c\right]F_{}\left[n,t\right]+z\left[k,t,c\right];& \left(3\right)\end{array}$

(21) where k is the observed signal sample index; t is the experiment index; c is the receiver channel index; A is a KNT spatial-spectral system model; f is an NC matrix of ancillary parameters; is an N1 vector or NT matrix of physiological parameters of interest; F.sub. is an NT physics-based signal model; and z is complex Gaussian noise. In magnetic resonance imaging applications, the system model, A, often corresponds to a partial discrete Fourier transform (DFT) as:
A[k,n,t]=e.sup.jx[n].Math.[k,t](4).

(22) Additionally, it can be assumed that the corrupting noise process is described by {z[k,t,c]}, {z[k,t,c]}N(0,.sup.2). The ancillary parameters, f, represent quantities that are not of primary interest. In quantitative MRI applications, for example, this may include coil sensitivity profiles used in accelerated imaging studies. As will be shown, as desired, these parameters can be marginalized out of the estimation problem. The physics-based signal model, F.sub., selected at process block 204 can be, for example, either a linear or nonlinear functional. For example, in the case of an MRI data acquisition using a spoiled gradient recalled echo (SPGR) pulse sequence, R.sub.1 relaxometry, is an N1 vector representing the relaxation values for different tissues and:

(23) $\begin{array}{cc}{F}_{}\left[n,t\right]=\sin \left(\left[t\right]\right)\frac{1-e^{-\left[n\right]}}{1-\cos \left(\left[t\right]\right)e^{-\left[n\right]}};& \left(5\right)\end{array}$

(24) where and are pulse sequence parameters.

(25) Thus, parameter estimation can begin. The signal model in Eqn. (3) can be equivalently expressed as a KTC matrix:

(26) $\begin{array}{cc}g=\left[\begin{array}{cccc}A\left[:,:,0\right]& 0& \mathrm{.Math.}& 0\\ 0& & & \mathrm{.Math.}\\ \mathrm{.Math.}& & & 0\\ 0& \mathrm{.Math.}& 0& A\left[:,:,T-1\right]\end{array}\right]\left[\begin{array}{c}\mathrm{diag}\left(F_{}\left[;,0\right]\right)\\ \mathrm{.Math.}\\ \mathrm{diag}\left(F_{}\left[:,T-1\right]\right)\end{array}\right]f+z=\mathrm{AG}\left(\right)f+z.& \left(6\right)\end{array}$

(27) Because z is Gaussian, maximum likelihood (ML) estimates of f and , which are both unknown parameters, are solutions to a separable nonlinear least squares problem. In practice, however, K<N, and the system in Eqn. (6) is underdetermined.

(28) At process block 206, the model selected at process block 204 is used to generate a quantitative estimate. As will be described, estimating quantitative estimates of physiological parameters from the acquired data can be performed using an iterative estimation in which a model of the medical imaging system is decoupled from a physics-based model of the acquired data. More particularly, the process can iteratively minimize a cost functional that includes the model of the medical imaging system and the physics-based model of the acquired data.

(29) In particular, an initial estimate may be derived from the acquired data and, for this estimation to be well-posed, auxiliary regularization or constraints may be used. To this end and continuing with the non-limiting example of MRI data, the following regularized nonlinear least squares model can be employed for physiological parameter estimation in undersampled quantitative MRI:
J(f,)=P()+AG()fg.sub.F.sup.2(7);

(30) where P(.Math.) is a penalty functional, which may not necessarily be smooth, and .Math..sub.F denotes the Frobenius matrix norm. Because f is ultimately not of interest, it can be left unregularized.

(31) Consider the special case of Eqn. (7) where A is the identity operator, where:
J(f,)=P()+G()fg.sub.F.sup.2(8).

(32) Although f and can be jointly estimated using standard optimization strategies, like a nonlinear conjugate gradient iteration, such strategies are often unstable, inefficient, or both. Moreover, much effort is spent estimating the nuisance variable, f, which is often subsequently discarded. Alternatively, because J(f,) is quadratic with respect to f, a closed-form expression for the minimizing value of this variable, as a function of , exists and can be used to marginalize out this variable from the cost functional.

(33) This approach can be referred to as the variable projection (VARPRO) strategy for separable nonlinear least squares problems. Specifically, for Eqn. (8), the optimizing expression for f is:
f=G.sup.()g(9);

(34) where (.).sup. denotes the left pseudo-inverse operator. Plugging this expression back into Eqn. (8) yields the following:
J()=P()+H()g.sub.F.sup.2(10);
where,
H()=IG()G.sup.()(11).

(35) For many quantitative problems, P() is a convex Markov penalty. In such cases, the single-variable problem in Eqn. (10) can be approximately minimized using discrete optimization strategies, such as graph cut techniques. For the general quantitative MRI problem, A is not the identity operator and VARPRO unfortunately cannot be used directly. As will be described, however, f can still be marginalized out of the estimation problem by constructing a specific alternating direction method of multipliers (ADMM) paradigm.

(36) It is noted that VARPRO cannot be easily used for the general problem in Eqn. (7) because of the action of A on G(). If these two functions can be separated, then a VARPRO-like operation can potentially be applied to the sub-problem involving only G(). To this end, instead of jointly minimizing Eqn. (7) over f and , the following equivalent constrained optimization problem can be considered:

(37) $\begin{array}{cc}\left[\hat{f},\hat{},\hat{u},\hat{v}\right]=\arg \underset{f,,u,v}{\min }{J}_{\mathrm{ext}}\left(f,,u,v\right)\mathrm{such}\mathrm{that}u=v=G\left(\right)f;& \left(12\right)\end{array}$

(38) where u and v are dummy variables and the extended cost functional is:
J.sub.ext(f,,u,v)=P()+Aug.sub.F.sup.2(13).

(39) Although standard ADMM uses only a single operator splitting, here a second level of splitting is employed. This second level of splitting allows the preemptive decoupling between A and G () during the VARPRO step. In other words, knowing that a step backwards will eventually be taken, two steps forward are taken up front so as to ensure that the algorithm is still only one step ahead in the end. The augmented Lagrangian (AL) functional for Eqn. (12) is:
J.sub.AL(f,,u,v,.sub.u,.sub.v)=P()+AUg.sub.F.sup.2+.sub.uuv.sub.u.sub.F.sup.2+.sub.vvG()f.sub.v.sub.F.sup.2 (14);

(40) with .sub.u,.sub.v>0, and wherein .sub.u and .sub.v are Lagrange multiplier vectors. An advantage of this relaxation is that only G() acts on f, and it can now be marginalized out using VARPRO,
J.sub.AL(,u,v,.sub.u,.sub.v)=P()+Aug.sub.F.sup.2+.sub.uuv.sub.vH()(v.sub.v).sub.F.sup.2(15).

(41) This problem can now be approached using an ADMM technique, which entails alternating minimization of J over u, v, and .

(42) J.sub.AL is quadratic with respect to u, and so the minimizing value can be found by solving the following linear system:

(43) $\begin{array}{cc}{}_{\stackrel{\_}{u}}{J}_{\mathrm{AL}}=A^{*}\left(A\hat{u}-g\right)+{}_u\left(\hat{u}-v-{}_u\right)=0\hat{u}={\left[A^{*}A+{}_{u}I\right]}^{-1}\left[A^{*}g+{}_u\left(v+{}_u\right)\right].& \left(16\right)\end{array}$

(44) In many cases, A is diagonalizable (e.g., via FFT) and so this problem can be readily solved. In non-ideal cases, such as in non-Cartesian sampling patters when using MRI data, this positive-definite linear system can be solved via conjugate gradient (CG) iteration. Thus, at decision block 208, a check is made against a threshold, which may be a boundary condition or optimization criteria. As for the base case of A being equal to the identity operator, approximate (i.e., inexact) minimization of J.sub.AL with respect to can be achieved using a discrete optimization method like graph cuts. Since f is defined as a function of in this setup, it is also implicitly regularized by constraints imposed onto . This leaves only the estimation of v. Like u, v can be determined by solving the linear system corresponding to .sub.vJ.sub.AL=0:

(45) $\begin{array}{cc}\hat{v}={\left[H^{*}\left(\right)H\left(\right)+\frac{{}_u}{{}_v}I\right]}^{-1}\left[H^{*}\left(\right)H\left(\right){}_v+\frac{{}_u}{{}_v}\left(u-{}_u\right)\right];& \left(17\right)\end{array}$

(46) which can be stated equivalently as:

(47) $\begin{array}{cc}\hat{v}={\left[H\left(\right)+\frac{{}_u}{{}_v}I\right]}^{-1}\left[H\left(\right){}_v+\frac{{}_u}{{}_v}\left(u-{}_u\right)\right];& \left(18\right)\end{array}$

(48) because H() is both self-adjoint and idempotent. The complexity of this problem is determined by the inversion of:

(49) $\begin{array}{cc}{\left[H\left(\right)+\left(-1\right)I\right]}^{-1}=\left[I-{G\left(\right)\left[G^{*}\left(\right)G\left(\right)\right]}^{-1}{G}^{*}\left(\right)\right];\mathrm{where}\text{:}& \left(19\right)\\ =\left(1+\frac{{}_u}{{}_v}\right).& \left(20\right)\end{array}$

(50) However, Woodbury's matrix identity asserts that:

(51) 0$\begin{array}{cc}{\left[H\left(\right)+\left(-1\right)I\right]}^{-1}={}^{-1}\left[\left(\frac{}{-1}\right)I-\left(\frac{1}{-1}\right)H\left(\right)\right];& \left(21\right)\end{array}$

(52) which does not contain an outer matrix inverse. Moreover, the matrix implicitly inverted inside of H() is diagonal and, thus, also trivially handled. Therefore, v can also be easily estimated.

(53) Putting all of these pieces together in the ADMM framework yields the following iterative scheme:

(54) TABLE-US-00001 input: g , maxIter ; define: k = 0, v = G(0)f , u = v , = 0, .sub.u = 0, .sub.v = 0, .sub.u > 0, .sub.v > 0; repeat | u = [A*A + .sub.uI].sup.1 [A*g + .sub.u (v + .sub.u)]; | $v={}^{-1}\left[\left(\frac{}{-1}\right)I-\left(\frac{1}{-1}\right)H\left(\right)\right]\left[H\left(\right){}_v+\frac{{}_u}{{}_v}\left(u-{}_u\right)\right];$ | = arg min{P() + .sub.vH()(v .sub.v).sub.F.sup.2}; | .sub.u = .sub.u (u v); | .sub.v = .sub.v (v G()f); | k = k + 1 k until k maxIter ;

(55) A review of the framework provided above reveals that it is of the classic ADMM form, where at each iteration a succession of independent variable estimations are performed followed by updating of the Lagrange multiplier vectors.

(56) By using the multi-stage splitting approach described above, the large-scale system operator, A, and the nonlinear signal model functional, H (), are actively decoupled and, thus, have no complex interactions. Moreover, the third sub-problem of the iterative routine is similar to the classic, non-accelerated NLLS problem for which robust and efficient computational routines already exist.

(57) Once the iterative process is completed such that the threshold at decision block 208 is met, at process block 210 the quantitative estimates can be registered with anatomical information from the data acquired at process block 202 or other anatomical information. At process block 212, a report is communicated, for example, as a set of images with registered quantitative estimates coupled with anatomical images. For example, the quantitative estimates may be overlaid on anatomical images, such as using a color coding or other display or report mechanism.

(58) Thus, the present disclosure may utilize a specific ADMM strategy, which is a form of augmented Lagrangian routine, to actively decouple the challenging components of an estimation model so they can be focused on independently. More specifically, the ADMM strategy makes use of a multi-level splitting routine that prospectively enables the use of a dimensionality reduction strategy called variable projection (VARPRO). This approach greatly reduces the overall complexity of the estimation problem while simultaneously providing a pathway to circumvent mathematical challenges associated with non-uniqueness problems, like phase wrapping.

(59) The present disclosure includes discussion of one or more preferred embodiments, and it should be appreciated that many equivalents, alternatives, variations, and modifications, aside from those expressly stated, are possible and within the scope of the invention.