Bethe Projections for Non-Local Inference

03/04/2015 ∙ by Luke Vilnis, et al. ∙ University of Massachusetts Amherst 0

Many inference problems in structured prediction are naturally solved by augmenting a tractable dependency structure with complex, non-local auxiliary objectives. This includes the mean field family of variational inference algorithms, soft- or hard-constrained inference using Lagrangian relaxation or linear programming, collective graphical models, and forms of semi-supervised learning such as posterior regularization. We present a method to discriminatively learn broad families of inference objectives, capturing powerful non-local statistics of the latent variables, while maintaining tractable and provably fast inference using non-Euclidean projected gradient descent with a distance-generating function given by the Bethe entropy. We demonstrate the performance and flexibility of our method by (1) extracting structured citations from research papers by learning soft global constraints, (2) achieving state-of-the-art results on a widely-used handwriting recognition task using a novel learned non-convex inference procedure, and (3) providing a fast and highly scalable algorithm for the challenging problem of inference in a collective graphical model applied to bird migration.



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

Structured prediction has shown great success in modeling problems with complex dependencies between output variables. Practitioners often use undirected graphical models, which encode conditional dependency relationships via a graph. However, the tractability of exact inference in these models is limited by the graph’s treewidth, often yielding a harsh tradeoff between model expressivity and tractability.

Graphical models are good at modeling local dependencies between variables, such as the importance of surrounding context in determining the meaning of words or phrases. However, their sensitivity to cyclic dependencies often renders them unsuitable for modeling preferences for certain globally consistent states. For example, in the canonical NLP task of part-of-speech tagging, there is no clear way to enforce the constraint that every sentence have at least one verb without increasing the likelihood that every token is predicted to be a verb.

Concretely, exact marginal inference in a discrete graphical model can be posed as the following optimization problem



is a concatenated vector of node and clique marginals,

is the entropy, is the marginal polytope, and are parameters. Here we face a tradeoff: adding long-range dependencies directly to the model increases the clique size and thus the complexity of the problem and size of , rendering inference intractable. However, the linear scoring function breaks down over cliques, preventing us from enforcing global regularities in any other way. In this work, we propose to augment the inference objective (1) and instead optimize


Here, is some arbitrary parametric function of the entire concatenated marginal vector, where may depend on input features. Since is non-linear, it can enforce many types of non-local properties. Interestingly, whenever is convex, and whenever inference is easy in the underlying model, i.e., solving (1) is tractable, we can solve (2) using non-Euclidean projected gradient methods using the Bethe entropy as a distance-generating function. Unlike many message-passing algorithms, our procedure maintains primal feasibility across iterations, allowing its use as an  anytime algorithm. Furthermore, for non-convex , we also show convergence to a local optimum of (2). Finally, we present algorithms for discriminative learning of the parameters . In a slight abuse of terminology, we call a non-local energy function.

Ours is not the first work to consider modeling global preferences by augmenting a tractable base inference objective with non-local terms. For example, generalized mean-field variational inference algorithms augment a tractable distribution (the distribution) with a non-linear, non-convex global energy function that scores terms in the full model (the distribution) using products of marginals of  (Wainwright & Jordan, 2008). This is one special case of our non-local inference framework, and we present algorithms for solving the problem for much more general , with compelling applications.

Additionally, the modeling utility provided by global preferences has motivated work in dual decomposition, where inference in loopy or globally-constrained models is decomposed into repeated calls to inference in tractable independent subproblems (Komodakis et al., 2007; Sontag et al., 2011). It has seen wide success due to its ease of implementation, since it reuses existing inference routines as black boxes. However, the technique is restricted to modeling linear constraints, imposed a priori. Similarly, these types of constraints have also been imposed on expectations of the posterior distribution for use in semi-supervised learning, as in posterior regularization and generalized expectation (Ganchev et al., 2010; Mann & McCallum, 2010). In contrast, our methods are designed to discriminatively learn expressive inference procedures, with minimal domain knowledge required, rather than regularizing inference and learning.

First, we provide efficient algorithms for solving the marginal inference problem (2) and performing MAP prediction in the associated distribution, for both convex and non-convex global energy functions. After that, we provide a learning algorithm for and the parametrized functions using an interpretation of (2) as approximate variational inference in a a probabilistic model. All of our algorithms are easy to implement and rely on simple wrappers around black-box inference subroutines.

Our experiments demonstrate the power and generality of our approach by achieving state-of-the-art results on several tasks. We extract accurate citations from research papers by learning discriminative global regularities of valid outputs, outperforming a strong dual decomposition-based baseline (Anzaroot et al., 2014). In a benchmark OCR task (Taskar et al., 2004), we achieve state-of-the-art results with a learned non-convex, non-local energy function, that guides output decodings to lie near dictionary words. Finally, our general algorithm for solving (2) provides large speed improvements for the challenging task of inference in chain-structured collective graphical models (CGMs), applied to bird migration (Sheldon & Dietterich, 2011).

2 Background

Let denote a set of discrete variables and be a collection of input features. We define the conditional distribution , where is a mapping from to a set of sufficient statistics, is a differentiable vector-valued mapping, and . Conditional random fields (CRFs) assume that are given a graph structure and maps to a 0-1 vector capturing joint settings of each clique (Lafferty et al., 2001). Going forward, we often suppress the explicit dependency of on . For fixed , the model is called a Markov random field (MRF).

Given a distribution , define the expected sufficient statistics operator . For the CRF statistics above, is a concatenated vector of node and clique marginals. Therefore, marginal inference, the task of finding the marginal distribution of over , is equivalent to computing the expectation .

