A Network Flow Method for Improved MR Field Map Estimation in the Presence of Water and Fat D. Hernando, P. Kellman, J. P. Haldar and Z.-P. Liang Abstract— Field map estimation is an important problem in MRI, with applications such as water/fat separation and cor- rection of fast acquisitions. However, it constitutes a nonlinear and severely ill-posed problem requiring regularization. In this paper, we introduce an improved method for regularized field map estimation, based on a statistically motivated formulation, as well as a novel algorithm for the solution of the corresponding optimization problem using a network flow approach. The proposed method provides theoretical guarantees (local optimality with respect to a large move), as well as an efficient implementation. It has been applied to the water/fat separation problem and tested on a number of challenging datasets, showing high-quality results. I. I NTRODUCTION In MRI, a very homogeneous main (B 0 ) magnetic field is desirable. However, inhomogeneities in the B 0 field are of- ten unavoidable, due to susceptibility differences introduced by the object being imaged, as well as magnet imperfec- tions. These inhomogeneities introduce undesired, spatially- varying phase shifts in the MR signal, which can be corrected given knowledge of the true B 0 field. Hence, estimation of the B 0 field inhomogeneity map (or “field map”) is an important problem in MRI, as it allows, e.g., effective water/fat separation, correction of EPI/spiral acquisitions, and automated shimming [1]–[3]. The field map can be estimated based on the phase evolu- tion of a sequence of images acquired at different echo times, t 1 , t 2 ,..., t N . In this work, we consider the presence of signal originating from water and fat, which further complicates the problem, as these two components have different phase behavior [4]. The signal at an individual voxel q can therefore be modeled as: s q (t n )= e i2π f B t n ( ρ W + ρ F e i2π f F t n ) (1) where t n is the echo time shift, f B (in Hz) is the local frequency shift due to B 0 field inhomogeneity, ρ W and ρ F are the intensities of the water and fat components, respectively, and f F (in Hz) is the frequency shift of fat, which is assumed known a priori. Field map estimation is a difficult problem due to the nonlinearity of the signal model and the presence of phase This work was supported in part by the following research grants: NIH- P41-EB03631-16 and NIH-R01-CA098717. D. Hernando, J. P. Haldar and Z.-P. Liang are with the Depart- ment of Electrical and Computer Engineering, University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA {dhernan2, haldar, z-liang} @uiuc.edu P. Kellman is with the Laboratory of Cardiac Energetics, National Heart, Lung and Blood Institute, National Institutes of Health, De- partment of Health and Human Services, Bethesda, MD 20892, USA kellmanp@nhlbi.nih.gov wraps (particularly in cases of high field inhomogeneity). To alleviate these problems, the estimated field map is typically regularized by imposing spatial smoothness. Most previously proposed methods resort to a two-step approach for estimat- ing the regularized field map: first, f B is estimated voxel- by-voxel using a maximum-likelihood (ML) criterion, and second, the resulting (noisy) field map is low-pass filtered to achieve the desired smoothness [5]. The main drawback of this method is that, while the low-pass filtering is generally effective in removing small noise-related perturbations in the field map, it is unable to correct the large errors due to the ill-posedness of the voxel-by-voxel estimation problem. Several extensions have been proposed to improve the initial voxel-by-voxel estimation [6]. In Ref. [7], a method is developed for directly estimating the regularized field map, assuming the presence of only water (i.e., ρ F = 0 in Eq. 1). This method formulates the estimation as a penalized ML (PML) problem, which is solved iteratively using conjugate gradients, producing a locally optimal solution. In this paper, we introduce a novel method for regularized field map estimation in the presence of water and fat, based on a PML formulation and an improved iterative optimization algorithm consisting on mapping each step to an equivalent network flow problem on a suitable graph. II. METHODS A. Problem formulation The signal model in Eq. 1 contains three unknown param- eters: {ρ W , ρ F , f B }. Under the assumption of white additive Gaussian noise, the ML estimate for {ρ W , ρ F , f B } is obtained by minimizing the following cost function at each voxel q: R 0 (ρ W , ρ F , f B ; s q )= N ∑ n=1   s q (t n ) − e i2π f B t n ( ρ W +ρ F e i2π f F t n )  2 (2) where N is the number of different echo times employed (typically N = 3), s q (t n ) is the measured signal at voxel q and echo time t n , and s q =[s q (t 1 ) ··· s q (t N )] T . Minimizing R 0 (ρ W , ρ F , f B ; s q ) is a separable nonlinear least- squares (NLLS) problem. As shown in [8], estimation of f B can be isolated using the variable projection (VARPRO) formulation, reducing to the minimization of: R( f B ; s q )=    I − Ψ( f B )Ψ † ( f B )  s q   2 2 (3)