A Continuous Max-Flow Approach to Cyclic Field Reconstruction

11/11/2015 ∙ by John S. H. Baxter, et al. ∙ 0

Reconstruction of an image from noisy data using Markov Random Field theory has been explored by both the graph-cuts and continuous max-flow community in the form of the Potts and Ishikawa models. However, neither model takes into account the particular cyclic topology of specific intensity types such as the hue in natural colour images, or the phase in complex valued MRI. This paper presents cyclic continuous max-flow image reconstruction which models the intensity being reconstructed as having a fundamentally cyclic topology. This model complements the Ishikawa model in that it is designed with image reconstruction in mind, having the topology of the intensity space inherent in the model while being readily extendable to an arbitrary intensity resolution.



There are no comments yet.


page 1

page 2

page 3

page 4

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

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 max-flow 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 complex-valued 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 max-flow (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 max-flow 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 complex-valued magnetic resonance images.

3 Previous Work

Previous work by Yuan et al. [6] has addressed the continuous binary min-cut problem:

as well as the convex relaxed continuous Potts Model[9]:

Both of these techniques both used a continuous max-flow model with augmented Lagrangian multipliers. In the case of the convex-relaxed 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 min-cut 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 max-flow 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 non-hierarchical, making them inaccessible to generalized hierarchical max-flow solvers[10] and directed acyclic graph max-flow models[7] of these orderings require quadratic time to handle linear increases in the intensity resolution.

4 Cyclic Topology in Continuous Max-Flow Segmentation

4.1 Continuous Max-Flow Model

The desired optimization model for cyclic field reconstruction would be:


in which is defined over a cyclic manifold to prevent modulation errors.

One can re-express the reconstructed intensity with respect to a series of fuzzy indicator functions:

Using these indicator functions, the data term can be re-written as:

and the smoothness term can be re-written as:

yielding the alternative minimization formulation:


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 max-flow optimization, specifically:


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.

Figure 1: Topology used in cyclic max-flow in which the spatial domain, is one dimensional. The resulting topology has one more dimensional than which encodes the intensity for the spatial flows, , and indicator functions, . The gradient magnitude and divergence operators over and include both the spatial and the dimensions.

4.3 Primal-Dual Formulation

The primal model can be converted to a primal-dual model through the use of Lagrangian multipliers on the flow conservation constraint . This yields the equation:


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 primal-dual 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 saddle-point. First, consider the primal-dual formulation with the sink flows, isolated:


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:


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:


This yields the smoothness component in Eq. (1) thus showing the equivalence of the two models.

4.5 Solution to Primal-Dual Formulation

The optimization problem being addressed is the augmented form of Eq. (4):


in which is a non-negative 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:

  1. Maximize (8) over :

    which is a Chambolle projection iteration with descent parameter .[13] Note that both the divergence and gradient operators are evaluated over both the linear spatial domain and the cyclic domain.

  2. Maximize (4) over the sink flows, , by:

  3. Maximize (4) over the source flow analytically by:

  4. Minimize (4) over analytically by:

The specific augmented Lagrangian based algorithm used is:

while not converged do
end while

4.6 Pseudo-Flow Formulation

The approach taken here for the use of proximal Bregman or pseudo-flow methods is derived from Baxter et al.[14]. The corresponding non-smooth pseudo-flow model can be written as:


which can be derived from Eq (4) by:

(as the equivalence of equations (1) and (4) guarantees the constraint )
(as the equivalence of equations (1) and (4) guarantees the constraint and
if is non-negative)
(since when
and )

Now that a pseudo-flow representation is developed, we can take advantage of Bregman proximal projections[15] to optimize this formula. Consider the distance between labeling functions and as:


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:


where is a positive constant. Using a Lagrangian multiplier on the constraint , we can solve for analytically as:


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 non-smooth equation Eq (9).

4.7 Solution to Pseudo-Flow 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 :


where is a positive gradient descent parameter. Thus, the proximal Bregman, or pseudo-flow, 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:

while not converged do
end while

5 Discussion and Conclusions

In this paper, we present two algorithms for addressing continuous max-flow 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 max-flow 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 open-source for both two- and three-dimensional images at www.advancedsegmentationtools.org.

Future work in continuous max-flow 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 less-constrained gradients in the direction orthogonal to said level set.

The authors would like to acknowledge Zahra Hosseini and Dr. Martin Rajchl for their invaluable discussion and editing.