For tree-structured graphical models, is a bijection, though this is not true for general graphs. Furthermore, for trees the entropy is equal to the Bethe entropy , defined, for example, in Wainwright & Jordan (2008). The marginal polytope is the set of that correspond to some .

As mentioned in the introduction, marginal inference can be posed as the optimization problem (1). MAP inference finds the joint setting

with maximum probability. For CRFs, this is equivalent to


For tree-structured CRFs, marginal and MAP inference can be performed efficiently using dynamic programming. Our experiments focus on such graphs. However, the inference algorithms we present can be extended to general graphs wherever marginal inference is tractable using a convex entropy approximation and a local polytope relaxation.

3 Marginal Inference With Non-Local Energies

We move beyond the standard inference objective (1), augmenting it with a non-local energy term as in (2):

Here, is some arbitrary parametrized function of the marginals, and may depend on input features .

Intuitively, we are augmenting the inference objective (1) by allowing it to optimize a broader set of tradeoffs – not only between expected node scores, clique scores, and entropy, but also global functions of the marginals. To be concrete, in our citation extraction experiments (Section 8.1), for example, we employ the simple structure:


Where each is a univariate convex function and each is constrained to be non-negative, in order to maintain the overall convexity of . We further employ


where encodes a ‘linear measurement’ of the marginals and is some univariate convex function.

4 Variational Interpretation and Map Prediction

We next provide two complementary interpretations of (2

) as variational inference in a class of tractable probability distributions over

. They yield precisely the same variational expression. However, both are useful because the first helps motivate a MAP prediction algorithm, while the second helps characterize our learning algorithm in Section 7 as (approximate) variational EM.

Proposition 1.

For fixed and , the output of inference in the augmented objective (2) is equivalent to the output of standard inference (1) in an MRF with the same clique structure as our base model, but with a modified parameter .


Forming a Lagrangian for (2), the stationarity conditions with respect to the variable are:


where are collected terms relating to the marginal polytope constraints. The proposition follows because (6) is the same as the stationarity conditions for


Therefore, we can characterize a joint distribution over

by first finding by solving (2) and then defining an MRF over with parameters . Even more conveniently, our inference technique in Section 6

iteratively estimates

on the fly, namely via the dual iterate in Algorithm 1.

Ultimately, in many prediction problems we seek a single output configuration rather than an inferred distribution over outputs. Proposition 1 suggests a simple prediction procedure: first, find the variational distribution over parametrized as an MRF with parameter . Then, perform MAP in this MRF. Assuming an available marginal inference routine for this MRF, we assume the tractability of MAP – for example using a dynamic program. We avoid predicting by locally maximizing nodes’ marginals, since this would not necessarily yield feasible outputs.

Instead of solving (2), we could have introduced global energy terms to the MAP objective (3) that act directly on values rather than on expectations , as in (2

). However, this yields a difficult combinatorial optimization problem for prediction and does not yield a natural way to learn the parametrization of the global energy. Section 

8.1 demonstrates that using energy terms defined on marginals, and performing MAP inference in the associated MRF, performs as well or better than an LP technique designed to directly perform MAP subject to global penalty terms.

Our second variational interpretation characterizes as a variational approximation to a complex joint distribution:


We assume that isolated marginal inference in is tractable, while is an alternative structured distribution over for which we do not have an efficient inference algorithm. Specifically, we assume that can be solved for . Furthemore, we assume that , where is a convex function, conditional on input features . Going forward, we will often surpress the dependence of on . Above, is the normalizing constant of the combined distribution. Note that if was linear, inference in both and would be tractable, since the distribution would decompose over the same cliques as .

Not surprisingly, (8) is intractable to reason about, due to the non-local terms in (2), so we approximate it with a variational distribution . The connection between this variational approximation and Proposition 1 is derived in Appendix A. Here, we assume no clique structure on , but show that minimizing a variational approximation of , for a given , yields a that is parametrized compactly as the MRF in Proposition 1. We discuss the relationship between this and general mean-field inference in Section 5.

Although the analysis of this section assumes convexity of , our inference techniques can be applied to non-convex , as discussed in Section 6.3, and our learning algorithm produces state-of-the-art results even in the non-convex regime for a benchmark OCR task.

5 Related Modeling Techniques

Mean field variational inference in undirected graphical models is a particular application of our inference framework, with a non-convex  (Wainwright & Jordan, 2008). The technique estimates marginal properties of a complex joint distribution using the clique marginals of some tractable base distribution , not necessarily fully factorized. This induces a partitioning of the cliques of into those represented directly by and those where we define clique marginals as a product distribution of the relevant nodes’ marginals in . To account for the energy terms of the full model involving cliques absent in the simple base model, the energy of the base model is augmented with an extra function of .


where is the set of cliques not included in the tractable sub-model, are the potentials of the original graphical model corresponding to the missing cliques, and

represents a repeated outer (tensor) product of the node marginals for the nodes in those cliques.

Note is non-linear and non-convex. Our work generalizes (9) by allowing arbitrary non-linear interaction terms between components of . This is very powerful – for example, in our citation extraction experiments in Section 8.1, expressing these global terms in a standard graphical model would require many factors touching all variables. Local coordinate ascent mean-field can be frustrated by these rigid global terms. Our gradient-based method avoids these issues by updating all marginals simultaneously.

