1 Introduction
Markov random fields, specifically in the form of the Ishikawa model[1]
, have previously been used to address the issue of reconstructing an image from noisy data using probabilistic priors and maximum a posteriori probability optimization. The benefit of such frameworks is their submodularity and therefore global optimizability via graph cuts.
[2] The primary improvement of Ishikawa models over previous graph cut approaches in reconstructing an image[3] or depth field[4] is that the topology of the intensity space, an ordered linear topology, is fundamentally built into the functional.Continuous maxflow approaches have traditionally been seen as continuous analogs to earlier discrete models. Along these lines, a continuous Ishikawa model [5] was proposed which maintains the same linear topology over the intensity space but using a continuous, rather than discrete domain. This is especially crucial in image reconstruction as the mitigation of discrete metrification artifacts achieved by the continuous formulation [6] has an immediate visual effect on the reconstruction.
Certain types of fields, however, do not have ordered linear intensities. For example, hue in HSL and HSV colour spaces or phase in a complexvalued image represents a cyclic, rather than linear, topology. These require a specialized reconstruction optimizer which, as the Ishikawa model was to linear intensity models, sensitive to the intensity’s cyclic topology. Although this could in theory be accomplished through a directed acyclic graph continous maxflow (DAGMF) model as described by Baxter et al.,[7] the quadratic growth of said models is undesirable, especially for medical applications such as phase processing in susceptibility weighted magnetic resonance imaging[8] in which high intensity resolution is necessary.
2 Contributions
The main contribution of this work is the development of a continuous maxflow model describing cyclic field image reconstruction, that is, reconstruction of an image in which the image intensity is defined over a cyclic topology rather than a linear topology as explored by Bae et al.[5]. This reconstruction algorithm is better suited for the reconstruction of cyclic fields such as colour/hue in natural images or phase in complexvalued magnetic resonance images.
3 Previous Work
Previous work by Yuan et al. [6] has addressed the continuous binary mincut problem:
as well as the convex relaxed continuous Potts Model[9]:
Both of these techniques both used a continuous maxflow model with augmented Lagrangian multipliers. In the case of the convexrelaxed continuous Potts model, the source flow had infinite capacity and the costs in the functional corresponded with constraints on the sink flows.
Bae et al. [5] extended the work on the continuous binary mincut problem to the continuous Ishikawa model:
using similar variational methods but a tiered continuous graph analogous to that used by Ishikawa [1] in the discrete case, that is, with finite capacities on intermediate flows between labels.
Both continuous Potts and Ishikawa models are important to consider as they were the first maxflow models to be used in image reconstruction. However, they are not topologically suitable for cyclic fields because the Potts model disallows any notion of label topology, and the Ishikawa model is defined only over linear orderings. Additionally, cyclic orderings are nonhierarchical, making them inaccessible to generalized hierarchical maxflow solvers[10] and directed acyclic graph maxflow models[7] of these orderings require quadratic time to handle linear increases in the intensity resolution.
4 Cyclic Topology in Continuous MaxFlow Segmentation
4.1 Continuous MaxFlow Model
The desired optimization model for cyclic field reconstruction would be:
(1) 
in which is defined over a cyclic manifold to prevent modulation errors.
One can reexpress the reconstructed intensity with respect to a series of fuzzy indicator functions:
Using these indicator functions, the data term can be rewritten as:
and the smoothness term can be rewritten as:
yielding the alternative minimization formulation:
(2)  
The benefit of this new formulation via indicator functions is that it resembles the previously explored continuous Potts model[9] with the exception that the gradient magnitude operator is defined over the Cartesian product of the spatial domain with a cyclic intensity domain, . Thus, the spatial domain augmented with the direction forms a high dimensional cylinder as shown in Figure 1.
4.2 Primal Formulation
As in the previous work, the primal formulation will be expressed as a maxflow optimization, specifically:
(3) 
subject to the capacity constraints:
and the flow conservation constraint:
, the spatial flow now has an additional dimension, the direction. The direction is equipped with a cyclic topology and therefore is not susceptible to modulation artifacts. Similar to , makes use of the spatial domain augmented with the cyclic domain of as shown in Figure 1.
4.3 PrimalDual Formulation
The primal model can be converted to a primaldual model through the use of Lagrangian multipliers on the flow conservation constraint . This yields the equation:
(4)  
To ensure that this function meets the criteria of the minimax theorem, we must ensure that it is convex with respect to , considering to be fixed, and concave with respect to with fixed. [11] For the first, it is sufficient to note that if are fixed, then is fixed as well, meaning that (4) is linear and therefore convex with respect to . It should also be noted that is a linear function of , implying the linearity (and concavity) of (4) with respect to . Thus, a saddle point must exist and the formulations are equivalent. [11]
4.4 Dual Formulation
To show the equivalence of the primaldual formulation (Eq. (4)) and the desired dual formulation (Eq. (1)), we must show that the original formula can be reconstructed from the definition of the saddlepoint. First, consider the primaldual formulation with the sink flows, isolated:
(5) 
which implies that as if , the minimization is unbounded as . This reconstructs the data term portion of the original formulation using the fuzzy labeling function.
The source flow, in Eq. 4 can be isolated as:
(6) 
at the saddle point defined by . This reconstructs the labeling constraint, that is, that only one intensity value is assigned to each voxel.
Considering coupling the spatial domain with the cyclic domain of via a simple Cartesian product with members taking the form . In this coupled space, the maximization of the spatial flow functions can be expressed [12] as:
(7)  
This yields the smoothness component in Eq. (1) thus showing the equivalence of the two models.
4.5 Solution to PrimalDual Formulation
The optimization problem being addressed is the augmented form of Eq. (4):
(8)  
in which is a nonnegative constant which has no effect on Eq. (4) when the constraint holds, but serves as an additional penalty. Optimizing Eq. (8) can be achieved by performing the following steps iteratively, optimizing Eq. (8) with respect to each individual variable:
The specific augmented Lagrangian based algorithm used is:
4.6 PseudoFlow Formulation
The approach taken here for the use of proximal Bregman or pseudoflow methods is derived from Baxter et al.[14]. The corresponding nonsmooth pseudoflow model can be written as:
(9) 
which can be derived from Eq (4) by:
Now that a pseudoflow representation is developed, we can take advantage of Bregman proximal projections[15] to optimize this formula. Consider the distance between labeling functions and as:
(10) 
which can be verified to be a Bregman distance (when ) using the entropy function:
If we consider a feasible labeling, , we can find another proximal labeling, , which has a lower energy by addressing the optimization:
(11) 
where is a positive constant. Using a Lagrangian multiplier on the constraint , we can solve for analytically as:
(12) 
noting that this answer fulfills the constraint provided the same holds for . By letting the distance weighting parameter approach 0, if . Using that fact,
illustrating that the Bregman proximal method is a smoothed version of the nonsmooth equation Eq (9).
4.7 Solution to PseudoFlow Formulation
By taking the gradient of equation (11) with respect to (with substituted by equation (12)) one can derive the appropriate Chambolle iteration scheme[13] for maximizing (11) with respect to the spatial flow variables :
(13) 
where is a positive gradient descent parameter. Thus, the proximal Bregman, or pseudoflow, algorithm relies on the alternation between the analytic optimization of the label functions, through Eq. (12) and the specified Chambolle iteration.
The proximal Bregman based algorithm proposed is therefore:
5 Discussion and Conclusions
In this paper, we present two algorithms for addressing continuous maxflow image reconstruction in which the intensity being reconstructed is fundamentally organized along a cyclic, as opposed to linear, topology. Such a reconstruction method is more efficient than that developed from previous continuous maxflow methods capable of incorporating said topology, specifically DAGMF[7], by displaying linear rather than quadratic growth in complexity in terms of the intensity resolution. This solver has been implemented using GPGPU acceleration due to its inherent parallelizability. This solver is available opensource for both two and threedimensional images at www.advancedsegmentationtools.org.
Future work in continuous maxflow based image reconstruction includes the robust automatic detection of gradient direction which may allow for anisotropic smoothness terms to be employed. These terms could be used to estimate level sets in the image, maintaining them while allowing for faster or lessconstrained gradients in the direction orthogonal to said level set.
Acknowledgements.
The authors would like to acknowledge Zahra Hosseini and Dr. Martin Rajchl for their invaluable discussion and editing.References
 [1] H. Ishikawa, “Exact optimization for markov random fields with convex priors,” Pattern Analysis and Machine Intelligence, IEEE Transactions on 25(10), pp. 1333–1336, 2003.
 [2] V. Kolmogorov and R. Zabin, “What energy functions can be minimized via graph cuts?,” Pattern Analysis and Machine Intelligence, IEEE Transactions on 26(2), pp. 147–159, 2004.
 [3] Y. Boykov, O. Veksler, and R. Zabih, “Fast approximate energy minimization via graph cuts,” Pattern Analysis and Machine Intelligence, IEEE Transactions on 23(11), pp. 1222–1239, 2001.
 [4] V. Kolmogorov and R. Zabih, “Multicamera scene reconstruction via graph cuts,” in Computer Vision ECCV 2002, pp. 82–96, Springer, 2002.
 [5] E. Bae, J. Yuan, X.C. Tai, and Y. Boykov, “A fast continuous maxflow approach to nonconvex multilabeling problems,” in Efficient Algorithms for Global Optimization Methods in Computer Vision, pp. 134–154, Springer, 2014.