Dual decomposition is a popular method for performing MAP inference in complex structured prediction models by leveraging repeated calls to MAP in tractable submodels (Komodakis et al., 2007; Sontag et al., 2011). The family of models solvable with dual decomposition is limited, however, because the terms that link the submodels must be expressible as linear constraints. Similar MAP techniques (Ravikumar et al., 2010; Martins et al., 2011; Fu & Banerjee, 2013) based on the alternating direction method of multipliers (ADMM) can be adapted for marginal inference, in problems where marginal inference in submodels is tractable. However, the non-local terms are defined as linear functions on settings of graphical model nodes, while our non-linear terms provide practitioners with an expressive means to learn and enforce regularities of the inference output.

Posterior regularization (PR) (Ganchev et al., 2010), learning from measurements (LFM) Liang et al. (2009) , and generalized expectations (GE) (Mann & McCallum, 2010), are a family of closely-related techniques for performing unsupervised or semi-supervised learning of a conditional distribution or a generative model

using expectation-maximization (EM), where the E-step for latent variables

does not come directly from inference in the model, but instead from projection onto a set of expectations obeying global regularity properties. In PR and GE, this yields a projection objective of the form (2), where the terms come from a Lagrangian relaxation of regularity constraints, and corresponds to dual variables. Originally, PR employed linear constraints on marginals, but He et al. (2013) extend the framework to arbitrary convex differentiable functions. Similarly, in LFM such an inference problem arises because we perform posterior inference assuming that the observations have been corrupted under some noise model. Tarlow & Zemel (2012) also present a method for learning with certain forms of non-local losses in a max-margin framework.

Our goals are very different than the above learning methods. We do not impose non-local terms in order to regularize our learning process or allow it to cope with minimal annotation. Instead, we use to increase the expressivity of our model, performing inference for every test example, using a different , since it depends on input features. Since we are effectively ‘learning the regularizer,’ on fully-labeled data, our learning approach in Section 7 differs from these methods. Finally, unlike these frameworks, we employ non-convex terms in some of our experiments. The algorithmic consequences of non-convexity are discussed in Section 6.3.

6 Optimizing the Non-Local Marginal Inference Objective

We now present an approach to solving (2) using non-Euclidean projected gradient methods, which require access to a procedure for marginal inference in the base distribution (which we term the marginal oracle), as well as access to the gradient of the energy function . We pose these algorithms in the composite minimization framework, which gives us access to a wide variety of algorithms that are discussed in the supplementary material.

6.1 Convex Optimization Background

Before presenting our algorithms, we review several definitions from convex analysis (Rockafellar, 1997).

We call a function -strongly convex with respect to a norm , if for all ,

Proposition 2 (e.g. Beck & Teboulle (2003)).

The negative entropy function is 1-strongly convex with respect to the 1-norm over the interior of the simplex (restricting to ).

Given a smooth and strongly convex function , we can also define an associated generalized (asymmetric) distance measure called the Bregman divergence (Bregman, 1967) generated by ,

For example, the KL divergence is the Bregman divergence associated to the negative entropy function, and the squared Euclidean distance is its own associated divergence.

Composite minimization (Passty, 1979) is a family of techniques for minimizing functions of the form , where we have an oracle that allows us to compute minimizations over in closed form (usually here takes the form of a regularizer). Problems of this form are often solved with an algorithm called proximal gradient, which minimizes over some convex set using:

for some decreasing sequence of learning rates . Note that because of the requirement , proximal gradient generalizes projected gradient descent – since unconstrained minimization might take us out of the feasible region , computing the update requires projecting onto .

But there is no reason to use the squared Euclidean distance when computing our updates and performing the projection. In fact, the squared term can be replaced by any Bregman divergence. This family of algorithms includes the mirror descent and dual averaging algorithms (Beck & Teboulle, 2003; Nesterov, 2009).

  Input: parameters , energy function
  set to prox-center
Algorithm 1 Bethe-RDA

We base our projected inference algorithms on regularized dual averaging (RDA) (Xiao, 2010). The updates are:


where is the average gradient of encountered so far. One benefit of RDA is that it does not require the use of a learning rate parameter () when using a strongly convex regularizer. RDA can be interpreted as doing a projection onto using the Bregman divergence generated by the strongly convex function .

6.2 Our Algorithm

These non-Euclidean proximal algorithms are especially helpful when we are unable to compute a projection in terms of Euclidean distance, but can do so using a different Bregman divergence. We will show that this is exactly the case for our problem of projected inference: the marginal oracle allows us to project in terms of KL divergence.

However, to maintain tractability we avoid using the entropy function on the exponentially-large simplex , and instead optimize over the structured, factorized marginal polytope and its corresponding structured Bethe entropy . For tree-structured models, and have identical values, but different inputs. It remains to show the strong convexity of so we can use it in RDA.

Proposition 3.

For trees with nodes, the negative Bethe entropy function is -strongly convex with respect to the 2-norm over the interior of the marginal polytope .


Consequence of Lemma 1 in Fu & Banerjee (2013). ∎

With these definitions in hand, we present Bethe-RDA projected inference Algorithm 1. This algorithm corresponds to instantiating (10) with and . Note the simplicity of the algorithm when choosing . It is intuitively appealing that the algorithm amounts to no more than calling our marginal inference oracle with iteratively modified parameters.

Proposition 4.

For convex energy functions and convex , the sequence of primal averages of Algorithm 1 converges to the optimum of the variational objective (2) with suboptimality of at time .


This follows from Theorem 3 of Xiao (2010) along with the strong convexity of . ∎

If we have more structure in the energy functions, specifically a Lipschitz-continuous gradient, we can modify the algorithm to use Nesterov’s acceleration technique and achieve a convergence of . Details can be found in Appendix D. Additionally, in practice these problems need not be solved to optimality and give stable results after a few iterations, as demonstrated in Figure 1.

6.3 Inference With Non-Convex, Non-Local Energies

An analogy can be made here to loopy belief propagation – even in the case of non-convex loss functions (and even non-convex entropy functions with associated inexact marginal oracles), the updates of our inference (and learning) algorithms are well-defined. Importantly, since one of our motivations for developing non-local inference was to generalize mean field inference, and the additional penalty terms are non-convex in that case, we would like our algorithms to work for the non-convex case as well.

Unlike loopy belief propagation, however, since we derive our algorithms in the framework of composition minimization, we have access to a wealth of theoretical guarantees. Based on results from the theory of optimization with first-order surrogate loss functions (Mairal, 2013), in Appendix C we propose a small modification to Algorithm 1 with an asymptotic convergence condition even for non-convex energies. In practice we find that the unmodified Algorithm 1 also works well for these problems, and experimentally in Section 8.2, we see good performance in both inference and learning with non-convex energy functions.

7 Learning Models With Non-Local Energies

  Input: examples and inference oracle
  for distributions with the clique structure of .
  Output: parameters for .
     for all  do
         (Algorithm 1) // using and MARG()
         (Proposition 5) // using
        // note is a CRF with potentials .
     end for
     //M-Step (gradient-based learning of CRF parameters)
         //standard CRF inference
     until converged
  until converged OR iter max_iters
Algorithm 2 Learning with non-local energies
  Input: examples and
  for distributions with the clique structure of .
  Output: parameters for .
     sample randomly
      (Algorithm 1)
  until converged OR iter max_iters
Algorithm 3 Doubly-stochastic learning with given by a sum of scalar functions of linear measurements (5).

We seek to learn the parameters and of the underlying CRF base model and , respectively. Let be training examples. Let be the variational distribution for resulting from applying Proposition 1. Namely, is an MRF with parameters


We employ the notation to highlight the role of : for a given pair, the family of variational distributions over is indexed by possible values of (recall we suppress the explicit dependence of and on ). Finally, define the shorthand .

interacts with the data in a complex manner that prevents us from using standard learning techniques for the exponential family. Namely, we can not easily differentiate a likelihood with respect to , since this requires differentiating the output of a convex optimization procedure, and the extra term in (2) prevents the use of conjugate duality relationships available for the exponential family. We could have used automatic methods to differentiate the iterative inference procedure (Stoyanov et al., 2011; Domke, 2012), but found our learning algorithm works well.

We employ a variational learning algorithm, presented in Algorithm 2, alternately updating the parameters of our tractable CRF-structured variational distributions, and updating the parameters assuming the following surrogate likelihood given by these CRF approximations:


Given and , we update using Algorithm 1. Given , we update and by taking a single step in the direction of the gradient of the surrogate likelihood (12). We avoid taking more than one gradient step, since the gradients for and depend on and an update to and will break the property that . Therefore, we recompute every time we update the parameters.

Overall, it remains to show how to compute gradients of (12). For , we have the standard CRF likelihood gradient (Sutton & McCallum, 2006):


For , we have:


From (11), is also and


Clearly, this depends on the structure of . Consider the parametrization (4). With this, we have:


Therefore, we have . For linear measurements (5), this amounts to


This has a simple interpretation: the gradient with respect to equals the gradient of the scalar loss at the current marginals times the difference in linear measurements between the ground truth labels and the inferred marginals.

Algorithm 2 has an expensive double-loop structure. In practice it is sufficient to employ a ‘doubly-stochastic’ version given in Algorithm 3, where we sample a training example and use this to only perform a single gradient step on and . To demonstrate the simplicity of implementing our learning algorithm, we avoid any abstract derivative notation in Algorithm 3 by specializing it to the case of (17). In our experiments, however, we sometimes do not use linear measurements. Overall, all our experiments use the fast doubly-stochastic approach of Algorithm 3 solely, since it performs well. In general, our learning algorithms are not guaranteed to converge because we approximate the complex interaction between and with alternating updates. In practice, however, terminating after a fixed number of iterations yields models that generalize well.

Finally, recall that the notation suppresses the potential dependence of on . We assume each is a differentiable function of features of . Therefore, in our experiments where depends on , we perform gradient updates for the parametrization of

via further application of the chain rule.

8 Experiments

8.1 Citation Extraction

Model F1
Our Baseline 94.47
Non-local Energies 95.47
Baseline (Anzaroot et al., 2014) 94.41
Soft-DD (Anzaroot et al., 2014) 95.39
Table 1: Comparison of F1 scores on Citation Extraction dataset. We compare MAP inference F1 scores of our non-local energy model and the specialized dual decomposition model of Anzaroot et al. (2014). Both variants learn global regularities that significantly improve performance.
Figure 1: Citation extraction F1 when limiting maximum number of test-time inference iterations. Most of our accuracy gain is captured within the first 5-10 iterations.

We first apply our algorithm to the NLP task of performing text field segmentation on the UMass citation dataset (Anzaroot & McCallum, 2013), which contains strings of citations from research papers, segmented into fields (author, title, etc.). Our modeling approach, closely follows Anzaroot et al. (2014), who extract segmentations using a linear-chain segmentation model, to which they add a large set of ‘soft’ linear global regularity constraints.

Let be a candidate labeling. Imagine, for example, that we constrain predicted segmentations to have no more predicted last names than first names. Then, the numbers of first and last names can be computed by linear measurements and , respectively. A hard constraint on would enforce . This is relaxed in Anzaroot et al. (2014) to a penalty term


that is added to the MAP inference objective, where is a hinge function. For multiple soft constraints, the overall prediction problem is