[6]
J. Yuan, E. Bae, and X.C. Tai, “A study on continuous maxflow and mincut
approaches,” in
Computer Vision and Pattern Recognition (CVPR), 2010 IEEE Conference on
, pp. 2217–2224, IEEE, 2010.  [7] J. S. Baxter, M. Rajchl, J. Yuan, and T. M. Peters, “A continuous maxflow approach to multilabeling problems under arbitrary region regularization,” arXiv preprint arXiv:1405.0892 , 2014.
 [8] E. M. Haacke, Y. Xu, Y.C. N. Cheng, and J. R. Reichenbach, “Susceptibility weighted imaging (swi),” Magnetic Resonance in Medicine 52(3), pp. 612–618, 2004.
 [9] J. Yuan, E. Bae, X.C. Tai, and Y. Boykov, “A continuous maxflow approach to potts model,” in Computer Vision–ECCV 2010, pp. 379–392, Springer, 2010.
 [10] J. S. Baxter, M. Rajchl, J. Yuan, and T. M. Peters, “A continuous maxflow approach to general hierarchical multilabeling problems,” arXiv preprint arXiv:1404.0336 , 2014.
 [11] I. Ekeland and R. Temam, “Convex analysis and 9 variational problems,” 1976.
 [12] E. Giusti, Minimal surfaces and functions of bounded variation, no. 80, Springer Science & Business Media, 1984.
 [13] A. Chambolle, “An algorithm for total variation minimization and applications,” Journal of Mathematical imaging and vision 20(12), pp. 89–97, 2004.
 [14] J. S. Baxter, M. Rajchl, J. Yuan, and T. M. Peters, “A proximal bregman projection approach to continuous maxflow problems using entropic distances,” arXiv preprint arXiv:1501.07844 , 2015.
 [15] L. M. Bregman, “The relaxation method of finding the common point of convex sets and its application to the solution of problems in convex programming,” USSR computational mathematics and mathematical physics 7(3), pp. 200–217, 1967.
Comments
There are no comments yet.