where are the parameters of the underlying linear-chain model. They use a dual decomposition style algorithm for solving (19), that crucially relies on the specific structure of the hinge terms . They learn the

for hundreds of ‘soft constraints’ using a perceptron-style algorithm.

We consider the same set of measurement vectors , but impose non-local terms that act on marginals rather than specific values . Further, we use smoothed hinge functions, which improve the convergence rate of inference (Rennie, 2005). We find the variational distribution by solving the marginal inference version of (19), an instance of our inference framework with linear measurements (5):


As in  Anzaroot et al. (2014), we first learn chain CRF parameters on the training set. Then, we learn the parameters on the development set, using Algorithm 3

, and tune hyperparameters for development set performance. At both train and test time, we ignore any terms in (

20) for which .

We present our results in Table 1, measuring segment-level F1. We can see that our baseline chain has slightly higher accuracy than the baseline approach of Anzaroot et al. (2014), possibly due to optimization differences. Our augmented model (Non-Local Energies) matches and very slightly beats their soft dual decomposition (Soft-DD) procedure. This is especially impressive because they employ a specialized linear-programming solver and learning algorithm adapted to the task of MAP inference under hinge-loss soft constraints, whereas we simply plug in our general learning and inference algorithms for non-local structured prediction – applicable to any set of energy functions.

Our comparable performance provides experimental evidence for our intuition that preferences about MAP configurations can be expressed (and “relaxed”) as functions of expectations. Anzaroot et al. (2014) solve a penalized MAP problem directly, while our prediction algorithm first finds a distribution satisfying these preferences, and then performs standard MAP inference in that distribution.

Finally, in Figure 1 we present results demonstrating that our algorithm’s high performance can be obtained using only 5-10 calls per test example to inference in the underlying chain model. In Section B, we analyze the empirical convergence behavior of Algorithm 1.

8.2 Handwriting Recognition

N-Grams 2 3 4 5 6
Accuracy 85.02 96.20 97.21 98.27 98.54
Table 2: Character-wise accuracy of Structured Prediction Cascades (Weiss et al., 2012) on OCR dataset.
SPC (Weiss et al., 2012) Accuracy
2-gram 85.02
3-gram 96.20
4-gram 97.21
5-gram 98.27
6-gram 98.54
Table 3:

Character-wise accuracy of our baselines, and models using learned non-local energies on Handwriting Recognition dataset. Note that word classifier baseline is also given in character-wise accuracy for comparison.

Model Accuracy
2-gram (base model) 84.93
(MM) 94.96
(MM) 98.83
55-Class Classifier (MM) 86.06
Table 4: Character-wise accuracy of our baselines, and models using learned non-local energies on Handwriting Recognition dataset. Note that word classifier baseline is also given in character-wise accuracy for comparison.

We next apply our algorithms to the widely-used handwriting recognition dataset of Taskar et al. (2004). We follow the setup of Weiss et al. (2012), splitting the data into 10 equally sized folds, using 9 for training and one to test. We report the cross-validation results across all 10 folds.

The structured prediction cascades of Weiss et al. (2012) achieve high performance on this dataset by using extremely high order cliques of characters (up to 6-grams), for which they consider only a small number of candidate outputs. Their state-of-the-art results are reproduced in Table 2. The excellent performance of these large-clique models is consequence of the fact that the data contains only 55 unique words, written by 150 different people. Once the model has access to enough higher-order context, the problem becomes much easier to solve.

With this in mind, we design two non-convex, non-local energy functions. These energies are intended to regularize our predictions to lie close to known elements of the vocabulary. Our base model is a standard linear-chain CRF with image features on the nodes, and no features on the bigram edge potentials. Let be a function that takes the concatenated vector of node and edge marginals and sums up all of the node marginals, giving the global unigram expected sufficient statistics. Let indicate the set of all such unique vectors when applying to the train set empirical sufficient statistics for each data case . Simply, this gives 55 vectors of length 26 containing the unigram counts for each unique word in the train set.

Our intuition is that we would like to be able to “nudge” the results of inference in our chain model by pulling the inferred to be close to one of these global statistics vectors. We add the following non-convex non-local energy function to the model:


We learn two variants of this model, which differently parametrize the dependence of on . The first has a single bias feature on the non-local energy. The second conditions on a global representation of the sequence: concretely, we approximate the RBF kernel mean map (MM) (Smola et al., 2007) using random Fourier features (RFF) (Rahimi & Recht, 2007)

. This simply involves multiplying each image feature vector in the sequence by a random matrix with

rows, applying a pointwise non-linearity, and taking to be a linear function of the average vector.

Results of these experiments can be seen in Table 4. Adding the non-local energy brings our performance well above the baseline bigram chain model, and our training procedure is able to give substantially better performance when depends on the above input features.

The energy , based on unigram sufficient statistics, is not able to capture the relative ordering of letters in the vocabulary words, which the structured prediction cascades models do capture. This motivates us to consider another energy function. Let be the set of unique vectors of concatenated node marginal statistics for the train set. This gives 55 vectors of length , where is the length of the th distinct train word. Next, we define a different energy function to add to our base chain model:


Once again we implement featurized and non-featurized versions of this model. As noted in structured prediction cascades, giving the model access to this level of high-order structure in the data makes the inference problem extremely easy. Our model outperforms the best structured prediction cascades results, and we note again an improvement from using the featurized over the non-featurized .

Of course, since the dataset has only 55 actual labels, and some of those are not valid for different input sequences due to length mismatches, this is arguably a classification problem as much as a structured prediction problem. To address this, we create another baseline, which is a constrained 55-class logistic regression classifier (constrained to only allow choosing output classes with appropriate lengths given the input). We use our same global mean-map features from the

variants of the structured model and report these results in Table 4. We also tune the number of random Fourier features as a hyperparameter to give the classifier as much expressive power as possible. As we can see, the performance is still significantly below the best structured models, indicating that the interplay between local and global structure is important.

8.3 Collective Graphical Models

625 10k 50k
Our Method 0.19 2.7 14
IP 2.8 93 690
Table 5: Comparison of runtime (in seconds, averaged over 10 trials) between the interior point solver (IP) of Sheldon et al. (2013) v.s. Algorithm 1 on different CGM problem sizes , the cardinality of the edge potentials in the underlying graphical model, where marginal inference is .

Next, we demonstrate that that our proximal gradient-based inference framework dramatically speeds up approximate inference in collective graphical models (CGMs) (Sheldon & Dietterich, 2011)

. CGMs are a method for structured learning and inference with noisy aggregate observation data. The large-scale dependency structure is represented via a graphical model, but the nodes represent not just single variables, but aggregate sufficient statistics of large sets of underlying variables, corrupted by some noise model. In previous work, CGMs have been successfully applied to modeling bird migration. Here, the base model is a linear chain representing a time series of bird locations. Each observed variable corresponds to counts from bird watchers in different locations. These observations are assumed to be Poisson distributed with rate proportional to the true count of birds present. The CGM MAP task is to infer the underlying migration patterns.

Sheldon et al. (2013) demonstrate that MAP in CGMs is NP-hard, even for trees, but that approximate MAP can be performed by solving a problem of the form (2):


where are (concave) Poisson log-likelihoods and each is an observed bird count.

For the case where the underlying CGM graph is a tree, the ‘hard EM’ learning algorithm of Sheldon et al. (2013) is the same as Algorithm 2 specialized to their model. Therefore, Sheldon et al. (2013) provide additional experimental evidence that our alternating surrogate-likelihood optimization works well in practice.

The learning procedure of Sheldon et al. (2013) is very computationally expensive because they solve instances of (23) using an interior-point solver in the inner loop. For the special case of trees, Algorithm  1 is directly applicable to (23). Using synthetic data and code obtained from the authors, we compare their generic solver to Algorithm  1 for solving instances of (23). In Table 5, we see that our method achieves a large speed-up with no loss in solution accuracy (since it solves the same convex problem).

9 Discussion and Future Work

Our results show that our inference and learning framework allows for tractable modeling of non-local dependency structures, resistant to traditional probabilistic formulations. By approaching structured modeling not via independence assumptions, but as arbitrary penalty functions on the marginal vectors , we open many new modeling possibilities. Additionally, our generic gradient-based inference method can achieve substantial speedups on pre-existing problems of interest. In future work, we will apply our framework to new problems and new domains.


This work was supported in part by the Center for Intelligent Information Retrieval, in part by DARPA under agreement number FA8750-13-2-0020, and in part by NSF grant #CNS-0958392. The U.S. Government is authorized to reproduce and distribute reprints for Governmental purposes notwithstanding any copyright notation thereon. Any opinions, findings and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect those of the sponsor.


  • Anzaroot & McCallum (2013) Anzaroot, Sam and McCallum, Andrew. A new dataset for fine-grained citation field extraction. In ICML Workshop on Peer Reviewing and Publishing Models, 2013.
  • Anzaroot et al. (2014) Anzaroot, Sam, Passos, Alexandre, Belanger, David, and McCallum, Andrew. Learning soft linear constraints with application to citation field extraction. In ACL, 2014.
  • Beck & Teboulle (2003) Beck, Amir and Teboulle, Marc. Mirror descent and nonlinear projected subgradient methods for convex optimization. Operations Research Letters, 31(3):167–175, 2003.
  • Bregman (1967) Bregman, Lev M. 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):200–217, 1967.
  • Domke (2012) Domke, Justin. Generic methods for optimization-based modeling. In AISTATS, 2012.
  • Duchi et al. (2010) Duchi, John, Shalev-Shwartz, Shai, Singer, Yoram, and Tewari, Ambuj. Composite objective mirror descent. In COLT, 2010.
  • Fu & Banerjee (2013) Fu, Qiang and Banerjee, Huahua Wang Arindam. Bethe-admm for tree decomposition based parallel map inference. In UAI, 2013.
  • Ganchev et al. (2010) Ganchev, Kuzman, Graça, Joao, Gillenwater, Jennifer, and Taskar, Ben. Posterior regularization for structured latent variable models. JMLR, 99:2001–2049, 2010.
  • He et al. (2013) He, L., Gillenwater, J., and Taskar, B. Graph-Based Posterior Regularization for Semi-Supervised Structured Prediction. In CoNLL, 2013.
  • Komodakis et al. (2007) Komodakis, Nikos, Paragios, Nikos, and Tziritas, Georgios. Mrf optimization via dual decomposition: Message-passing revisited. In IEEE ICCV, 2007.
  • Lafferty et al. (2001) Lafferty, John, McCallum, Andrew, and Pereira, Fernando CN. Conditional random fields: Probabilistic models for segmenting and labeling sequence data. In ICML, 2001.
  • Liang et al. (2009) Liang, Percy, Jordan, Michael I, and Klein, Dan. Learning from measurements in exponential families. In ICML, 2009.
  • Mairal (2013) Mairal, Julien. Optimization with first-order surrogate functions. In ICML, 2013.
  • Mann & McCallum (2010) Mann, Gideon S and McCallum, Andrew. Generalized expectation criteria for semi-supervised learning with weakly labeled data. JMLR, 11:955–984, 2010.
  • Martins et al. (2011) Martins, André, Figueiredo, Mário, Aguiar, Pedro, Smith, Noah A, and Xing, Eric P. An augmented lagrangian approach to constrained map inference. In ICML, 2011.
  • Nesterov (2009) Nesterov, Yurii. Primal-dual subgradient methods for convex problems. Mathematical programming, 120(1):221–259, 2009.
  • Passty (1979) Passty, Gregory B. Ergodic convergence to a zero of the sum of monotone operators in hilbert space. Journal of Mathematical Analysis and Applications, 72(2):383 – 390, 1979.
  • Rahimi & Recht (2007) Rahimi, Ali and Recht, Benjamin. Random features for large-scale kernel machines. In NIPS, 2007.
  • Ravikumar et al. (2010) Ravikumar, Pradeep, Agarwal, Alekh, and Wainwright, Martin J. Message-passing for graph-structured linear programs: Proximal methods and rounding schemes. JMLR, 11:1043–1080, 2010.
  • Rennie (2005) Rennie, Jason DM. Smooth hinge classification, 2005.
  • Rockafellar (1997) Rockafellar, R Tyrell. Convex Analysis, volume 28. Princeton University Press, 1997.
  • Sheldon et al. (2013) Sheldon, Daniel, Sun, Tao, Kumar, Akshat, and Dietterich, Thomas G. Approximate inference in collective graphical models. In ICML, 2013.
  • Sheldon & Dietterich (2011) Sheldon, Daniel R and Dietterich, Thomas G. Collective graphical models. In NIPS, 2011.
  • Smola et al. (2007) Smola, Alex, Gretton, Arthur, Song, Le, and Schölkopf, Bernhard. A hilbert space embedding for distributions. In Algorithmic Learning Theory, pp. 13–31. Springer, 2007.
  • Sontag et al. (2011) Sontag, David, Globerson, Amir, and Jaakkola, Tommi. Introduction to dual decomposition for inference.

    Optimization for Machine Learning

    , 1:219–254, 2011.
  • Stoyanov et al. (2011) Stoyanov, Veselin, Ropson, Alexander, and Eisner, Jason. Empirical risk minimization of graphical model parameters given approximate inference, decoding, and model structure. In AISTATS, 2011.
  • Sutton & McCallum (2006) Sutton, Charles and McCallum, Andrew. An introduction to conditional random fields for relational learning. Introduction to statistical relational learning, pp. 93–128, 2006.
  • Tarlow & Zemel (2012) Tarlow, Daniel and Zemel, Richard S. Structured output learning with high order loss functions. In AISTATS, 2012.
  • Taskar et al. (2004) Taskar, Ben, Carlos, Guestrin, and Koller, Daphne. Max-margin markov networks. In NIPS, 2004.
  • Wainwright & Jordan (2008) Wainwright, Martin J and Jordan, Michael I. Graphical models, exponential families, and variational inference. Foundations and Trends in Machine Learning, (1-2):1–305, 2008.
  • Weiss et al. (2012) Weiss, D., Sapp, B., and Taskar, B. Structured Prediction Cascades. ArXiv e-prints, August 2012.
  • Xiao (2010) Xiao, Lin. Dual averaging methods for regularized stochastic learning and online optimization. JMLR, 11:2543–2596, 2010.

Appendix A Variational Approximation

During learning, reasoning about in (8) is difficult, due to the intractability of . In response, we approximate it with a variational distribution:




Given , , and , we select by minimizing the approximation (26). Note that the surrogate we minimize is a lower bound to (25), as , by Jensen’s inequality and the convexity of . This differs from many mean-field variational inference approaches that minimize an upper bound.

So far, we have not assumed any structure on . Next, we show that the minimizer of (26) is a MRF with the same clique structure as . This provides an alternative derivation of the techniques in Section 4.

Let denote the probability under of a given joint configuration . There are exponentially-many such , and is the entropy on the simplex . Since minimizes (26), we have the following stationarity condition for every :


Here, is a dual variable for the constraint . Rearranging, we have:


where is a normalizing constant.

Proposition 5.

There exists a vector such that the quantity for all . Furthermore, is a simple, closed-form function of .


We have , since . Therefore, . ∎

Corollary 1.

Since , Proposition 5 implies is an MRF with the same clique decomposition as .

So far, is implicitly defined in terms of its own marginals . Since we assume and have the same sufficient statistics , we can use the Bethe entropy representation . This transforms (26) to the augmented inference problem (2). Therefore, we can directly solve for , which can then be used to provide a closed-form expression for the CRF distribution .

Appendix B Additional Experiments

Figure 2: The number of iterations taken for inference to converge on test set citations, as a percentage of the total number of test cases. Number of iterations is capped at 40. We can see that the distribution is long tailed. Inference converges within 40 iterations for of examples, and each example takes an average of iterations to converge.

In Figure 2, we examine the convergence behavior of our algorithm on the citation dataset. This demonstrates that our inference procedure converges quite quickly except for a small number of difficult cases, where the global energy and the local evidence are in significant disagreement.

Appendix C Non-Convex Energies and Composite Mirror Descent

We introduce a small modification of Algorithm 1, along with a rough proof sketch of its convergence even in the case of non-convex energy functions. Because it leans heavily on significant prior work in optimization, it is hard to give a self-contained proof of the results in this section, and our argument takes the form of a proof sketch that appeals to these other works. However, the basic argument simply combines the strong convexity of and its associated Bregman divergence, along with the results of Mairal (2013) for the case of composite minimization of non-convex functions using the Euclidean Bregman divergence, and the fact that the local updates performed using entropy as a distance-generating function have a log-barrier function for the constraint set , effectively bounding the norm of the gradient of when restricted to the set of iterates actually visited during optimization.

While Algorithm 1 was built on the framework of regularized dual averaging (RDA), we introduce a slightly different formulation based on composite mirror descent (COMID) (Duchi et al., 2010). Like RDA, COMID is a gradient method for minimizing functions of the form . At each time step , COMID makes the update


where is some strongly convex function and is its associated Bregman divergence. In Algorithm 4, we present an instantiation of composite mirror descent for our inference problem.

At first glance, this seems significantly different from our original Algorithm 1, but remembering that because of conjugate duality of the exponential family, we can see that it actually only a corresponds to a slight re-weighting of the iterates of Algorithm 1.

  Input: parameters , energy function , learning rate sequence
  set to prox-center
Algorithm 4 Bethe-MD

First, we give Algorithm 4 similar guarantees in the convex setting as we did for Algorithm 1.

Proposition 6.

For convex energy functions and convex , given the learning rate sequence , where is the strong convexity of , the sequence of primal averages of Algorithm 4 converges to the optimum of the variational objective (2) with suboptimality of at time .


This follows from a standard online-to-batch conversion, along with the strong convexity of and Theorem 7 of Duchi et al. (2010). ∎

Now, having introduced composite mirror descent in  (30), will lean heavily on the framework for optimization with first-order surrogate losses of Mairal (2013) to show that these types of algorithms should converge even in the non-convex case. We now recall a few definitions from that work.

First, we define the asymptotic stationary point condition, which gives us a notion of convergence in the non-convex optimization case.

Definition 1 (Asymptotic Stationary Point (Mairal, 2013)).

For a sequence , and differentiable function , we say it satisfies an asymptotic stationary point condition if

We call a function -strongly smooth if

is a bound on the largest eigenvalue of the Hessian – this tells us how the norm of the gradient changes. This is also known as a

-Lipschitz continuous gradient. Now we recall the notion of a majorant first-order surrogate function.

Definition 2 (Majorant First-Order Surrogate (Mairal, 2013)).

A function g: is a majorant first-order surrogate of near when the following conditions are satisfied

  • Majorant: we have .

  • Smoothness: the approximation error is differentiable, and its gradient is -Lipschitz continuous, moreover, we have and

We denote by the set of such surrogates.

Now we recall the majorant first-order surrogate property for the composite minimization step in the case of Euclidean Bregman divergence (Euclidean distance).

Proposition 7 (Proximal Gradient Surrogates (Mairal, 2013)).

Assume that where is differentiable with an -Lipschitz gradient. Then, admits the following majorant surrogate in :


We can use this result to establish a majorant property for the composite mirror descent surrogate (30) given a strongly convex and strongly smooth Bregman divergence.

Proposition 8 (Composite Mirror Descent Surrogates).

Assume that where is differentiable with an -Lipschitz gradient, is a -strongly convex and -strongly smooth function, and is its Bregman divergence. Then, admits the following majorant surrogate in :


By the definition of strong convexity and the Bregman divergence, (32) upper bounds (31), so it is a majorant of . Additionally, by the additive property of strong smoothness, we get the strong smoothness constant for the surrogate. ∎

However, small technical conditions keep Proposition 8 from applying directly to our case. The Bethe entropy , and thus its associated Bregman divergence, is not strongly smooth – its gradient norm is unbounded as we approach the corners of the marginal polytope. However, it is locally Lipschitz – every point in the domain has a neighborhood for which the function is Lipschitz. In practice, since the mirror descent updates have a barrier function for the constraint set , our iterative algorithm will never get too close to the boundary of the polytope and it is effectively strongly smooth for purposes of our minimization algorithm. This is not a rigorous argument, but is both intuitively plausible and born out in experiments.

Proposition 9.

The sequence of iterates from Algorithm 4, when bounded away from the corners of the marginal polytope constraint set , and for appropriate choice of learning rates , convex , and -strongly smooth (but possibly non-convex) energy function , satisfies an asymptotic stationary point condition.


This follows from application of Proposition 8, and noting that Algorithm 4 corresponds to the generalized surrogate-minimization scheme in Algorithm 1 of Mairal (2013). The asymptotic stationary point condition then follows from Proposition 2.1 of Mairal (2013). The appropriate learning rates must be chosen by the Lipschitz constant of the gradient of , as well as the effective Lipschitz constant of the gradient of , given how far we are bounded from the edge of the constraint set (this effective smoothness constant is determined by the norm of our parameter vector ). ∎

In this section we have given a rough proof sketch for the asymptotic convergence of our inference algorithms even in the case of non-convex energies. Our heuristic argument for the effective smoothness of the entropy

is the most pressing avenue for future work, but we believe it could be made rigorous by examining the norm of the parameter vector and how it contributes to the “sharpness” of the barrier function for the mirror descent iterates.

Appendix D Accelerated Bethe-RDA

  Input: parameters , energy function
  set to prox-center
Algorithm 5 Accelerated Bethe-RDA

If we have -strongly smooth losses ( is a bound on the largest eigenvalue of the Hessian), we can use an accelerated dual averaging procedure to obtain an even faster convergence rate of . Let be the diameter of the marginal polytope as measured by the strongly convex distance-generating function (using its associated Bregman divergence.) Then Algorithm 5 gives us a convergence rate of by Corollary 7 of Xiao (2010).