# On the Convexity and Reliability of
the Bethe Free Energy Approximation

Harald Leisenberger<EMAIL_ADDRESS>
Signal Processing and Speech Communication Laboratory
Graz University of Technology
Graz, Austria Christian Knoll<EMAIL_ADDRESS>
Signal Processing and Speech Communication Laboratory
Graz University of Technology
Graz, Austria Franz Pernkopf<EMAIL_ADDRESS>
Signal Processing and Speech Communication Laboratory
Graz University of Technology
Graz, Austria

###### Abstract

The Bethe free energy approximation provides an effective way for relaxing NP-
hard problems of probabilistic inference. However, its accuracy depends on the
model parameters and particularly degrades if a phase transition in the model
occurs. In this work, we analyze when the Bethe approximation is reliable and
how this can be verified. We argue and show by experiment that it is mostly
accurate if it is convex on a submanifold of its domain, the ’Bethe box’. For
verifying its convexity, we derive two sufficient conditions that are based on
the definiteness properties of the Bethe Hessian matrix: the first uses the
concept of diagonal dominance, and the second decomposes the Bethe Hessian
matrix into a sum of sparse matrices and characterizes the definiteness
properties of the individual matrices in that sum. These theoretical results
provide a simple way to estimate the critical phase transition temperature of
a model. As a practical contribution we propose BETHE-MIN, a projected quasi-
Newton method to efficiently find a minimum of the Bethe free energy.

Keywords: Variational Inference, Bethe Free Energy Approximation,
Probabilistic Graphical Models, Quantum Mechanics, Mathematical Foundations of
Machine Learning.

##  1 Introduction

Statistical models from quantum mechanics have experienced an increasing
popularity in real-world applications [MacKay, 2003, Zdeborová and Krzakala,
2016]. Many of the practical challenges can be modeled by a graph and require
the solution of fundamental probabilistic problems, such as the computation of
the partition function or marginal probabilities. However, an exact solution
to these problems is only available for restricted model classes; e.g., the
one- and two-dimensional Ising model [Ising, 1925, Onsager, 1944] or small
systems [Lauritzen and Spiegelhalter, 1988]. Generally, one needs to apply
approximation methods whose specific choice is often a compromise between
efficiency and accuracy.

A large class of approximation algorithms relies on variational inference and
the Gibbs variational principle [Mezard and Montanari, 2009, Wainwright et
al., 2008]. First, probabilistic inference is formulated in terms of an
optimization problem with respect to the so-called variational Gibbs free
energy; then, to relax the problem complexity, one optimizes a simpler
objective that approximates the solutions to the inference problem. This
framework has its origin in quantum mechanics and mean-field theory [Parisi,
1988], and has been successfully adopted to computer science [Jordan et al.,
1999, Yedidia et al., 2005].

One important aspect is the choice of the auxilary objective to be optimized.
On the one hand, it should be simple enough to entail a true relaxation of the
NP-hard inference problem. On the other hand, a too simple objective might not
capture the global behavior of the model well enough. Over time, different
types of free energy approximations have been proposed, with the Bethe (or
Bethe-Peierls) approximation being one popular representative that often
proves to be superior to other constructions in terms of a tradeoff between
efficiency and accuracy [Bethe, 1935, Peierls, 1936]. It only considers
effects between interacting variables, while ignoring long range correlations
at the cost of an information loss; still, it often yields accurate results
111Often, the quality of the Bethe approximation is analyzed in terms of the
loopy belief propagation (LBP) algorithm; a so-called message passing
algorithm whose relationship to statistical mechanics and the Bethe free
energy has been established by Yedidia et al. [2001]. [Frey and MacKay, 1997,
Murphy et al., 1999].

The quality of the Bethe approximation strongly varies between different
models and particularly suffers from multiple and strong interactions in the
graph. From a physisist’s perspective this is due to a low model temperature;
we may thus interpret a sudden decrease in its accuracy at some critical
temperature as a phase transition in the model. This unwanted behavior of the
Bethe approximation is well known and has often been observed in the
literature [Meshi et al., 2009, Weller et al., 2014, Leisenberger et al.,
2022]; however, a precise quantification of its reliability (e.g., in terms of
some critical temperature) is an open problem 222Sometimes the success of the
Bethe approximation is attributed to the uniqueness of its minimum [Mooij and
Kappen, 2004]. However, we will show that this can be an inaccurate
assumption.. Another prohibitive factor is the practical minimization of the
Bethe free energy which – despite being a relaxation to the Gibbs free energy
– can be difficult.

In this work, we analyze the reliability of the Bethe approximation. We show
when can expect its provided estimates of the partition function and marginals
to be accurate. For that purpose, we distinguish between different ’stages’ of
the Bethe free energy: convexity, non-convexity but uniqueness of a minimum,
and multiple minima. Based on this differentiation, we assess the quality of
the Bethe approximation with respect to its stage.

The intuition for this split lies in the Gibbs variational principle: While
the Gibbs free energy is convex on its domain, this is not guaranteed for the
Bethe free energy. We aim to show that the Bethe approximation is accurate,
whenever it is convex as well. Moreover, we examine if it is still accurate if
the model temperature decreases, i.e., if the Bethe free energy changes its
stage – from convex to non-convex (but uniqueness of a minimum), and later to
multiple minima. To perform this analysis, we must determine the stage of the
Bethe free energy in a certain state of the model. While there exist methods
to prove the uniqueness of a minimum (e.g., Mooij and Kappen [2007]), this is
not the case for its convexity. We address this open issue by deriving two
sufficient conditions for convexity of the Bethe free energy on an
appropriately defined submanifold of its domain, the ’Bethe box’. Both results
rely on a profound second-order analysis of the Bethe free energy, by
characterizing the definiteness properties of its Hessian. They also provide
estimates for a phase transition (i.e., the temperature at which the accuracy
of the Bethe approximation begins to degrade).

Our first result includes conditions under which the Bethe Hessian is
diagonally dominant and hence positive definite (Sec. 3.2). The starting point
of this analysis is to reformulate diagonal dominance as a set of optimization
problems, each defined on a local neighborhood of the graph. By analyzing
these problems, we reduce the question of diagonal dominance of the Bethe
Hessian to the positivity of a certain set of one-dimensional polynomials.

Our second result decomposes the Bethe Hessian into a sum of simpler matrices
and characterizes the definiteness properties of the individual matrices in
that sum (Sec. 3.3). We consider a natural decomposition into ’edge-specific
Hessians’ that represent the curvature of the Bethe free energy on individual
edges and their end nodes. This specific choice induces a closed-form
characterization of positive definiteness that is based on their determinant.

We conclude the theoretical part of this work by comparing our two conditions
to similar results from the literature (Sec. 3.4). Specifically, we consider a
complete graph with a perfectly symmetric potential configuration; a well
understood model whose simplicity facilitates such a direct comparison. On
this model, our second condition is exact in the sense that it correctly
predicts the critical temperature at which a phase transition occurs.

After having established the theoretical framework, we perform an experimental
analysis in which we evaluate the performance of the Bethe approximation on
various graphical models. Specifically, we consider spin glasses (also known
as general or mixed models) on graphs that either have a lattice or a random
structure. For minimizing the Bethe free energy, we propose a ’projected
quasi-Newton’ algorithm (named BETHE-MIN, Sec. 4).

Our experiments show that the Bethe approximation is mostly accurate if it is
convex. If the temperature decreases, the Bethe free energy changes its state
while it also becomes much less accurate at some point. We demonstrate that
the estimated temperature at which the Bethe free energy becomes non-convex
often provides a good estimate of this phase transition. These experiments
help to assess, whether the Bethe approximation is reliable or not.

This paper is structured as follows: Sec. 2 provides all relevant notations
and background. Sec. 3 includes the theoretical part of this work in which
introduce all details on the Bethe free energy and discuss its characteristics
that are known from the literature. Afterwards, we derive our two sufficient
conditions for convexity of the Bethe free energy, before we discuss their
properties on the symmetric model. Sec. 4 includes the practical part of this
work. We first explain our algorithm BETHE-MIN for minimizing the Bethe free
energy, which we then use in our experiments to analyze the quality of the
Bethe approximation with respect to its stage. In Sec. 5 we conclude our work.

##  2 Background

In this section, we introduce the physical and mathematical concepts and
notation that are relevant for this work. This includes matrix algebra and
convex functions (Sec. 2.1), models on graphs and the Boltzmann distribution
(Sec. 2.2), and the Gibbs and Bethe variational principle for probabilistic
inference (Sec. 2.3). We also give an overview of the literature related to
our specific research questions (Sec. 2.4).

###  2.1 Mathematical Preliminaries

A real matrix $\textbf{M}=(m_{ij})_{i,j=1,\dots,n}$ of size $n\times n$ is
positive definite if $\bm{y}^{T}\textbf{M}\bm{y}>0$ for all vectors $\bm{y}$
in $\mathbb{R}^{n}$ whenever $\bm{y}\neq\bm{0}$, and positive semidefinite if
$\bm{y}^{T}\textbf{M}\bm{y}\geq 0$ for all $\bm{y}$ in $\mathbb{R}^{n}$.
Equivalently, M is positive definite if all its leading principal minors 333A
leading principal minor is the determinant of a $(k\times k)$ \- submatrix of
M, consisting of the first $k$ rows and columns of M (for $k$ in
$\\{1,\dots,n\\}$). are positive (Sylvester’s Criterion). The sum of positive
definite matrices is positive definite. Furthermore, M is diagonally dominant
if, for each row $i$, the magnitude of the corresponding diagonal entry is
larger than the sum of the magnitudes of all other entries in that row, i.e.,

$\displaystyle|m_{ii}|>\sum\limits_{j\neq i}|m_{ij}|.$ (1)

A diagonally dominant matrix with positive diagonal entries is positive
definite (e.g., [Golub and Van Loan, 2013]).

Let $f$ be a real-valued and twice continuously differentiable function
defined on an open, convex set $\mathcal{D}\subseteq\mathbb{R}^{n}$. Let
$\nabla f(\bm{p})=\big{(}\frac{\partial f}{\partial
p_{1}}(\bm{p}),\dots,\frac{\partial f}{\partial p_{n}}(\bm{p})\big{)}^{T}$ be
the gradient of $f$ and let
$\textbf{H}_{f}(\bm{p})\coloneqq\nabla^{2}f(\bm{p})=\big{(}\frac{\partial^{2}f}{\partial
p_{i}p_{j}}(\bm{p})\big{)}_{i,j=1,\dots,n}$ be the Hessian matrix of $f$, each
evaluated in some point $\bm{p}$ of $\mathcal{D}$. If $\nabla
f(\bm{p})=\bm{0}$, then $\bm{p}$ is a stationary point of $f$. A stationary
point $\bm{p}$ is a local minimum (resp., maximum) of $f$ if
$\textbf{H}_{f}(\bm{p})$ (resp., $-\textbf{H}_{f}(\bm{p})$ ) is positive
definite, and it is a saddle point if neither $\textbf{H}_{f}(\bm{p})$ nor
$-\textbf{H}_{f}(\bm{p})$ are positive definite. We define $f$ to be strictly
convex on $\mathcal{D}$ if $\textbf{H}_{f}$ is positive definite 444This
criterion is often easier to work with than the more common yet equivalent
definition of strict convexity: $f$ is strictly convex on $\mathcal{D}$ if,
for any pair of distinct points $\bm{p_{1}},\bm{p_{2}}$ in $\mathcal{D}$, the
straight line that connects $f(\bm{p_{1}})$ and $f(\bm{p_{2}})$ lies above the
curve of $f$ between $\bm{p_{1}}$ and $\bm{p_{2}}$, i.e., $\gamma
f(\bm{p_{1}})+(1-\gamma)f(\bm{p_{2}})>f(\gamma
f(\bm{p_{1}})+(1-\gamma)f(\bm{p_{2}}))$ for $0<\gamma<1$. for all $\bm{p}$ in
$\mathcal{D}$. In that case, $f$ has at most one local (and thus global)
minimum in $\mathcal{D}$ (e.g., Rockafellar [1970]).

###  2.2 Model Specification

We consider a finite system of $N$ variables $X_{1},\dots,X_{N}$ with an
identical configuration space $\mathcal{X}_{i}=\mathcal{X}$ describing the
possible states of individual variables. A configuration
$\bm{x}=(x_{1},\dots,x_{N})$ in the product space
$\mathcal{X}^{N}=\mathcal{X}\times\dots\times\mathcal{X}$ thus represents one
possible state of the entire system. A real-valued energy function $E(\bm{x})$
characterizes the behavior of the system by assigning each system state to its
energy. Then the Boltzmann distribution is defined as

$\displaystyle p_{\beta}(\bm{x})=\frac{1}{Z(\beta)}e^{-\beta E(\bm{x})},$ (2)

where $\beta=\frac{1}{T}$ is the inverse temperature 555Usually the
temperature is scaled by $k_{B}$, Boltzmann’s constant. Without loss of
generality we set $k_{B}=1$, hence measure temperature in units of energy. of
the system and $Z(\beta)$ is the (temperature dependent) normalization
constant, known as partition function. Hence, $p_{\beta}(\bm{x})$ gives the
probability that the system is found in state $\bm{x}$. We focus on an
important class of models, satisfying the following conditions:

  * •

Variables are described in terms of their spin which takes one of two states
(an _Ising spin_). The configuration space is thus $\mathcal{X}=\\{+1,-1\\}$.

  * •

Variables interact in pairs only. Specifically, we distinguish between
_ferromagnetic_ (the spins are aligned) and _antiferromagnetic_ (the spins are
anti-aligned) interactions.

  * •

An external (not necessarily uniform) magnetic field $\theta_{i}$ additionally
influences the spin of variables.

  * •

A graph $\mathbf{G}=(\mathbf{X},\mathbf{E})$ models the structure of the
system. Its node set 666With a slight abuse of notation, we usually identify a
variable $X_{i}$ with its representing node $i$ in the graph instead of making
an explicit distinction between them. $\mathbf{X}=\\{X_{1},\dots,X_{N}\\}$
represents the variables and its edge set $\mathbf{E}$ represents all
interacting pairs of variables. Depending on the kind of interaction, we call
edges ferromagnetic or antiferromagnetic. Furthermore, let $\mathcal{N}(i)$ be
the neighborhood of $X_{i}$ in the graph (i.e., the set of all variables
interacting with $X_{i}$) and $d_{i}=|\mathcal{N}(i)|$ be the degree of
$X_{i}$.

Under these assumptions, the energy function $E(\bm{x})$ decomposes into a sum
of pairwise and single terms:

$\displaystyle
E(\bm{x})=\sum\limits_{(i,j)\in\mathbf{E}}\\!\\!E_{ij}(x_{i},x_{j})+\sum\limits_{i\in\mathbf{X}}E_{i}(x_{i})=-\\!\\!\sum\limits_{(i,j)\in\mathbf{E}}\\!\\!J_{ij}x_{i}x_{j}-\sum\limits_{i\in\mathbf{X}}\theta_{i}x_{i}.$
(3)

The _couplings_ $J_{ij}$ are positive if $X_{i}$ and $X_{j}$ share a
ferromagnetic edge, and negative if they share an antiferromagnetic edge. The
_field_ $\theta_{i}$ biases the spin of variable $X_{i}$ towards $+1$ if
$\theta_{i}>0$, and towards $-1$ if $\theta_{i}<0$. A model with purely
positive (negative) couplings is ferromagnetic (antiferromagnetic); a model
with both positive and negative couplings is a _spin glass_. A spin glass that
is defined on a graph which contains cycles with an odd number of
antiferromagnetic edges (implying that the local constraints on the involved
spin pairs cannot be satisfied simultanously) is called _frustrated_.

Two fundamentally important problems of probabilistic inference 777Note that
these problems are closely connected, since marginals can be written as a
ratio between a ’partial’ and the ’true’ partition function [Weller and
Jebara, 2014b]. are:

  1. (P1)

The computation of the partition function:

$\displaystyle Z(\beta)=\sum\limits_{\bm{x}\in\mathcal{X}^{N}}e^{-\beta
E(\bm{x})}$ (4)

  2. (P2)

The computation of marginal probabilities of the Boltzmann distribution,
primarily of single variables and pairs of variables:

$\displaystyle p_{i}(x_{i})=$
$\displaystyle\sum\limits_{\bm{x}^{\prime}\in\mathcal{X}^{N}:x_{i}^{\prime}=x_{i}}p_{\beta}(\bm{x}^{\prime})$
(5) $\displaystyle p_{ij}(x_{i},x_{j})=$
$\displaystyle\sum\limits_{\bm{x}^{\prime}\in\mathcal{X}^{N}:x_{i}^{\prime}=x_{i},x_{j}^{\prime}=x_{j}}p_{\beta}(\bm{x}^{\prime})$
(6)

These problems can be solved efficiently on tree-structured (i.e., acyclic)
graphs by using dynamical programming [Gallager, 1962, Pearl, 1988], but are
generally NP-hard on _loopy_ (i.e., cyclic) graphs [Valiant, 1979, Cooper,
1990]. In that case one needs to approximate the solution. One large class of
approximation algorithms is based on variational inference. This approach
converts the inference problem into an (also intractable) optimization
problem, and afterwards optimizes an auxiliary objective that should be close
to the original one. The solution to the auxiliary problem is then used as an
estimate of the solution to the inference problem [Jordan et al., 1999,
Wainwright et al., 2008]. In Sec. 2.3 e present further details on the (Gibbs)
variational principle.

###  2.3 Gibbs Variational Principle and Bethe Approximation

We can address the problems (P1) and (P2) stated in Sec. 2.2 by reformulating
them in terms of the Gibbs variational principle [Mezard and Montanari, 2009,
Wainwright et al., 2008]. The idea is to interpret the partition function and
marginals as minima of the variational Gibbs free energy, which is then – due
to its intractability – approximated by an ’auxiliary free energy’ that is
easier to optimize. Examples are the Bethe free energy [Bethe, 1935, Peierls,
1936], the tree-reweighted free energies [Wainwright et al., 2005], and the
Kikuchi free energies [Kikuchi, 1951, Pelizzola, 2005, Yedidia et al., 2005],
each being defined over different clusters of variables. Their minima can be
obtained by numerical optimization and serve as estimates of the solutions to
(P1) and (P2).

Instead of directly minimizing an auxiliary free energy, one can alternatively
use so-called _message passing algorithms_. These are iterative procedures on
the graph whose fixed points are identical to minima of a certain kind of free
energy approximation 888More precisely, each class of message passing
algorithm corresponds to a specific type of variational free energy
approximation and vice versa.. The best-known such relationship is the duality
between the Bethe approximation and the _loopy belief propagation_ (LBP)
algorithm [Yedidia et al., 2001]. We do not primarily deal with LBP in this
work, but draw conclusions about its behavior by analyzing the Bethe free
energy (Sec. 3).

Let us introduce the Gibbs variational principle and Bethe free energy
approximation: Let $q(\bm{x})$ be any probability distribution over the
product space $\mathcal{X}^{N}$. Then we define the _variational Gibbs free
energy_ as the functional

$\displaystyle\mathcal{F}(q)=\mathcal{U}(q)-\frac{1}{\beta}\mathcal{S}(q)$ (7)

with the average energy (or expected energy)

$\displaystyle\mathcal{U}(q)=\mathbb{E}_{q}(E(\bm{x}))=\sum\limits_{\bm{x}\in\mathcal{X}^{N}}q(\bm{x})E(\bm{x})$
(8)

and the entropy

$\displaystyle\mathcal{S}(q)=-\sum\limits_{\bm{x}\in\mathcal{X}^{N}}q(\bm{x})\log
q(\bm{x}).$ (9)

Splitting the Gibbs free energy into the Kullback-Leibler divergence 999The
Kullback-Leibler divergence [Cover and Thomas, 1991] measures the ’distance’
between two probability distributions $p_{1}$ and $p_{2}$ over the same
configuration space $\mathcal{Y}$ and is – in the discrete case – defined as
$D_{KL}(p_{1}||p_{2})=-\sum\limits_{\bm{y}\in\mathcal{Y}}p_{1}(\bm{y})\log\Big{(}\frac{p_{2}(\bm{y})}{p_{1}(\bm{y})}\Big{)}$.
between $q$ and the Boltzmann distribution $p_{\beta}$, and the partition
function $Z(\beta)$ as

$\displaystyle\mathcal{F}(q)=\frac{1}{\beta}D(q||p_{\beta})-\frac{1}{\beta}\log
Z(\beta)$ (10)

shows that $\mathcal{F}(q)$ is convex and has a unique minimum if $q$ is
chosen as $p_{\beta}$ with the functional value
$\mathcal{F}(p_{\beta})=-\frac{1}{\beta}\log Z(\beta)$ (known as Helmholtz
free energy).

While minimizing $\mathcal{F}$ over all distributions $q$ would leave us with
an even increased complexity, the idea is to exploit the variational framework
and replace $\mathcal{F}$ by a simpler objective. One sophisticated approach
is to introduce the _Bethe free energy_ that makes two approximations: First,
it relaxes the space of feasible distributions $q$ to a larger space of
’pseudo-distributions’ $\tilde{p}$ that must only satisfy local instead of
global probability constraints. More precisely, $\tilde{p}$ induces pairwise
pseudo-marginals $\tilde{p}_{ij}$ and singleton pseudo-marginals
$\tilde{p}_{i}$ (all normalized) such that, for all pairs
$(i,j)\in\mathbf{E}$, summing $\tilde{p}_{ij}$ over all states of variable $j$
(or $i$) yields $\tilde{p}_{i}$ (or $\tilde{p}_{j}$); on the other hand,
global constraints of the form
$\sum\limits_{\bm{x}^{\prime}\in\mathcal{X}^{N}:x_{i}^{\prime}=x_{i},x_{j}^{\prime}=x_{j}}\tilde{p}(x_{1}^{\prime},\dots,x_{N}^{\prime})=\tilde{p}_{ij}(x_{i},x_{j})$
and $\sum\limits_{\bm{x}\in\mathcal{X}^{N}}\tilde{p}(\bm{x})=1$ need not be
satisfied. Formally, we describe the set $\mathbb{L}$ of all pseudo-
distributions – sometimes called the _local polytope_ 101010Note that
$\mathbb{L}$ is a convex set that is bounded by
$\mathcal{O}(|\mathbf{X}|+|\mathbf{E}|)$ linear constraints [Wainwright et
al., 2008]. – as

$\displaystyle\begin{split}\mathbb{L}=\\{\,\tilde{p}:\mathcal{X}^{N}\rightarrow\mathbb{R}_{>0}\,\,|&\,\,\,\sum\limits_{x_{j}\in\mathcal{X}}\tilde{p}_{ij}(x_{i},x_{j})=\tilde{p}_{i}(x_{i}),\\\
&\sum\limits_{x_{i},x_{j}\in\mathcal{X}}\tilde{p}_{ij}(x_{i},x_{j})=1,\sum\limits_{x_{i}\in\mathcal{X}}\tilde{p}_{i}(x_{i})=1,\\\
&\,\,\text{for
all}\,\,\,(i,j)\in\mathbf{E}\,\,\,\text{and}\,\,\,i\in\mathbf{X}\\}.\end{split}$
(11)

Second, the Bethe free energy approximates the entropy $\mathcal{S}(q)$ by a
simpler expression, the _Bethe entropy_ , that only takes pairwise and local
statistics into account. For $\tilde{p}\in\mathbb{L}$, we define the
variational Bethe free energy as

$\displaystyle\mathcal{F}_{B}(\tilde{p})=\mathcal{U}_{{B}}(\tilde{p})-\frac{1}{\beta}\mathcal{S}_{{B}}(\tilde{p})$
(12)

with the _Bethe average energy_

$\displaystyle\begin{split}\mathcal{U}_{{B}}(\tilde{p})=&\sum_{(i,j)\in\mathbf{E}}\sum_{x_{i},x_{j}\in\mathcal{X}}\tilde{p}_{ij}(x_{i},x_{j})\,E_{ij}(x_{i},x_{j})\,+\,\sum_{i\in\mathbf{X}}\sum\limits_{x_{i}\in\mathcal{X}}\,\tilde{p}_{i}(x_{i})\,E_{i}(x_{i})\end{split}$
(13)

and the _Bethe entropy_

$\displaystyle\begin{split}\mathcal{S}_{{B}}(\tilde{p})=&\sum_{(i,j)\in\mathbf{E}}\overbrace{\Big{(}-\sum_{x_{i},x_{j}\in\mathcal{X}}\tilde{p}_{ij}(x_{i},x_{j})\log\tilde{p}_{ij}(x_{i},x_{j})\Big{)}}^{\mathcal{S}_{ij}(\tilde{p})}\\\
&-\sum_{i\in\mathbf{X}}\,(d_{i}-1)\underbrace{\Big{(}-\sum_{x_{i}\in\mathcal{X}}\tilde{p}_{i}(x_{i})\log\tilde{p}_{i}(x_{i})\Big{)}}_{\mathcal{S}_{i}(\tilde{p})}.\end{split}$
(14)

The Bethe entropy $\mathcal{S}_{{B}}$ uses the sum over all pairwise entropies
$\mathcal{S}_{ij}$ to approximate $\mathcal{S}(q)$. This simplified
computation results in an overcounting of entropy contributions from
individual nodes, which is compensated by reducing the sum by a
$(d_{i}-1)$-multiple of the local entropies $\mathcal{S}_{i}$. Note that
$\mathcal{F}_{B}$ depends on $\tilde{p}$ only via the singleton and pairwise
pseudo-marginals $\tilde{p}_{i},\tilde{p}_{ij}$ (for each node and edge). This
allows us, with a slight abuse of notation, to identify any element
$\tilde{p}\in\mathbb{L}$ with the collection of all its associated pseudo-
marginals $\\{\tilde{p}_{i},\tilde{p}_{ij}\\}$.

By minimizing $\mathcal{F}_{B}$ over $\mathbb{L}$, one usually tries to
approximate the partition function $Z(\beta)$ and/or the marginals
$p_{i},p_{ij}$ induced by $p_{\beta}$ according to

$\displaystyle\min\limits_{\tilde{p}_{i},\tilde{p}_{ij}\,\in\,\mathbb{L}}\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{F}_{B}\,\,$
$\displaystyle\approx\,\,-\frac{1}{\beta}\log Z(\beta)$ (15)
$\displaystyle\,\,\operatorname*{argmin}\limits_{\tilde{p}_{i},\tilde{p}_{ij}\,\in\,\mathbb{L}}\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\mathcal{F}_{B}\,\,$
$\displaystyle\approx\,\,\\{\,p_{i},p_{ij}\,\,|\,\,(i,j)\in\mathbf{E}\,\,\,\text{and}\,\,\,i\in\mathbf{X}\\}.$
(16)

The expressions on the left-hand side are well-defined, because
$\mathcal{F}_{B}$ is bounded from below and its derivatives tend to infinity
as $\mathcal{F}_{B}$ approaches the boundary of $\mathbb{L}$; i.e.,
$\mathcal{F}_{B}$ attains its global minimum in an interior point of
$\mathbb{L}$ [Yedidia et al., 2005, Weller and Jebara, 2013].

The quality of the approximations made in (15) and (16) degrades with a high
connectivity and strong correlations (or equivalently, a low temperature) in
the graph [Meshi et al., 2009, Weller et al., 2014, Leisenberger et al.,
2022]. This can be explained by the shape of the Bethe free energy: Whenever
the graph contains at least two cycles, $\mathcal{F}_{B}$ is – in contrast to
the Gibbs free energy – non-convex and may lose the favorable property of
having a unique global minimum [Pakzad and Anantharam, 2005, Heskes, 2006,
Watanabe and Fukumizu, 2009]. However, this is not the complete story and we
will argue in Sec. 3 that the Bethe approximation is still accurate if it is
convex on a specific submanifold of $\mathbb{L}$.

###  2.4 Related Work

Statistical mechanics deals with infinite systems and their properties in
thermodynamic equilibrium [Reif, 1965, Huang, 1987]. On the basis of
simplified mathematical models, it explains a system’s macroscopic behavior in
terms of its microscopic constitutents. The statistical properties of the
system are fully described by the partition function and its derivatives;
however, their analytical form is only known for a few models, e.g., the one-
and two-dimensional Ising model, [Ising, 1925, Onsager, 1944]. Of particular
interest is the phenomenon of phase transitions, i.e., sudden changes in the
qualitative behavior of the system [Yeomans, 1992]. These are manifested as
singularities of the thermodynamic potentials (or the partition function) with
respect to the temperature.

The theory of Gibbs measures [Georgii, 1988] provides an elegant framework to
study phase transitions directly in infinite systems, instead of computing the
thermodynamic limit of finite systems. A Gibbs measure represents an
equilibrium state of an infinite system and we may interpret a phase
transition as the non-uniqueness of a Gibbs measure 111111Formally, a
probability distribution $\mu$ on an infinite system $\mathcal{S}$ is a Gibbs
measure if, for all finite subsystems $\Lambda\subset\mathcal{S}$, the
conditional probabilities induced by $\mu$ on $\Lambda$ are consistent with
the canonical ensemble in $\Lambda$, given that all degrees of freedom outside
$\Lambda$ are frozen.. Conversely, the existence of a unique Gibbs measure
predicts an absence of phase transitions in the system. As a consequence,
researchers have made effort to derive conditions for uniqueness of Gibbs
measures [Dobrushin, 1968, Simon, 1979].

In computer science, many applications use probabilistic graphical models
(PGMs) to represent a probability distribution [Koller and Friedman, 2009].
These are principally identical to models from statistical mechanics, with the
major difference that the involved graphs are often relatively small and do
not necesarrily exhibit a lattice structure. In particular, finite systems do
not exhibit phase transitions in the proper sense, because the thermodynamic
potentials are always analytical functions. [Brezin and Zinn-Justin, 1985,
Binder, 1987].

Researchers in this field strongly benefited from statistical mechanics and
found a surprising connection between concepts that seem entirely unrelated:
the loopy belief propagation (LBP) algorithm on the one hand, which is an
iterative method to approximate marginal probabilities in PGMs; and the Bethe-
Peierls approximation on the other hand, which has its origins in the theory
of superlattices [Bethe, 1935, Peierls, 1936]. More precisely, Yedidia et al.
[2001] have proven that there exists a bijective mapping between fixed points
121212LBP iteratively exchanges statistical information in form of ’messages’
between neighboring nodes in the graph. If this process converges, then we say
that LBP has converged to a fixed point. of LBP and stationary points of the
Bethe free energy. This was later complemented by the result that each stable
fixed point of LBP corresponds to a local minimum of the Bethe free energy
[Heskes, 2003]. A successful application of the Bethe approximation (or loopy
belief propagation) is defined in terms of its quality towards estimating
marginal probabilities and the partition function. Due to the lack of general
performance guarantees (Sec. 2.3), one requires alternative ways to assess the
reliability of LBP and the Bethe approximation.

One line of research analyzes when LBP converges to a unique fixed point. This
is often considered as an indicator that the LBP solution (and thus the Bethe
approximation) is accurate. As is known, LBP always converges to a unique
fixed point, if the graph is a tree or contains precisely one cycle [Pearl,
1988, Weiss, 2000]; for other graphs, the number of fixed points and the
convergence behavior depend on the potential configuration (or temperature).
Results on this topic include the following: Tatikonda and Jordan [2002] have
shown that LBP converges to a unique fixed point, if there exists a unique
Gibbs measure on the computation tree 131313The computation tree is an
’unwrapping’ of the graph, starting from some node $i$ and consisting of all
non-backtracking paths. E.g., the Cayley tree is the computation tree of a
complete graph. The concept of Gibbs measures is well defined on the
computation tree, because it consists of infinitely many nodes.; i.e., they
established a direct connection to the theory of Gibbs measures. Heskes [2004]
formulated a minimax problem on the Bethe free energy whose concavity is
sufficient for the uniqueness of an LBP fixed point. Ihler et al. [2005],
Mooij and Kappen [2007] consider LBP as a dynamical system and use contraction
arguments to prove convergence to a unique fixed point. Further, Watanabe and
Fukumizu [2009] use the Ihara zeta function from graph theory to prove
uniqueness of an LBP fixed point. Figure 1 shows the relations between the
different criteria that are described above.

Alternatively, one can try to quantify the approximation error induced by the
Bethe free energy or LBP. Various results bound its approximation error with
respect to the marginal probabilities, either in closed form [Tatikonda, 2003,
Taga and Mase, 2006a] or recursively [Ihler, 2007, Mooij and Kappen, 2008,
Weller and Jebara, 2013]. Others use loop series [Chertkov and Chernyak, 2006]
or graph covers [Vontobel, 2013] to establish a relationship between the Gibbs
and the Bethe minimum, which induces bounds on the approximation error with
respect to the partition function. In ferromagnetic models, the Bethe minimum
bounds the Gibbs minimum from above [Willsky et al., 2007, Ruozzi, 2013] and
can be approximated in polynomial time [Weller and Jebara, 2014a].
Unfortunately, these bounds can be arbitrarily loose and are sometimes
difficult to compute.

Intuitively, as the Bethe free energy is used to approximate a convex
functional (the Gibbs free energy, Sec. 2.3), one would expect it to provide
accurate estimates whenever it is convex as well. However, Watanabe and
Fukumizu [2009] have proven that the Bethe free energy is convex, if and only
if the graph is either a cycle or a tree. This might imply that the concept of
convexity might not be appropriate for inferring about the quality of the
Bethe approximation (because it cannot be convex except for these special
graphs anyway).

In this work, we will have a closer look at the convexity of the Bethe free
energy. We show how it can be made a useful concept for assessing the quality
of the Bethe approximation.

Figure 1: Summary of the current state of knowledge about the relationships
between the Bethe free energy and loopy belief propagation. The white
rectangles contain various properties of the Bethe free energy and LBP. The
main results of this work are contained in the blue ellipse and serve as
sufficient conditions for the convexity of the Bethe free energy. Sufficient
conditions from previous research are included in the green ellipses. Other
important results are shown as red label next to the logical arrows between
properties.

##  3 Convexity of the Bethe Free Energy: Theory

This section includes our main theoretical contributions in this work.
Particularly, we address the question under which conditions $\mathcal{F}_{B}$
is convex and how this can be verified. Our theory is built on the results of
previous work and exploits the current state of knowledge as a starting point.
In detail, we proceed as follows: First, we summarize the most important
properties of the Bethe free energy and specify the problem to be tackled
(Sec. 3.1). Next, we perform a second-order analysis of the Bethe free energy
and derive two sufficient conditions for its convexity, one of which is based
on the criterion of diagonal dominance (Sec. 3.2), and the other on a sum
decomposition of the Hessian matrix (Sec. 3.3). Finally, we compare our
results to similar statements from the literature with respect to a perfectly
symmetric model that facilitates a direct comparison (Sec. 3.4). Our notation
in this section closely follows that of Welling and Teh [2001], Weller and
Jebara [2013], Leisenberger et al. [2022].

###  3.1 Properties of the Bethe Free Energy

We first make a range of simplifications regarding the dimensionality of the
local polytope (11) and the Bethe free energy (12) for binary variables. Let
$X_{i},X_{j}$ be any pair of connected variables and let
$\tilde{p}_{i},\tilde{p}_{j},\tilde{p}_{ij}$ be a set of associated pseudo-
marginal vectors 141414More precisely,
$\tilde{p}_{i}=\begin{pmatrix}\tilde{p}_{i}(X_{i}=+1)\\\
\tilde{p}_{i}(X_{i}=-1)\end{pmatrix}$ and
$\tilde{p}_{ij}=\begin{pmatrix}\tilde{p}_{ij}(X_{i}=+1,X_{j}=+1)\\\
\tilde{p}_{ij}(X_{i}=+1,X_{j}=-1)\\\ \tilde{p}_{ij}(X_{i}=-1,X_{j}=+1)\\\
\tilde{p}_{ij}(X_{i}=-1,X_{j}=-1)\end{pmatrix}$. from the local polytope
$\mathbb{L}$ (i.e., with all entries satisfying the constraints in the
definition of $\mathbb{L}$). Let us denote the individual pseudo-marginal
probabilities for $X_{i}$ and $X_{j}$ being in state $+1$ by
$q_{i}\coloneqq\tilde{p}_{i}(X_{i}=+1)$ and
$q_{j}\coloneqq\tilde{p}_{j}(X_{j}=+1)$ and let us further denote the joint
pseudo-marginal probability for $X_{i},X_{j}$ being both in state $+1$ by
$\xi_{ij}\coloneqq\tilde{p}_{ij}(X_{i}=+1,X_{j}=+1)$. Then, as
$\tilde{p}_{i},\tilde{p}_{j},\tilde{p}_{ij}$ must satisfy the constraints of
$\mathbb{L}$, all five remaining entries of these vectors can already be
expressed in terms of $q_{i},q_{j},\xi_{ij}$ according to the following table:

Table 1: Joint probability table of two binary variables $X_{i}$ and $X_{j}$. $\tilde{p}_{ij}(X_{i},X_{j})$ | $X_{j}=+1$ | $X_{j}=-1$ |
---|---|---|---
$X_{i}=+1$ | $\xi_{ij}$ | $q_{i}-\xi_{ij}$ | $q_{i}$
$X_{i}=-1$ | $q_{j}-\xi_{ij}$ | $1+\xi_{ij}-q_{i}-q_{j}$ | $1-q_{i}$
| $q_{j}$ | $1-q_{j}$ |

Additionally, we have constraints on the independent parameters
$q_{i},q_{j},\xi_{ij}$ of the form

$\displaystyle 0<$ $\displaystyle\,\,q_{i}<1,$ (17) $\displaystyle 0<$
$\displaystyle\,\,q_{j}<1,$ (18) $\displaystyle\max(0,q_{i}+q_{j}-1)<$
$\displaystyle\,\,\xi_{ij}<\min(q_{i},q_{j}).$ (19)

This step reparameterizes the local polytope by exploiting probabilistic
dependencies between variables and removing a considerable number of redundant
parameters. If we collect all node-specific variables $q_{i}$ in a vector
$\bm{q}$ and all edge-specific parameters $\xi_{ij}$ in a vector $\bm{\xi}$,
we can compactly rewrite the local polytope (11) as the set

$\displaystyle\begin{split}\mathbb{L}=\\{&(\bm{q};\bm{\xi})\in\mathbb{R}^{\lvert\mathbf{X}\rvert+\lvert\mathbf{E}\rvert}:0<q_{i}<1,i\in\mathbf{X};\\\
&\max(0,q_{i}+q_{j}-1)<\xi_{ij}<\min(q_{i},q_{j}),(i,j)\in\mathbf{E}\\}.\end{split}$
(20)

Figure 2 shows the local polytope on a graph with two vertices and a single
edge.

(a)  (b)

Figure 2: Local polytope $\mathbb{L}$ (left) on a single-edge graph (right).

This reparameterization of the local polytope also reparameterizes the Bethe
free energy: By substituting the expressions from Table 1 into the
corresponding terms of (12) - (14), $\mathcal{F}_{B}$ takes the simpler form

$\displaystyle\begin{split}\mathcal{F}_{B}(\bm{q};\bm{\xi})=&\overbrace{-\sum_{(i,j)\in\mathbf{E}}\,(1+2\;(2\,\xi_{ij}-q_{i}-q_{j}))\,J_{ij}+\,\,\sum_{i\in\mathbf{X}}(1-2q_{i})\,\theta_{i}}^{\mathcal{U}_{{B}}(\bm{q};\bm{\xi})}\\\
&-\frac{1}{\beta}\underbrace{\Big{(}\sum_{(i,j)\in\mathbf{E}}\mathcal{S}_{ij}\;-\;\sum_{i\in\mathbf{X}}(d_{i}-1)\,\mathcal{S}_{i}\Big{)}}_{\mathcal{S}_{{B}}(\bm{q};\bm{\xi})},\end{split}$
(21)

with the pairwise entropies

$\displaystyle\begin{split}\mathcal{S}_{ij}=&-\xi_{ij}\log\xi_{ij}-(1+\xi_{ij}-q_{i}-q_{j})\log(1+\xi_{ij}-q_{i}-q_{j})\\\
&-(q_{i}-\xi_{ij})\log(q_{i}-\xi_{ij})-(q_{j}-\xi_{ij})\log(q_{j}-\xi_{ij})\end{split}$
(22)

and the local entropies

$\displaystyle\begin{split}\mathcal{S}_{i}=-q_{i}\log
q_{i}-(1-q_{i})\log(1-q_{i}).\end{split}$ (23)

This and similar reparameterizations of $\mathbb{L}$ and $\mathcal{F}_{B}$ are
not new and have been used for different purposes [Welling and Teh, 2001,
Mooij and Kappen, 2005]. One result from Watanabe and Fukumizu [2009] that is
particularly relevant for our work is the following:

######  Proposition 1 (Corollary 2 in Watanabe and Fukumizu [2009])

$\mathcal{F}_{B}$ is convex on $\mathbb{L}$, if and only if the graph
$\mathbf{G}$ is either a tree or contains precisely one cycle.

In principle, this statement fully characterizes the convexity of
$\mathcal{F}_{B}$. As it is only applicable for two types of graphs, convexity
might not seem to be a suitable concept for explaining the success (or
failure) of the Bethe approximation 151515In particular, there are popular
applications involving graphs with multiple cycles where the Bethe
approxiation (resp., LBP) is highly accurate [Frey and MacKay, 1997, Yanover
and Weiss, 2002].. However, we will show how we can refine Proposition 1 by
defining an appropriate subset $\mathbb{B}$ of $\mathbb{L}$, on which
$\mathcal{F}_{B}$ can preserve its convexity for general graphs. Moreover, we
claim that $\mathcal{F}_{B}$ provides accurate solutions if is convex on that
specific subset $\mathbb{B}$ of $\mathbb{L}$. From this, three questions
arise:

  1. (i)

How can we characterize $\mathbb{B}$ such that $\mathcal{F}_{B}$ can preserve
its convexity for general graphs?

  2. (ii)

How can we determine if $\mathcal{F}_{B}$ is convex on $\mathbb{B}$?

  3. (iii)

Can we explain the reliability of the Bethe approximation through its
convexity on $\mathbb{B}$?

For the rest of this work, we focus on answering these three questions: (i)
will be addressed in the remainder of the current subsection; (ii) will be
subject of Sec. 3.2 and Sec. 3.3; and (iii) will be theoretically discussed in
Sec. 3.4 and experimentally evaluated in Sec. 4.

First we identify a subset $\mathbb{B}$ of $\mathbb{L}$ on which
$\mathcal{F}_{B}$ can preserve its convexity for general graphs. For that
purpose, recall that we are primarily interested in minimizing the Bethe free
energy to approximate the minimum of the Gibbs free energy (Sec. 2.3). If we
exclude regions of $\mathbb{L}$ where no minimum of $\mathcal{F}_{B}$ can lie,
then we lose no relevant information on its accuracy. Or, in other words, we
expect that $\mathcal{F}_{B}$ should primarily be a convex approximation to
the Gibbs free energy on a subset $\mathbb{B}$ of $\mathbb{L}$, that includes
all minima of $\mathcal{F}_{B}$. Following this idea, we adopt an approach
made by Welling and Teh [2001], Shin [2012], Weller and Jebara [2013] – albeit
in their work, to reduce the computational complexity of Bethe minimization:

Let $(\bm{q};\bm{\xi})$ be a point in $\mathbb{L}$ where all partial
derivatives of $\mathcal{F}_{B}$ with respect to the edge-specific variables
$\xi_{ij}$ (not necessarily with respect to the node-specific variables
$q_{i}$) are zero.; i.e., the gradient of $\mathcal{F}_{B}$ evaluated in that
point is
$\nabla\mathcal{F}_{B}(\bm{q};\bm{\xi})=\big{(}\frac{\partial\mathcal{F}_{B}}{\partial
q_{1}}(\bm{q};\bm{\xi}),\dots,\frac{\partial\mathcal{F}_{B}}{\partial
q_{N}}(\bm{q};\bm{\xi});0,\dots,0\big{)}^{T}$. This imposes a necessary
condition on any $(\bm{q};\bm{\xi})\in\mathbb{L}$ for being a stationary point
(and thus a minimum) of $\mathcal{F}_{B}$ (Sec. 2.1). In Welling and Teh
[2001] it has been shown that, if we choose an arbitrary (but fixed) $\bm{q}$,
there exists a unique $\bm{\xi}^{\ast}(\bm{q})$ – i.e., depending on the
choice of $\bm{q}$ – such that $\big{(}\bm{q};\bm{\xi}^{\ast}(\bm{q})\big{)}$
is a candidate for a stationary point of $\mathcal{F}_{B}$; moreover, the
explicit form of $\bm{\xi}^{\ast}(\bm{q})$ has been computed with its
individual components given by

$\displaystyle\begin{split}&\xi_{ij}^{\ast}(q_{i},q_{j})=\frac{1}{2\alpha_{ij}}\Big{(}Q_{ij}-\sqrt[\
]{Q_{ij}^{2}-4\alpha_{ij}(1+\alpha_{ij})q_{i}q_{j}}\,\,\Big{)},\\\
&\text{where}\,\,\quad\alpha_{ij}=e^{4\beta
J_{ij}}-1\quad\,\,\text{and}\quad\,\,Q_{ij}=1+\alpha_{ij}(q_{i}+q_{j}).\end{split}$
(24)

In other words, for any $\bm{q}$ there exists at most one
$\bm{\xi}^{\ast}(\bm{q})$ (defined in (24)) such that
$\big{(}\bm{q};\bm{\xi}^{\ast}(\bm{q})\big{)}$ can be a minimum of
$\mathcal{F}_{B}$. Consequently, we address question (i) according to our
previously described requirements on $\mathbb{B}$, by defining a suitable
subset $\mathbb{B}$ of $\mathbb{L}$ as follows:

######  Theorem 2 (Adopted from Welling and Teh [2001])

Every minimum of $\mathcal{F}_{B}$ is located on an $\lvert\mathbf{X}\rvert$\-
dimensional submanifold $\mathbb{B}$ of $\mathbb{L}$, the Bethe box, which is
defined as

$\displaystyle\begin{split}\mathbb{B}\coloneqq\\{&(\bm{q};\bm{\xi}^{\ast}(\bm{q}))\in\mathbb{L}:\,\,0<q_{i}<1,i\in\mathbf{X};\\\
&\,\,\xi_{ij}^{\ast}(q_{i},q_{j})\,\,\,\text{given
by}\,\,\,\eqref{eq:xi_optimal},\,(i,j)\in\mathbf{E}\\}.\end{split}$ (25)

Note that $\mathbb{B}$ is parameterized in terms of the $q_{i}$ only, while
the $\xi_{ij}$ have lost their role as independent parameters but depend on
the $q_{i}$ through the bijective mapping (24). This implies that $\mathbb{B}$
is isomorphic to the open unit cube $(0,1)^{\lvert\mathbf{X}\rvert}$ of
dimension $\lvert\mathbf{X}\rvert$, which explains the naming of $\mathbb{B}$.
Figure 3 shows the Bethe box with respect to a graph on two vertices and a
single edge for different values of the associated coupling $J_{ij}$ and
inverse temperature $\beta$.

Figure 3: Bethe box $\mathbb{B}\subset\mathbb{L}$ for different values of
coupling $J_{ij}$ and inverse temperature $\beta$ on a single-edge graph. For
a ferromagnetic edge ($J_{ij}>0$; first row), $\mathbb{B}$ approaches the
upper boundary of $\mathbb{L}$ as $\beta$ increases; for an antiferromagnetic
edge ($J_{ij}<0$; second row), $\mathbb{B}$ approaches the lower boundary of
$\mathbb{L}$ as $\beta$ increases.

The Bethe box $\mathbb{B}$ satisfies our intuitive requirements on a subset of
$\mathbb{L}$, on which $\mathcal{F}_{B}$ can preserve its convexity for
general graphs (in fact, we prove this in Sec. 3.2 and Sec. 3.3). Also, the
definition (25) of $\mathbb{B}$ once again reduces the number of independent
variables $\mathcal{F}_{B}$, as we can define $\mathcal{F}_{B}$ directly on
the space of its potential minima, by substituting
$\xi_{ij}^{\ast}(q_{i},q_{j})$ from (24) into $\xi_{ij}$ in (21). An immediate
implication of Theorem 2 is:

###### Corollary 3

If $\mathcal{F}_{B}$ is strictly convex on $\mathbb{B}$, then it has a unique
global minimum in $\mathbb{L}$.

In Sec. 3.4, we will have a closer look at the connection between convexity of
$\mathcal{F}_{B}$ uniqueness of a minimum. Finally, we require the first-order
and second-order derivatives of $\mathcal{F}_{B}$ on $\mathbb{B}$:

######  Lemma 4 ((a): Welling and Teh [2001]; (b): Lemma 9 in Weller and
Jebara [2013])

  1. (a)

The first-order partial derivatives of $\mathcal{F}_{B}$ on $\mathbb{B}$ are

$\displaystyle\begin{split}\frac{\partial\mathcal{F}_{B}}{\partial
q_{i}}=-2\,\theta_{i}+2\sum\limits_{j\in\mathcal{N}(i)}J_{ij}+\frac{1}{\beta}\log\Big{(}\frac{(1-q_{i})^{d_{i}-1}}{q_{i}^{d_{i}-1}}\prod\limits_{j\in\mathcal{N}(i)}\frac{q_{i}-\xi_{ij}^{\ast}}{1+\xi_{ij}^{\ast}-q_{i}-q_{j}}\Big{)}.\end{split}$
(26)

  2. (b)

The second-order partial derivatives of $\mathcal{F}_{B}$ on $\mathbb{B}$ are

$\displaystyle\begin{split}&\frac{\partial^{2}\mathcal{F}_{B}}{\partial
q_{i}\partial
q_{j}}=\begin{dcases}\frac{1}{\beta}\,\Big{(}\\!-\frac{d_{i}-1}{q_{i}(1-q_{i})}+\sum\limits_{j\in\mathcal{N}(i)}\frac{q_{j}(1-q_{j})}{T_{ij}}\Big{)}&\quad
i=j\\\
\frac{1}{\beta}\,\Big{(}\frac{q_{i}q_{j}-\xi_{ij}^{\ast}}{T_{ij}}\Big{)}&\quad
j\in\mathcal{N}(i)\\\ 0&\quad\text{otherwise},\end{dcases}\end{split}$ (27)
$\displaystyle\text{where}\quad T_{ij}(q_{i},q_{j})\coloneqq
q_{i}q_{j}(1-q_{i})(1-q_{j})-(\xi_{ij}^{\ast}(q_{i},q_{j})-q_{i}q_{j})^{2}>0.$
(28)

Also, we have that $\frac{\partial^{2}\mathcal{F}_{B}}{\partial q_{i}^{2}}$ is
always $>0$.

A simple, but remarkable consequence from Lemma 4 is:

###### Corollary 5

The Hessian matrix of $\mathcal{F}_{B}$ is independent from $\theta_{i}$. In
other words, the curvature of the Bethe free energy is not influenced by the
local potentials.

Next, we deal with question (ii): How can we determine, if $\mathcal{F}_{B}$
is convex on the Bethe box $\mathbb{B}$? To address this problem, we propose
two different methodologies – the first one in Sec. 3.2, the second one in
Sec. 3.3 – that both rely on a profound analysis of the Hessian matrix, i.e.,
the curvature of the Bethe free energy. The basic idea of both approaches is
the same: to apply the formal definition of strict convexity (Sec. 2.1) to the
Bethe free energy: Let
$\textbf{H}_{B}(\bm{q}):\mathbb{B}\to\mathbb{R}^{\lvert\mathbf{X}\rvert\times\lvert\mathbf{X}\rvert}$
be the Hessian matrix of $\mathcal{F}_{B}$ (the _Bethe Hessian_), considered
as a matrix-valued function on the Bethe box 161616To simplify the notation,
we denote the Bethe Hessian by $\textbf{H}_{B}$ instead of
$\textbf{H}_{\mathcal{F}_{B}}$, as formally defined in Sec. 2.1. Also, we
mostly omit the dependency of $\textbf{H}_{B}$ on $\bm{q}$, since this is
usually clear from the context.. Then $\mathcal{F}_{B}$ is strictly convex on
$\mathbb{B}$, if $\textbf{H}_{B}(\bm{q})$ is positive definite for all
$\bm{q}\in\mathbb{B}$. Unfortunately, the complex analytical form of
$\textbf{H}_{B}$ (Lemma 4, (b)) prevents a direct application of this
definition and requires relaxations to the problem.

###  3.2 Node-Specific Convexity Criterion: Diagonal Dominance of Bethe
Hessian

Our first approach to relax the convexity of $\mathcal{F}_{B}$ uses the
concept of diagonal dominance, a sufficient condition for positive
definiteness of a matrix (Sec. 2.1). The main result is Theorem 15, which
reduces diagonal dominance of the Bethe Hessian to the positivity of a set of
one-dimensional polynomials. All proofs for this section are collected in
Appendix A.

Let us compactly write the $(i,j)$-th entry of $\textbf{H}_{B}$ (i.e.,
$\frac{\partial^{2}\mathcal{F}_{B}}{\partial q_{i}\partial q_{j}}$) as
$h_{ij}$. Then the diagonal dominance criterion (equation (1)) says that
$\textbf{H}_{B}$ is positive definite in $\bm{q}\in\mathbb{B}$, if

$\displaystyle h_{ii}-\sum\limits_{j\in\mathcal{N}(i)}|h_{ij}|>0$ (29)

for each row of $\textbf{H}_{B}$, i.e., for all $i\in\mathbf{X}$. Note that
$h_{ii}>0$, and $h_{ij}=0$ if neither $i=j$ nor $j\in\mathcal{N}(i)$. If we
insert the expressions from (LABEL:eq:Bethe_second_derivatives) into (29),
this is equivalent to

$\displaystyle\frac{1}{\beta}\,\Big{(}\\!-\frac{d_{i}-1}{q_{i}(1-q_{i})}+\sum\limits_{j\in\mathcal{N}(i)}\frac{q_{j}(1-q_{j})}{T_{ij}}-\sum\limits_{j\in\mathcal{N}(i)}\,\frac{|q_{i}q_{j}-\xi_{ij}^{\ast}|}{T_{ij}}\Big{)}>0,$
(30)

because $T_{ij}>0$ (Lemma 4, (b)). To compute the absolute values, we can use
the following result:

######  Lemma 6 (Lemma 2 in Weller and Jebara [2013])

If $(i,j)$ is a ferromagnetic edge ($J_{ij}>0$), then
$\xi_{ij}^{\ast}>q_{i}q_{j}$. If $(i,j)$ is an antiferromagnetic edge
($J_{ij}<0$), then $\xi_{ij}^{\ast}<q_{i}q_{j}$.

If we distinguish between the two types of edges, inequality (30) becomes

$\displaystyle-\frac{d_{i}-1}{q_{i}(1-q_{i})}+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i)\\\
J_{ij}>0\end{subarray}}\frac{q_{j}(1-q_{j})-\xi_{ij}^{\ast}+q_{i}q_{j}}{T_{ij}}+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i)\\\
J_{ij}<0\end{subarray}}\frac{q_{j}(1-q_{j})-q_{i}q_{j}+\xi_{ij}^{\ast}}{T_{ij}}>0.$
(31)

So far, we have only considered one specific point $\bm{q}\in\mathbb{B}$;
however, to guarantee the convexity of $\mathcal{F}_{B}$ on the entire Bethe
box $\mathbb{B}$, inequality (31) must be valid for _all_
$\bm{q}\in\mathbb{B}$ or, equivalently, for the infimum 171717Strictly
speaking, the infimum could also be equal to zero, while positivity holds for
all interior points of $\mathbb{B}$. However, we will ignore this subtlety as
it is does not affect our later results. over all $\bm{q}\in\mathbb{B}$. In
summary, we arrive at the following problem:

###### Proposition 7

Let us denote the expression on the left-hand side of inequality (31) by
$R_{i}(\bm{q})$. Then $\mathcal{F}_{B}$ is convex on the Bethe box
$\mathbb{B}$ if, for all nodes $i\in\mathbf{X}$,

$\displaystyle\inf\limits_{\bm{q}\in\mathbb{B}}\,\,R_{i}(\bm{q})>0.$ (32)

In Proposition 7 we have reformulated diagonal dominance of the Bethe Hessian
in terms of an optimization problem. The infimum of $R_{i}(\bm{q})$ measures
to what extent the diagonal dominance criterion is satisfied (or violated) in
node $i$. As $R_{i}$ is a nonlinear function of potentially many variables
181818More precisely, $R_{i}$ depends on $(d_{i}+1)$ variables: on $q_{i}$
that is associated to node $i$, and additionally on all $q_{j}$ that are
associated to neighboring nodes $j\in\mathcal{N}(i)$., its solution (i.e., the
verification of inequality (32)) is not trivial.

In the remainder of this section, we analyze the left-hand side of the
inequality (32) and finally present a condition that is equivalent to (32) but
much easier to verify (Theorem 15). Our first observation is that
$R_{i}(\bm{q})$ has the special form

$\displaystyle R_{i}(\bm{q})$
$\displaystyle=r_{i}(q_{i})+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i)\\\
J_{ij}>0\end{subarray}}r_{ij}^{+}(q_{i},q_{j})+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i)\\\
J_{ij}<0\end{subarray}}r_{ij}^{-}(q_{i},q_{j}),\quad\text{where}$ (33)
$\displaystyle
r_{i}(q_{i})\coloneqq-\frac{d_{i}-1}{q_{i}(1-q_{i})}\quad\text{depends only
on}\,\,q_{i},$ (34) $\displaystyle
r_{ij}^{+}(q_{i},q_{j})\coloneqq\frac{q_{j}(1-q_{j})-\xi_{ij}^{\ast}+q_{i}q_{j}}{T_{ij}}\quad\text{depends
on}\,\,q_{i}\,\,\text{and}\,\,q_{j},\,\,\text{and}$ (35) $\displaystyle
r_{ij}^{-}(q_{i},q_{j})\coloneqq\frac{q_{j}(1-q_{j})-q_{i}q_{j}+\xi_{ij}^{\ast}}{T_{ij}}\quad\text{depends
on}\,\,q_{i}\,\,\text{and}\,\,q_{j}.$ (36)

Note that, for each ferromagnetic edge we have one term of the form (35), and
for each antiferromagnetic edge we have one term of the form (36) (with
$q_{i}$ and $q_{j}$ being the variables associated to the end nodes of a
specific edge $(i,j)$). In particular, $q_{i}$ (associated to node $i$) is the
only variable that appears in each individual sum term, while any other
variable $q_{j}$ (for $j\in\mathcal{N}(i)$) appears in precisely one of the
pairwise terms $r_{ij}^{+}(q_{i},q_{j})$ or $r_{ij}^{-}(q_{i},q_{j})$
(depending on the type of the edge that connects $i$ and $j$). Figure 4
illustrates this situation.

Figure 4: Dependencies and components of the function $R_{i}(\bm{q})$
associated to some node $i$. Ferromagnetic edges are drawn in green and
antiferromagnetic edges are drawn in red.

We can exploit the specific structure (33) of $R_{i}(\bm{q})$ and decompose
the computation of $\inf\limits_{\bm{q}\in\mathbb{B}}R_{i}(\bm{q})$ into a
collection of simpler sub-problems, each involving only a single variable.
Concretely, we apply a two-step procedure: In the first step, we fix $q_{i}$
and determine

$\displaystyle\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{+}(q_{i},q_{j})\,\,\quad\,\,\text{and}\quad\,\,\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{-}(q_{i},q_{j})$
(37)

for each of the pairwise terms
$r_{ij}^{+}(q_{i},q_{j}),r_{ij}^{-}(q_{i},q_{j})$ separately (with respect to
$q_{j}$ as a function of $q_{i}$). In the second step, we fix the solutions to
these sub-problems and compute the infimum of the overall function
$R_{i}(\bm{q})$ with respect to $q_{i}$. Altogether, we determine

$\displaystyle\tilde{R}_{i}\coloneqq\inf\limits_{q_{i}\in(0,1)}\,\,\Big{[}\,\,\,r_{i}(q_{i})+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i):\\\
J_{ij}>0\end{subarray}}\,\,\Big{(}\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{+}(q_{i},q_{j})\Big{)}+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i):\\\
J_{ij}<0\end{subarray}}\,\,\Big{(}\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{-}(q_{i},q_{j})\Big{)}\,\,\,\Big{]}$
(38)

which is equal to the infimum of $R_{i}(\bm{q})$ as the following Lemma shows:

###### Lemma 8

Let $R_{i}(\bm{q})$ be the expression 191919We remark that our proof of Lemma
8 does not require the specific definition of the functions
$r_{i},r_{ij}^{+},r_{ij}^{-}$, but only relies on the structure of the
decomposition (33) of $R_{i}(\bm{q})$.  in (33) and let $\tilde{R}_{i}$ be the
expression in (38). Then

$\displaystyle\inf\limits_{\bm{q}\in\mathbb{B}}\,\,R_{i}(\bm{q})=\tilde{R}_{i}.$
(39)

Lemma 8 justifies the described two-step procedure for evaluating
$\inf\limits_{\bm{q}\in\mathbb{B}}R_{i}(\bm{q})$ from ’inside’ (fix
$q_{i}\in(0,1)$ and compute the infima of $r_{ij}^{+}(q_{i},q_{j})$ and
$r_{ij}^{-}(q_{i},q_{j})$ with respect to $q_{j}$ as a function of $q_{i}$) to
’outside’ (given the solutions to these sub-problems, compute the infimum of
$R_{i}(\bm{q})$ with respect to $q_{i}$). Hence, we must first solve the
’inner’ sub-problems (37) whose solutions are required in the ’outer’ problem
(38). We have two types of sub-problems to be solved: one for ferromagnetic
edges and one for antiferromagnetic edges. Instead of treating them
separately, we show that they can be transformed into each other:

###### Lemma 9

Let $J_{ij}>0$, and let $r_{ij}^{+}(q_{i},q_{j})$ and
$r_{ij}^{-}(q_{i},q_{j})$ be the two types of pairwise functions in the sum
decomposition (33) of $R_{i}(\bm{q})$, with $r_{ij}^{+}(q_{i},q_{j})$
implicitly depending on $J_{ij}$ (corresponding to a ferromagnetic edge) and
$r_{ij}^{-}(q_{i},q_{j})$ implicitly depending on $-J_{ij}$ (corresponding to
an antiferromagnetic edge). Then we have

$\displaystyle
r_{ij}^{-}(q_{i},q_{j})=r_{ij}^{+}(q_{i},1-q_{j})=r_{ij}^{+}(1-q_{i},q_{j}).$
(40)

In other words, if we reflect $r_{ij}^{+}(q_{i},q_{j})$ in its domain
$(0,1)\times(0,1)$ either across the line $q_{i}=0.5$ or the line $q_{j}=0.5$,
we obtain $r_{ij}^{-}(q_{i},q_{j})$ (and vice versa).

Lemma 9 establishes a direct relationship between $r_{ij}^{-}(q_{i},q_{j})$
and $r_{ij}^{+}(q_{i},q_{j})$. It therefore facilitates the determination of
$\tilde{R}_{i}$ in (38) by focusing on a single edge type, instead of
distinguishing between two types of edges. Without loss of generality, we may
replace each pairwise term $r_{ij}^{-}(q_{i},q_{j})$ (which implicitly depends
on $J_{ij}<0$) in (38) with $r_{ij}^{+}(1-q_{i},q_{j})$ and the sign-inverted
coupling $|J_{ij}|>0$. Under this assumption, (38) becomes

$\displaystyle\tilde{R}_{i}\coloneqq\inf\limits_{q_{i}\in(0,1)}\,\,\Big{[}\,\,\,r_{i}(q_{i})+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i):\\\
J_{ij}>0\end{subarray}}\,\,\Big{(}\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{+}(q_{i},q_{j})\Big{)}+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i):\\\
J_{ij}<0\end{subarray}}\,\,\Big{(}\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{+}(1-q_{i},q_{j})\Big{)}\,\,\,\Big{]}.$
(41)

We now return to the ’inner’ sub-problems stated in (37) with our focus on the
first problem type: the determination of
$\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{+}(q_{i},q_{j})$ with $q_{i}\in(0,1)$
arbitrary but fixed 202020By Lemma 9, the determination of
$\inf\limits_{q_{j}\in(0,1)}\,\,r_{ij}^{-}(q_{i},q_{j})$ can be reduced to the
ferromagnetic case.. In particular, we consider $r_{ij}^{+}(q_{i},q_{j})$ as a
one-dimensional function of $q_{j}$ defined on the interval $(0,1)$
(representing the line segment that connects $(q_{i},0)$ and $(q_{i},1)$ in
the domain $(0,1)\times(0,1)$). To find the greatest lower bound on
$r_{ij}^{+}$ in that interval, we first show that $r_{ij}^{+}$ has a unique
stationary point in $(0,1)$ (Theorem 11) which, however, is a maximum (Theorem
13). This implies that the infimum is not taken by $r_{ij}^{+}$ but exists as
a limit value at the boundary of $(0,1)$. In Lemma 12 and Corollary 14, we
derive the explicit formula for the infimum, in dependence of the chosen
$q_{i}$. Figure 5 shows sketches of $r_{ij}^{+}(q_{i},q_{j})$ for different
values of $q_{i}$.

Figure 5: $r_{ij}^{+}(q_{i},q_{j})$ as a function of $q_{j}\in(0,1)$ and for
different fixed values of $q_{i}$. The green line $q_{j}=1-q_{i}$ represents
the space of all stationary points of $r_{ij}^{+}(q_{i},q_{j})$ (exactly one
for each $q_{i}\in(0,1)$, Theorem 11), which are always maxima (drawn as red
points, Theorem 13). The blue points at the boundary represent the infima of
$r_{ij}^{+}(q_{i},q_{j})$ with respect to $q_{j}$, which are never taken by
$r_{ij}^{+}(q_{i},q_{j})$ (Lemma 12 and Corollary 14). Note that $q_{i}=0.5$
represents the only scenario where the infimum exists at both boundary points
$0$ and $1$.

For proving Theorem 11, we require the first-order partial derivatives of the
functions $\xi_{ij}^{\ast},T_{ij},$ and $r_{ij}^{+}$ with respect to $q_{i}$
and $q_{j}$. Moreover, we prepare the second-order partial derivatives of
$\xi_{ij}^{\ast}$ and $T_{ij}$, that we will require in the proof of Theorem
21 (Sec. 3.3).

###### Lemma 10

  * (a)

The first-order partial derivatives of
$\xi_{ij}^{\ast}(q_{i},q_{j}),T_{ij}(q_{i},q_{j}),r_{ij}^{+}(q_{i},q_{j})$ are

$\displaystyle\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{i}}=\frac{\alpha_{ij}(q_{j}-\xi_{ij}^{\ast})+q_{j}}{1+\alpha_{ij}(q_{i}+q_{j}-2\xi_{ij}^{\ast})},$
(42) $\displaystyle\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}=\frac{\alpha_{ij}(q_{i}-\xi_{ij}^{\ast})+q_{i}}{1+\alpha_{ij}(q_{i}+q_{j}-2\xi_{ij}^{\ast})},$
(43) $\displaystyle\frac{\partial T_{ij}}{\partial
q_{i}}=(1-2q_{i})q_{j}(1-q_{j})-2(\xi_{ij}^{\ast}-q_{i}q_{j})(\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{i}}-q_{j}),$ (44) $\displaystyle\frac{\partial T_{ij}}{\partial
q_{j}}=(1-2q_{j})q_{i}(1-q_{i})-2(\xi_{ij}^{\ast}-q_{i}q_{j})(\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}-q_{i}),$ (45) $\displaystyle\frac{\partial r_{ij}^{+}}{\partial
q_{i}}=\frac{(-\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{i}}+q_{j})T_{ij}-(q_{j}(1-q_{j})-\xi_{ij}^{\ast}+q_{i}q_{j})\frac{\partial
T_{ij}}{\partial q_{i}}}{T_{ij}^{2}},$ (46) $\displaystyle\frac{\partial
r_{ij}^{+}}{\partial
q_{j}}=\frac{(1-2q_{j}-\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}+q_{i})T_{ij}-(q_{j}(1-q_{j})-\xi_{ij}^{\ast}+q_{i}q_{j})\frac{\partial
T_{ij}}{\partial q_{j}}}{T_{ij}^{2}}.$ (47)

  * (b)

The second-order partial derivatives of $\xi_{ij}^{\ast}(q_{i},q_{j})$ and
$T_{ij}(q_{i},q_{j})$ are

$\displaystyle\frac{\partial^{2}\xi_{ij}^{\ast}}{\partial
q_{i}^{2}}=\alpha_{ij}\frac{\big{(}-\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{i}}\big{)}\cdot
A_{ij}-\big{(}\alpha_{ij}(q_{j}-\xi_{ij}^{\ast})+q_{j}\big{)}\big{(}1-2\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{i}}\big{)}}{A_{ij}^{2}},$ (48)
$\displaystyle\frac{\partial^{2}\xi_{ij}^{\ast}}{\partial q_{i}\partial
q_{j}}=\frac{\alpha_{ij}\big{[}\big{(}1-\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}+\frac{1}{\alpha_{ij}}\big{)}\cdot
A_{ij}-\big{(}\alpha_{ij}(q_{j}-\xi_{ij}^{\ast})+q_{j}\big{)}\big{(}1-2\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}\big{)}\big{]}}{A_{ij}^{2}},$ (49)
$\displaystyle\frac{\partial^{2}\xi_{ij}^{\ast}}{\partial
q_{j}^{2}}=\alpha_{ij}\frac{\big{(}-\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}\big{)}\cdot
A_{ij}-\big{(}\alpha_{ij}(q_{i}-\xi_{ij}^{\ast})+q_{i}\big{)}\big{(}1-2\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}\big{)}}{A_{ij}^{2}},$ (50)
$\displaystyle\frac{\partial^{2}T_{ij}}{\partial
q_{i}^{2}}=-2q_{j}(1-q_{j})-2\Big{(}\big{(}\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{i}}-q_{j}\big{)}^{2}+(\xi_{ij}^{\ast}-q_{i}q_{j})\big{(}\frac{\partial^{2}\xi_{ij}^{\ast}}{\partial
q_{i}^{2}}\big{)}\Big{)},$ (51)
$\displaystyle\frac{\partial^{2}T_{ij}}{\partial q_{i}\partial
q_{j}}=(1-2q_{j})(1-2q_{i})-2\Big{(}\big{(}\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{i}}-q_{j}\big{)}\big{(}\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}-q_{i}\big{)}+(\xi_{ij}^{\ast}-q_{i}q_{j})\big{(}\frac{\partial^{2}\xi_{ij}^{\ast}}{\partial
q_{i}\partial q_{j}}-1\big{)}\Big{)},$ (52)
$\displaystyle\frac{\partial^{2}T_{ij}}{\partial
q_{j}^{2}}=-2q_{i}(1-q_{i})-2\Big{(}\big{(}\frac{\partial\xi_{ij}^{\ast}}{\partial
q_{j}}-q_{i}\big{)}^{2}+(\xi_{ij}^{\ast}-q_{i}q_{j})\big{(}\frac{\partial^{2}\xi_{ij}^{\ast}}{\partial
q_{j}^{2}}\big{)}\Big{)},$ (53)

where we have defined
$A_{ij}\coloneqq\big{(}1+\alpha_{ij}(q_{i}+q_{j}-2\xi_{ij}^{\ast})\big{)}$.

###### Theorem 11

Let $q_{i}$ be an arbitrary (but fixed) value in $(0,1)$ and consider
$r_{ij}^{+}(q_{i},q_{j})$ as a function of $q_{j}$ in the interval $(0,1)$.
The only stationary point of $r_{ij}^{+}(q_{i},q_{j})$ in that interval is
$q_{j}^{\ast}=1-q_{i}$. In other words, the nonlinear equation

$\frac{\partial r_{ij}^{+}}{\partial q_{j}}=0$ (54)

has the unique solution $q_{j}^{\ast}=1-q_{i}$.

Next, we analyze what kind of stationary point $r_{ij}^{+}(q_{i},q_{j})$ has
in $q_{j}^{\ast}=1-q_{i}$ and if this depends on the choice of $q_{i}$. We
answer this in Theorem 13, by comparing the value of $r_{ij}^{+}$ in
$q_{j}^{\ast}$ to its limit values at the boundary points. These are computed
in the next Lemma:

###### Lemma 12

Let $q_{i}\in(0,1)$ be arbitrary (but fixed)- and consider
$r_{ij}^{+}(q_{i},q_{j})$ as a function of $q_{j}$ in $(0,1)$.

  1. (a)

In its unique stationary point $q_{j}^{\ast}=1-q_{i}$, $r_{ij}^{+}$ takes the
value

$\displaystyle
r_{ij}^{+}(q_{i},1-q_{i})=\frac{1}{\xi_{ij}^{\ast}(q_{i},1-q_{i})}.$ (55)

  2. (b)

If $r_{ij}^{+}$ approaches the boundary points of $(0,1)$, its limit values
are

$\displaystyle\lim\limits_{q_{j}\to
0}r_{ij}^{+}(q_{i},q_{j})=\frac{1+\alpha_{ij}q_{i}^{2}}{(1+\alpha_{ij}q_{i})q_{i}(1-q_{i})},\quad\text{and}$
(56) $\displaystyle\lim\limits_{q_{j}\to
1}r_{ij}^{+}(q_{i},q_{j})=\frac{1+\alpha_{ij}(1-q_{i})^{2}}{\big{(}1+\alpha_{ij}(1-q_{i})\big{)}q_{i}(1-q_{i})}.$
(57)

This further implies the relation

$\displaystyle\lim\limits_{q_{j}\to
1}r_{ij}^{+}(q_{i},q_{j})=\lim\limits_{q_{j}\to 0}r_{ij}^{+}(1-q_{i},q_{j}).$
(58)

For a more compact notation, let us define

$\displaystyle r_{ij}^{(+,0)}(q_{i})\coloneqq\lim\limits_{q_{j}\to
0}r_{ij}^{+}(q_{i},q_{j}),\quad\text{and}$ (59) $\displaystyle
r_{ij}^{(+,1)}(q_{i})\coloneqq\lim\limits_{q_{j}\to
1}r_{ij}^{+}(q_{i},q_{j}).$ (60)

Note that $r_{ij}^{(+,0)}$ and $r_{ij}^{(+,1)}$ are not functions of $q_{j}$
any longer but only depend on $q_{i}$.

The next Theorem shows that $q_{j}^{\ast}=1-q_{i}$ is always a maximum of
$r_{ij}^{+}(q_{i},q_{j})$ (still considered as a function of $q_{j}$).

###### Theorem 13

We have the following relations between the values of
$r_{ij}^{+}(q_{i},q_{j})$ in its stationary point $q_{j}^{\ast}=1-q_{i}$ and
its limit values at the boundary:

  * •

If $0<q_{i}<0.5$, then
$r_{ij}^{(+,0)}(q_{i})<r_{ij}^{(+,1)}(q_{i})<r_{ij}^{+}(q_{i},1-q_{i})$.

  * •

If $0.5<q_{i}<1$, then
$r_{ij}^{(+,1)}(q_{i})<r_{ij}^{(+,0)}(q_{i})<r_{ij}^{+}(q_{i},1-q_{i})$.

  * •

If $q_{i}=0.5$, then
$r_{ij}^{(+,0)}(q_{i})=r_{ij}^{(+,1)}(q_{i})<r_{ij}^{+}(q_{i},1-q_{i})$.

In particular, for any choice of $q_{i}\in(0,1)$, the unique stationary point
$q_{j}^{\ast}=1-q_{i}$ is a maximum of $r_{ij}^{+}(q_{i},q_{j})$.

Let us recap the insights from the previous statements: To find the infimum of
$r_{ij}^{+}(q_{i},q_{j})$ with respect to $q_{j}\in(0,1)$ as the variable and
$q_{i}\in(0,1)$ arbitrary but fixed, we have first shown that this function
has the unique stationary point $q_{j}^{\ast}=1-q_{i}$ (Theorem 11). We argued
that this must be a maximum (Theorem 13), which implies that the infimum is
not taken by $r_{ij}^{+}(q_{i},q_{j})$ in the interior of $(0,1)$, but
approached at one of the boundary points (i.e, for $q_{j}\to 0$ or for
$q_{j}\to 1$). We have derived an analytical formula for the limit values of
$r_{ij}^{+}(q_{i},q_{j})$ at these boundary points as a function of $q_{i}$
(Lemma 12), and determined which of the two limit values is less than the
other, depending on the choice of $q_{i}\in(0,1)$. Together with the other
observations, we conclude that this must then be the infimum.

###### Corollary 14

For $q_{i}\in(0,1)$ arbitrary (but fixed), we have

$\displaystyle\inf\limits_{q_{j}\in(0,1)}r_{ij}^{+}(q_{i},q_{j})=\begin{dcases}\frac{1+\alpha_{ij}q_{i}^{2}}{(1+\alpha_{ij}q_{i})q_{i}(1-q_{i})},\quad&\text{if}\quad
0<q_{i}\leq 0.5,\quad\text{and}\\\
\frac{1+\alpha_{ij}(1-q_{i})^{2}}{\big{(}1+\alpha_{ij}(1-q_{i})\big{)}q_{i}(1-q_{i})},\quad&\text{if}\quad
0.5<q_{i}<1.\end{dcases}$

Figure 6 illustrates $\inf\limits_{q_{j}\in(0,1)}r_{ij}^{+}(q_{i},q_{j})$ as a
function of $q_{i}$ (drawn in solid blue curves).

Figure 6: Visualization of the limit values of $r_{ij}^{+}(q_{i},q_{j})$ at
the boundary points $q_{j}=0$ and $q_{j}=1$ as a function of $q_{i}$. The
solid blue curves indicate the position of the infimum, which is approached
for $q_{j}\to 0$ if $q_{i}\leq 0.5$, and for $q_{j}\to 1$ if $q_{i}\geq 0.5$
(Corollary 14). The dashed blue lines represent the level of the infimum for
different fixed values of $q_{i}$.

Corollary 14 provides the solution to the ’inner’ sub-problems stated in (37).
Hence, we return to the ’outer’ problem (41) which is the evaluation of
$\tilde{R_{i}}$. Before we substitute the formulas derived in Corollary 14 and
the expression for $r_{i}(q_{i})$ from (34) into (41), we observe that, for
all $q_{i}\in(0,1)$,

$\displaystyle\inf\limits_{q_{j}\in(0,1)}r_{ij}^{+}(q_{i},q_{j})=\inf\limits_{q_{j}\in(0,1)}r_{ij}^{+}(1-q_{i},q_{j}).$
(61)

This follows directly from the relation (58) from Lemma 12. Consequently,
$\inf\limits_{q_{j}\in(0,1)}r_{ij}^{+}(q_{i},q_{j})$ (as a function of
$q_{i}$) is symmetric across $q_{i}=0.5$, which makes it sufficient to
evaluate $\tilde{R_{i}}$ for $q_{i}\in(0,0.5]$. If we assume that, for each
edge, we only take the coupling strength $|J_{ij}|$ but not the sign of
$J_{ij}$ into account (which is justified by Lemma 9), equation (41) becomes

$\displaystyle\tilde{R}_{i}$
$\displaystyle=\inf\limits_{q_{i}\in(0,0.5)}\,\,\Big{[}\frac{-(d_{i}-1)}{q_{i}(1-q_{i})}+\sum_{\begin{subarray}{c}j\in\mathcal{N}(i)\end{subarray}}\,\,\frac{1+\alpha_{ij}q_{i}^{2}}{(1+\alpha_{ij}q_{i})q_{i}(1-q_{i})}\Big{]}$
(62)
$\displaystyle=\inf\limits_{q_{i}\in(0,0.5)}\,\,\Big{[}\frac{-(d_{i}-1)\prod\limits_{j\in\mathcal{N}(i)}(1+\alpha_{ij}q_{i})+\sum\limits_{j\in\mathcal{N}(i)}\big{(}(1+\alpha_{ij}q_{i}^{2})\prod\limits_{k\in\mathcal{N}(i)\setminus
j}(1+\alpha_{ik}q_{i})\big{)}}{q_{i}(1-q_{i})\prod\limits_{j\in\mathcal{N}(i)}(1+\alpha_{ij}q_{i})}\Big{]}.$
(63)

To check the diagonal dominance criterion (32) stated in Proposition 7, we do
not need to compute the exact value of this infimum (which, by Lemma 8, is
identical to $\inf\limits_{\bm{q}\in\mathbb{B}}\,\,R_{i}(\bm{q})$), but only
have to determine if it is positive. One simple observation further facilitate
this task: All terms in the denominator of (63) are positive; i.e., positivity
of the infimum depends only on the numerator, which is a polynomial in a
single variable $q_{i}$ and of degree $d_{i}+1$.
Altogether, we present the final Theorem of this subsection:

###### Theorem 15

The convexity criterion stated in Proposition 7 is equivalent to the following
condition: The Bethe free energy $\mathcal{F}_{B}$ is convex on the Bethe box
$\mathbb{B}$ if, for all nodes $i\in\mathbf{X}$, the one-dimensional
polynomial

$\displaystyle\Psi_{i}(q_{i})\coloneqq-(d_{i}-1)\prod\limits_{j\in\mathcal{N}(i)}(1+\alpha_{ij}q_{i})+\sum\limits_{j\in\mathcal{N}(i)}\Big{(}(1+\alpha_{ij}q_{i}^{2})\prod\limits_{k\in\mathcal{N}(i)\setminus
j}(1+\alpha_{ik}q_{i})\Big{)}$ (64)

is strictly positive on the interval $(0,0.5]$.

Theorem 15 is the central result of this subsection. It provides a suffient
condition for the convexity of the Bethe free energy $\mathcal{F}_{B}$ on the
Bethe box $\mathbb{B}$. Figure 7 shows the behavior of the polynomial
$\Psi_{i}(q_{i})$ in dependence of the inverse temperature $\beta$.

Figure 7: Development of the polynomial $\Psi_{i}(q_{i})$ for increasing
values of the inverse temperature $\beta$, associated to a node $i$ with three
neighbors ($d_{i}=3$), and a coupling $J_{ij}=J=0.5$ that is shared by all
three incident edges. For $\beta=0.25$ and $\beta=0.55$, condition (64) from
Theorem 15 is satisfied with respect to node $i$ (blue). $\beta^{\ast}\approx
0.693$ is the critical threshold at which (64) is violated for the first time,
as $\Psi_{i}(q_{i})$ touches the horizontal axis and takes root in the
interval $(0,0.5)$ (orange); for higher values of the inverse temperature
(e.g., $\beta=0.75$), $\Psi_{i}(q_{i})$ has multiple roots in $(0,0.5)$ which
causes (64) to be violated (red).

###  3.3 Edge-Specific Convexity Criterion: Sum Decomposition of Bethe Hessian

Our second approach to relax the convexity of $\mathcal{F}_{B}$ decomposes the
Bethe Hessian into a sum of sparse matrices, and derives a sufficient
condition for the positive definiteness of the individual matrices in that
sum. The main result is Theorem 23, a closed formula for the critical
threshold of the inverse temperature, at which convexity of the Bethe free
energy is not guaranteed any longer. All proofs for this section are collected
in Appendix B.

The starting point is the same as in Sec. 3.2, namely to ask under which
conditions the Bethe Hessian $\textbf{H}_{B}(\bm{q})$ is positive definite for
all $\bm{q}$ in the Bethe box $\mathbb{B}$ (the definition of convexity, Sec.
2.1). However, we use an alternative approach to relax this problem: While, in
Sec. 3.2, we have applied the diagonal dominance criterion to the entire Bethe
Hessian, we now decompose $\textbf{H}_{B}$ into a sum of sparse matrices
$\textbf{H}_{B}^{\prime}$ (i.e., by writing
$\textbf{H}_{B}=\sum\textbf{H}_{B}^{\prime}$) and then characterize the
positive definiteness of the individual matrices in that sum. In case that all
$\textbf{H}_{B}^{\prime}$ are positive definite over $\mathbb{B}$, their sum
$\textbf{H}_{B}$ is positive definite over $\mathbb{B}$ as well. This induces
a sufficient (but not necessary 212121Note that $\textbf{H}_{B}$ might still
be positive definite if some (or even all) matrices $\textbf{H}_{B}^{\prime}$
in the sum decomposition are not. This explains why our approach is merely a
relaxation to the problem of characterizing the positive definiteness of
$\textbf{H}_{B}$.) condition for the convexity of $\mathcal{F}_{B}$ on the
Bethe box $\mathbb{B}$. Further differences between the our approaches will be
discussed in Sec. 3.4.

The first problem we are facing is how to define the matrices
$\textbf{H}_{B}^{\prime}$ in a sum decomposition of $\textbf{H}_{B}$. On the
one hand, they should be simple enough to allow for an exact characterization
of their definiteness properties. On the other hand, too many (simple)
matrices in the sum may induce a relatively high gap between the ’true’ and
the predicted convexity (as a sufficient criterion for convexity requires
positive definiteness of each individual matrix in the sum). Instead of trying
to minimize this gap, we propose an edge-specific class of sum decompositions,
in which the individual matrices are naturally related to the edges in the
graph (i.e., for each edge $(i,j)$ we have one associated matrix
$\textbf{H}_{B}^{(i,j)}$ in the sum). This particular choice can be motivated
by the properties of the Bethe free energy, which is itself a sum of edge- and
node-specific statistics (21). While it might be tempting to write the Bethe
Hessian as a sum of Hessians corresponding to these statistics (i.e.,
computing the Hessians of the individual terms in (21) and summing them up),
we must consider one aspect: the local entropies $\mathcal{S}_{i}$ are always
concave and their Hessians can never be positive definite (i.e., considering
them as separate matrices in the sum can never result in a sufficient
criterion for convexity). Also, these are the only node-specific statistics
contributing to the Bethe Hessian, because $\textbf{H}_{B}$ is independent
from the fields $\theta_{i}$ (Corollary 5). Thus, instead of defining separate
’node-specific Hessians’ $\textbf{H}_{B}^{i}$ for each node $i$, it is
preferable to absorb the node-specific statistics from $\textbf{H}_{B}$ into
the ’edge-specific Hessians’ $\textbf{H}_{B}^{(i,j)}$.

We consider the following class of sum decompositions: Let us write the Bethe
Hessian as

$\displaystyle\textbf{H}_{B}=\sum\limits_{(i,j)\in\mathbf{E}}\textbf{H}_{B}^{(i,j)},$
(65)

with the entries of $\textbf{H}_{B}^{(i,j)}$ – for an arbitrary edge $(i,j)$
and indices $k,l\in\\{1,\dots,|\mathbf{X}|\\}$ – being defined as

$\displaystyle\begin{split}\big{(}\textbf{H}_{B}^{(i,j)}\big{)}_{k,l}=\begin{dcases}\frac{1}{\beta}\,\Big{(}\\!-s_{ij}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{T_{ij}}\Big{)}&\quad
k=l=i\\\
\frac{1}{\beta}\,\Big{(}\\!-s_{ji}\,\frac{(d_{j}-1)}{q_{j}(1-q_{j})}+\frac{q_{i}(1-q_{i})}{T_{ij}}\Big{)}&\quad
k=l=j\\\
\frac{1}{\beta}\,\Big{(}\frac{q_{i}q_{j}-\xi_{ij}^{\ast}}{T_{ij}}\Big{)}&\quad
k,j\in\\{i,j\\},\,k\neq l\\\ 0&\quad\text{otherwise}.\end{dcases}\end{split}$
(66)

This definition requires some explanation. First, note that we have introduced
parameters $s_{ij}$ and $s_{ji}$ (i.e, two parameters per edge) whose purpose
is to ’distribute’ the second-order local entropy contributions from
$\textbf{H}_{B}$ (using the equations from
(LABEL:eq:Bethe_second_derivatives), (b)) over the edge-specific Hessians
$\textbf{H}_{B}^{(i,j)}$. More precisely, $s_{ij}$ represents the relative
proportion of the second-order local entropy contribution related to node $i$
(i.e, of $\frac{\partial^{2}}{\partial
q_{i}^{2}}\big{(}-\frac{1}{\beta}(d_{i}-1)\mathcal{S}_{i}\big{)}=-\frac{1}{\beta}\frac{(d_{i}-1)}{q_{i}(1-q_{i})}$)
that we absorb from $\textbf{H}_{B}$ into the edge-specific Hessian
$\textbf{H}_{B}^{(i,j)}$. Analogously, $s_{ji}$ does the same with the second-
order local entropy contribution related to node $j$ (i.e., with
$-\frac{1}{\beta}\frac{(d_{j}-1)}{q_{j}(1-q_{j})}$)). In other words, for each
node $i$, we distribute its second-order local entropy contributions from
$\textbf{H}_{B}$ over all edge-specific Hessians $\textbf{H}_{B}^{(i,j)}$ with
$j\in\mathcal{N}(i)$. To make this a reasonable approach, where each $s_{ij}$
represents a certain proportion of the corresponding second-order local
entropy contribution, we assume that all parameters satisfy

$\displaystyle 0\leq s_{ij},s_{ji}\leq 1.$ (67)

Then, to guarantee that the sum decomposition (65) holds, we introduce linear
constraints

$\displaystyle\sum\limits_{j\in\mathcal{N}(i)}s_{ij}=1$ (68)

for all nodes. By this, the second-order local entropy contributions are
exactly distributed over the edge-specific Hessians. The other terms in the
definition (66) of $\textbf{H}_{B}^{(i,j)}$ (i.e., the main diagonal entries
$\frac{1}{\beta}\big{(}\frac{q_{j}(1-q_{j})}{T_{ij}}\big{)}$ resp.
$\frac{1}{\beta}\big{(}\frac{q_{i}(1-q_{i})}{T_{ij}}\big{)}$, and the off-
diagonal entries
$\frac{1}{\beta}\,\Big{(}\frac{q_{i}q_{j}-\xi_{ij}^{\ast}}{T_{ij}}\Big{)}$)
stem from the second-order pairwise entropy contributions
$\frac{\partial^{2}}{\partial
q_{i}^{2}}\big{(}-\frac{1}{\beta}\big{(}\mathcal{S}_{ij}\big{)}\big{)}$ in
$\textbf{H}_{B}$.

One favorable property of this specific class of sum decompositions is the
sparseness of its involved matrices. From (66), we observe that each edge-
specific Hessian $\textbf{H}_{B}^{(i,j)}$ has precisely four non-zero entries
at the index positions $[i,i],[i,j],[j,i]$, and $[j,j]$, with identical
entries at $[i,j]$ and $[j,i]$. More generally, we deal with a sum of
quadratic matrices of the form

$\displaystyle\textbf{M}^{(i,j)}\coloneqq\begin{pmatrix}0&\cdots&&0&\cdots&&0\\\
\vdots&\ddots&&\vdots&&&\vdots\\\ &&a^{(i,j)}&\cdots&c^{(i,j)}&&\\\
0&\cdots&\vdots&&\vdots&\cdots&0\\\
\vdots&&c^{(i,j)}&\cdots&b^{(i,j)}&&\vdots\\\ &&&\vdots&&\ddots&\\\
0&\cdots&&0&\cdots&&0\end{pmatrix}$ (69)

where at most one matrix in the sum can have non-zero entries at a specific
pair of off-diagonal positions $[i,j]$ and $[j,i]$ (which is precisely the
case if $(i,j)$ is an edge in the graph). Here is another subtlety to be
addressed: Strictly speaking, this kind of matrix can never be positive
definite – but only positive semidefinite – as one can always find vectors
$\bm{y}\neq\bm{0}$ such that $\bm{y}^{T}\textbf{M}^{(i,j)}\bm{y}=\bm{0}$.
However, one can show that a sum of positive semidefinite matrices of the form
(69) is positive definite, if only the $(2\times 2)$ \- submatrices

$\displaystyle\textbf{M}_{2\times
2}^{(i,j)}\coloneqq\begin{pmatrix}a^{(i,j)}&c^{(i,j)}\\\
c^{(i,j)}&b^{(i,j)}\end{pmatrix}$ (70)

are positive definite for each individual matrix $\textbf{M}^{(i,j)}$ in the
sum 222222To prove this, write
$\bm{y}^{T}\big{(}\sum\textbf{M}^{(i,j)}\big{)}\bm{y}$ equivalently as a sum
of quadratic forms
$\sum\begin{pmatrix}y_{i}&y_{j}\end{pmatrix}\textbf{M}_{2\times
2}^{(i,j)}\begin{pmatrix}y_{i}\\\ y_{j}\end{pmatrix}$ (by exploiting the
structure (69) of $\textbf{M}^{(i,j)}$) and use the positive definiteness of
all submatrices $\textbf{M}_{2\times 2}^{(i,j)}$.. This facilitates our
analysis, as we can focus on studying the definiteness properties of $(2\times
2)$ \- submatrices of the edge-specific Hessians $\textbf{H}_{B}^{(i,j)}$
defined as

$\displaystyle\textbf{H}_{2\times
2}^{(i,j)}\coloneqq\frac{1}{\beta}\begin{pmatrix}-s_{ij}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{T_{ij}}&\frac{q_{i}q_{j}-\xi_{ij}^{\ast}}{T_{ij}}\\\
\frac{q_{i}q_{j}-\xi_{ij}^{\ast}}{T_{ij}}&-s_{ji}\,\frac{(d_{j}-1)}{q_{j}(1-q_{j})}+\frac{q_{i}(1-q_{i})}{T_{ij}}\end{pmatrix}.$
(71)

So far, we have explained our approach to deriving a sufficient criterion for
convexity of the Bethe free energy $\mathcal{F}_{B}$ on the Bethe box
$\mathbb{B}$ by writing its Hessian $\textbf{H}_{B}$ (Lemma 4, (b)) as a sum
of matrices and proving positive definiteness of the individual matrices in
that sum. Moreover, we have introduced a class of edge-specific sum
decompositions (65) and shown that an edge-specific Hessian (66) is positive
semidefinite, if and only if a $(2\times 2$) - submatrix $\textbf{H}_{2\times
2}^{(i,j)}$ of $\textbf{H}_{B}^{(i,j)}$, defined by (71), is positive
definite. We have also argued that positive semidefiniteness of all edge-
specific Hessians $\textbf{H}_{B}^{(i,j)}$ is sufficient for positive
definiteness of the Bethe Hessian $\textbf{H}_{B}$. In the remainder of this
section, we analyze when the submatrices $\textbf{H}_{2\times 2}^{(i,j)}$ are
positive definite 232323Note that $\textbf{H}_{B}$ and the matrices in its sum
decomposition are functions defined over the Bethe box, i.e., depend on
$\bm{q}\in\mathbb{B}$. This implies that the edge-specific Hessians and their
submatrices must be positive (semi-)definite for all pairs $(q_{i},q_{j})$ in
$(0,1)\times(0,1)$ (as $\textbf{H}_{B}^{(i,j)}$ depends only on the associated
variables $q_{i}$ and $q_{j}$) to guarantee convexity of $\mathcal{F}_{B}$.;
i.e., in dependence of its involved parameters (which are the independent
variables $q_{i},q_{j}$; the model-specific parameters $d_{i}$, $d_{j}$,
$J_{ij}$; the constrained parameters $s_{ij},s_{ji}$; and the inverse
temperature $\beta$).

According to Sylvester’s criterion (Sec. 2.1), a quadratic real matrix is
positive definite if and only if all its leading principle minors are
positive. In the case of a $(2\times 2)$ \- matrix, these are the upper left
element of the matrix and the determinant of the entire matrix. This implies
that a submatrix $\textbf{H}_{2\times 2}^{(i,j)}$ of $\textbf{H}_{B}^{(i,j)}$
is positive definite in an arbitrary point $(q_{i},q_{j})$ in
$(0,1)\times(0,1)$, if and only if the inequalities

$\displaystyle-
s_{ij}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{T_{ij}}>0\quad\text{and}$
(72)
$\displaystyle\begin{split}&\Big{(}-s_{ij}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{T_{ij}}\Big{)}\cdot\Big{(}-s_{ji}\,\frac{(d_{j}-1)}{q_{j}(1-q_{j})}+\frac{q_{i}(1-q_{i})}{T_{ij}}\Big{)}\\\
&-\Big{(}\frac{q_{i}q_{j}-\xi_{ij}^{\ast}}{T_{ij}}\Big{)}^{2}>0\end{split}$
(73)

are satisfied (for some prespecified parameters $s_{ij},s_{ji}$). The left-
hand side of (73) is the determinant of $\textbf{H}_{2\times 2}^{(i,j)}$,
where we have omitted the positive scale factor $\frac{1}{\beta^{2}}$ as it
does not influence the sign of the determinant (still, we have implicit
dependence on the inverse temperature $\beta$ via $\xi_{ij}$ and $T_{ij}$).
These formulas are analytically complex, as $\xi_{ij}$ and $T_{ij}$ are
themselves functions of $q_{i},q_{j}$ (defined in (24) and (28)). Also, both
inequalites must be satisfied for all $(q_{i},q_{j})$ in $(0,1)\times(0,1)$
(or, equivalently, for the infimum over all such pairs) which complicates
their verification. To simplify the analysis, we further specify our
definition of edge-specific Hessians. Recall that there is some freedom in
selecting the constrained parameters $s_{ij},s_{ji}$, as long as they satisfy
the constraints (67) and (68). A particularly natural choice arises if we
distribute the second-order local entropy contributions related to node $i$
uniformly over the edge-speficic Hessians $\textbf{H}_{B}^{(i,j)}$ (for $j$ in
$\mathcal{N}(i)$); i.e., if we take $s_{ij}=\frac{1}{d_{i}}$ and
$s_{ji}=\frac{1}{d_{j}}$. Clearly, this is a valid choice since
$0<\frac{1}{d_{i}}<1$ and
$\sum\limits_{j\in\mathcal{N}(i)}\frac{1}{d_{i}}=d_{i}\cdot\frac{1}{d_{i}}=1$.
Note that other choices of $s_{ij},s_{ji}$ might induce more powerful criteria
for convexity 242424One possible alternative would be to distribute the
second-order local entropy contribution according to the strength of couplings
related to the edges between neighboring nodes; i.e., take
$s_{ij}=\frac{|J_{ij}|}{\sum\limits_{k\in\mathcal{N}(i)}|J_{ik}|}$ . Also, it
remains an open problem if we can ’optimize’ the choice of $s_{ij},s_{ji}$ to
minimize the gap between the convexity induced by a certain criterion and the
true convexity. (as this might be the case for different kind of sum
decompositions than the edge-specific one); however, our choice entails the
positive aspect that inequality (72) is always satisfied:

$\displaystyle-
s_{ij}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{T_{ij}}$
$\displaystyle=\,$
$\displaystyle-\frac{1}{d_{i}}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{T_{ij}}$
$\displaystyle\stackrel{{\scriptstyle\eqref{eq:Tij_optimal}}}{{=}}\,$
$\displaystyle-\frac{1}{d_{i}}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{q_{i}q_{j}(1-q_{i})(1-q_{j})-(\xi_{ij}^{\ast}-q_{i}q_{j})^{2}}$
$\displaystyle>\,$
$\displaystyle-\frac{1}{d_{i}}\,\frac{(d_{i}-1)}{q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{q_{i}q_{j}(1-q_{i})(1-q_{j})}$
$\displaystyle=\,$
$\displaystyle\frac{-(d_{i}-1)+d_{i}}{d_{i}q_{i}(1-q_{i})}\,=\,\frac{1}{d_{i}q_{i}(1-q_{i})}\,>\,0$

Using this ’uniform’ choice of $s_{ij},s_{ji}$, we can thus focus on
inequality (73) whose left-hand side is the determinant of
$\textbf{H}_{2\times 2}^{(i,j)}$, defined as

$\displaystyle\begin{split}&\det\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})\\\
&\coloneqq\Big{(}\frac{-(d_{i}-1)}{d_{i}q_{i}(1-q_{i})}+\frac{q_{j}(1-q_{j})}{T_{ij}}\Big{)}\cdot\Big{(}\frac{-(d_{j}-1)}{d_{j}q_{j}(1-q_{j})}+\frac{q_{i}(1-q_{i})}{T_{ij}}\Big{)}-\Big{(}\frac{q_{i}q_{j}-\xi_{ij}^{\ast}}{T_{ij}}\Big{)}^{2}.\end{split}$
(74)

Before we proceed, we simplify the expression in (74) and work with a more
compact formula of $\det\textbf{H}_{2\times 2}^{(i,j)}$. Namely, (74) is the
same as

$\displaystyle\frac{d_{i}-1}{d_{i}}\cdot\frac{d_{j}-1}{d_{j}}\cdot\frac{1}{q_{i}q_{j}(1-q_{i})(1-q_{j})}-\frac{1}{T_{ij}}\cdot\Big{(}\frac{d_{i}-1}{d_{i}}+\frac{d_{j}-1}{d_{j}}\Big{)}$
$\displaystyle+\,$
$\displaystyle\frac{\overbrace{q_{i}q_{j}(1-q_{i})(1-q_{j})-(q_{i}q_{j}-\xi_{ij}^{\ast})^{2}}^{=T_{ij}}}{T_{ij}^{2}}$
$\displaystyle=\,$
$\displaystyle\frac{d_{i}-1}{d_{i}}\cdot\frac{d_{j}-1}{d_{j}}\cdot\frac{1}{q_{i}q_{j}(1-q_{i})(1-q_{j})}+\frac{1}{T_{ij}}\Big{(}1-\frac{d_{i}-1}{d_{i}}-\frac{d_{j}-1}{d_{j}}\Big{)}.$
(75)

Let us summarize the essentials of our present analysis:

###### Proposition 16

Let $\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})$ be the submatrices of the
edge-specific Hessians $\textbf{H}_{B}^{(i,j)}(q_{i},q_{j})$ defined in (71)
resp. (66) together with a set of parameters $\\{s_{ij},s_{ji}\\}$ for each
edge $(i,j)$ satisfying the constraints (67) and (68). Then the Bethe free
energy $\mathcal{F}_{B}$ is convex on the Bethe box $\mathbb{B}$, if all
submatrices $\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})$ are positive
definite for all $(q_{i},q_{j})$ in $(0,1)\times(0,1)$, i.e., if the
inequalities (72) and (73) are satisfied. If we set $s_{ij}=\frac{1}{d_{i}}$,
then $\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})$ is positive definite if and
only if the determinant of $\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})$
defined in (75) is positive.

To use these qualitative statements for proving the convexity of
$\mathcal{F}_{B}$, we need to characterize the behavior of
$\det\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})$ (as a function of $q_{i}$
and $q_{j}$) in dependence of its parameters. In particular, we must determine
whether $\det\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})$ is positive for all
points $(q_{i},q_{j})$ in $(0,1)\times(0,1)$, or if there exists a point where
it is non-positive (in which case positive definiteness of
$\textbf{H}_{2\times 2}^{(i,j)}$ would be violated). In principle, we could
once again formulate an optimization problem in which we try to compute the
infimum of $\det\textbf{H}_{2\times 2}^{(i,j)}$ and determine whether it is
positive or not (analogously to problem (32) in Sec. 3.2); however, we follow
a slightly different path in this analysis. Let us consider an edge $(i,j)$ in
the graph with fixed degrees $d_{i}$, $d_{j}$ of its end nodes and with an
associated coupling $J_{ij}$, and let us assume that we vary the inverse
temperature $\beta$. First of all, we observe that $\det\textbf{H}_{2\times
2}^{(i,j)}$ is positive on its domain if $\beta$ approaches zero (or in other
words, in the infinite temperature limit $T\to\infty$):

###### Lemma 17

In the infinite temperature limit, we have for all $(q_{i},q_{j})$ in
$(0,1)\times(0,1):$

$\displaystyle\lim\limits_{\beta\to\,0}\big{[}\det\textbf{H}_{2\times
2}^{(i,j)}(q_{i},q_{j})\big{]}=\frac{1}{d_{i}d_{j}}\cdot\frac{1}{q_{i}q_{j}(1-q_{i})(1-q_{j})}>0$
(76)

Lemma 17 characterizes the behavior of $\det\textbf{H}_{2\times 2}^{(i,j)}$ in
the infinite temperature limit. It also states that the determinant is
positive on its domain, if the inverse temperature $\beta$ is sufficiently
small (because $\det\textbf{H}_{2\times 2}^{(i,j)}$ is a continuous function
of $q_{i}$ and $q_{j}$). This proves the existence 252525We remark that
$\beta^{\ast}$ might take the value $\infty$, i.e., $\mathcal{F}_{B}$ might be
convex for any temperature. Also, if we additionally require that
$\beta^{\ast}$ is the greatest positive number (including $\infty$) satisfying
the described property, then it is even unique. of some critical threshold
$\beta^{\ast}>0$ such that, for all values $0<\beta<\beta^{\ast}$, Proposition
16 is satisfied 262626More precisely, Lemma 17 is only an edge-specific
result, where we have considered the submatrix $\textbf{H}_{2\times
2}^{(i,j)}$ associated to some arbitrary edge $(i,j)$ in the graph. However,
the statement in the text follows immediately, if we have one such ’edge-
specific’ threshold $(\beta^{\ast})^{(i,j)}$ for each edge $(i,j)$ and then
identify their minimum as the critical threshold $\beta^{\ast}$., and
therefore justifies our approach in this section. Moverover, Proposition 16
and Lemma 17 prove a fundamental result:

###### Corollary 18

For a sufficiently high temperature, the Bethe free energy $\mathcal{F}_{B}$
is convex on the Bethe box $\mathbb{B}$.

Regarding our previous analysis, this statement is not surprising; however, to
the best of our knowledge, this is the first time that this result has been
proven in the literature.

The next question is how $\det\textbf{H}_{2\times 2}^{(i,j)}$ behaves, if we
make an infinitesimal change in the (inverse) temperature. Intuitively, we
expect the determinant to decrease if we reduce the temperature, as the Bethe
approximation often deteriorates with stronger interactions in the model
(thus, the conditions in Proposition 16 would be expected to be violated at
some point, Meshi et al. [2009], Weller et al. [2014], Knoll and Pernkopf
[2020]). This is confirmed by the following Lemma, in which we prove an even
stronger result, namely that $\det\textbf{H}_{2\times 2}^{(i,j)}$ is a
monotonically decreasing function of $\beta$:

###### Lemma 19

For all $(q_{i},q_{j})$ in $(0,1)\times(0,1)$, the determinant of
$\textbf{H}_{2\times 2}^{(i,j)}(q_{i},q_{j})$ is monotonically decreasing, as
$\beta$ increases, i.e.,

$\displaystyle\frac{\partial}{\partial\beta}\big{[}\det\textbf{H}_{2\times
2}^{(i,j)}(q_{i},q_{j})\big{]}<0.$ (77)

Lemma 19 makes an important statement as it says the following: If the
determinant of $\textbf{H}_{2\times 2}^{(i,j)}$ is negative for some value of
the inverse temperature $\beta$ (by which positive definiteness of
$\textbf{H}_{2\times 2}^{(i,j)}$ and thus the condition in Proposition 16
would be violated), then it is also negative for all larger values of $\beta$.
On the other hand, we know by Lemma 17 that $\det\textbf{H}_{2\times
2}^{(i,j)}$ is always positive if $\beta$ is sufficiently small. To complete
the image, we will next show that there exists a real (and unique) value
$\beta^{\ast}$, for which $\det\textbf{H}_{2\times 2}^{(i,j)}$ equals zero in
some point $(q_{i},q_{j})$ in $(0,1)\times(0,1)$, and which is consequently
the smallest $\beta$ for which the determinant is not positive. Moreover, we
show that the point where this scenario occurs is precisely the ’center’ of
the domain, i.e., the point $(0.5,0.5)$. Figure 27 illustrates this behavior
of the determinant. Regarding our analysis, this implies that we can
explicitly compute a critical value $\beta^{\ast}$ such that Proposition 16 is
either satisfied or violated.

Figure 8: Behavior of $\det\textbf{H}_{2\times 2}^{(i,j)}$ on its domain
$(0,1)\times(0,1)$ for $d_{i}=d_{j}=3$ and $J_{ij}=0.5$ (corresponding to a
ferromagnetic edge whose end nodes have both three neighbors), and for varying
values of the inverse temperature $\beta$.
Top row: 3D plots of $\det\textbf{H}_{2\times 2}^{(i,j)}$ for increasing
values of $\beta$. The green planes represent the $q_{i}q_{j}$ \- plane (i.e,
the zero level of the determinant). If the graph of $\det\textbf{H}_{2\times
2}^{(i,j)}$ lies above this plane, then $\textbf{H}_{2\times 2}^{(i,j)}$ is
positive definite for the associated value of $\beta$; if the graph cuts (or
touches) the plane, it is not. According to Theorem 23,
$\det\textbf{H}_{2\times 2}^{(i,j)}$ is positive on its domain if and only if
it is positive in the center point $(0.5,0.5)$. This is, e.g., the case for
$\beta=0.5$, as is indicated by the green circle in the left-hand figure. In
the middle figure, the graph of $\det\textbf{H}_{2\times 2}^{(i,j)}$ touches
the green plane in the center point. This occurs precisely for the critical
value $\beta^{\ast}$, which is the smallest $\beta$ such that
$\det\textbf{H}_{2\times 2}^{(i,j)}$ is not longer positive on its domain.
This is because $\det\textbf{H}_{2\times 2}^{(i,j)}$ takes the value zero in
the center point (Lemma 20), as is indicated by a red circle. Since
$\det\textbf{H}_{2\times 2}^{(i,j)}$ is monotonically decreasing as $\beta$
increases (Lemma 19), it is negative in $(0.5,0.5)$ for all larger values of
$\beta$. In the right-hand figure ($\beta=2$), the graph of
$\det\textbf{H}_{2\times 2}^{(i,j)}$ cuts the green plane in multiple points.
Bottom row: The same graphs as in the top row, but from a side perspective at
the level of the $q_{i}q_{j}$ \- plane (with the $q_{i}$ \- axis located in
front of the viewpoint and the $q_{j}$ \- axis hidden behind). The area around
the center point has been magnified, so that the described behavior is better
visible272727Due to numerical reasons, an exact visualization of the situation
for $\beta=\beta^{\ast}$ was not possible..

Before we address the global properties of the determinant, we focus on the
center point $(0.5,0.5)$. The following result provides the critical value
$\beta^{\ast}$ of the inverse temperature, at which the determinant is zero in
that point:

###### Lemma 20

Let $d_{i},d_{j}>2$ and let

$\displaystyle\beta^{\ast}\coloneqq\frac{1}{2|J_{ij}|}\operatorname{arccosh}\big{(}1+\frac{2}{d_{i}d_{j}-d_{i}-d_{j}}\big{)}$
(78)

be the ’critical’ value of the inverse temperature. Then

$\displaystyle\det\textbf{H}_{2\times
2}^{(i,j)}(0.5,0.5)\begin{cases}&>0\quad\text{for}\quad\beta<\beta^{\ast},\\\
&=0\quad\text{if}\quad\,\,\,\,\beta=\beta^{\ast},\quad\text{and}\\\
&<0\quad\text{for}\quad\beta>\beta^{\ast}.\end{cases}$ (79)

Lemma 20 provides the critical threshold $\beta^{\ast}$ of the inverse
temperature at which $\det\textbf{H}_{2\times 2}^{(i,j)}$, evaluated in the
center point $(0.5,0.5)$, changes its sign from positive to negative (i.e., at
which $\det\textbf{H}_{2\times 2}^{(i,j)}(0.5,0.5)$ becomes indefinite). This
may not seem especially important as there might exist other points
$(q_{i},q_{j})$ in the domain, in which the determinant is already negative
for smaller values of $\beta$ than $\beta^{\ast}$ (while we are searching for
the greatest value of $\beta$ at which the determinant is positive on its
entire domain). However, as already indicated ahead of Lemma 20, the center
point is precisely the point where $\det\textbf{H}_{2\times 2}^{(i,j)}$ first
changes its sign if the inverse temperature steadily increases (when starting
from close to zero).

###### Theorem 21

Let $d_{i},d_{j}>2$, let $\beta^{\ast}$ be the critical threshold of the
inverse temperature derived in Lemma 20, and let $\det\textbf{H}_{2\times
2}^{(i,j)}$ (considered as a function of $(q_{i},q_{j})$) implicitly depend on
$\beta^{\ast}$. Then $\det\textbf{H}_{2\times 2}^{(i,j)}$ has a unique root in
the center point $(0.5,0.5)$ and is positive for all other points
$(q_{i},q_{j})$ in $(0,1)\times(0,1)$.

Theorem 21 is a strong result that fully characterizes the global behavior of
$\det\textbf{H}_{2\times 2}^{(i,j)}$ on its domain. It says that, if the
inverse temperature of the model takes the critical value $\beta^{\ast}$
(derived in Lemma 20), there exists precisely one point in the domain – the
center point – where the determinant is non-positive (while it is positive for
all other points). As a consequence, the conditions for convexity of the Bethe
free energy from Proposition 16 are not satisfied if $\beta=\beta^{\ast}$. If
we combine this knowledge with the insights gained in Lemma 19 (the
determinant monotonically decreases as $\beta$ increases), we can conclude as
follows: On the one hand, $\det\textbf{H}_{2\times 2}^{(i,j)}$ must be
positive on the entire domain whenever the inverse temperature $\beta$ is
below that critical threshold $\beta^{\ast}$; on the other hand,
$\det\textbf{H}_{2\times 2}^{(i,j)}$ is negative in $(0.5,0.5)$ whenever the
inverse temperature $\beta$ is beyond that critical threshold $\beta^{\ast}$.
Hence, the complete information about the global definiteness properties of
$\textbf{H}_{2\times 2}^{(i,j)}$ lies already in the center point: If
$\textbf{H}_{2\times 2}^{(i,j)}$ is positive definite in $(0.5,0.5)$, then it
is positive definite on its entire domain. This induces a particularly simple
criterion for checking the positive definiteness of $\textbf{H}_{2\times
2}^{(i,j)}$ on its domain, as one simply needs to compute the critical value
of the inverse temperature $\beta^{\ast}$ derived in Lemma 20. If
$\beta<\beta^{\ast}$, then $\textbf{H}_{2\times 2}^{(i,j)}$ is positive
definite; otherwise it is indefinite (because at least the main diagonal
entries are always positive). Lemma 19, Lemma 20, and Theorem 21 therefore
imply the following corollary:

###### Corollary 22

If $\det\textbf{H}_{2\times 2}^{(i,j)}$ is positive in the center point
$(0.5,0.5)$, then it is positive for all $(q_{i},q_{j})$ in
$(0,1)\times(0,1)$. This is precisely the case, if the inverse temperature of
the model lies below the critical threshold $\beta^{\ast}$. Equivalently, if
$\textbf{H}_{2\times 2}^{(i,j)}$ is positive definite in the center point
$(0.5,0.5)$, then it is positive definite 282828More precisely, the second
statement is true if we have defined the constrained parameters
$s_{ij},s_{ji}$ in the definition of the edge-specfic Hessians (66) as
$s_{ij}=\frac{1}{d_{i}},s_{ji}=\frac{1}{d_{j}}$; not necessarily, however, for
all parameters $s_{ij},s_{ji}$ satisfying the constraints (67), (68). for all
$(q_{i},q_{j})$ in $(0,1)\times(0,1)$.

Up to this point, we have analyzed the properties of a single edge-specific
Hessian $\textbf{H}_{B}^{(i,j)}$ in the sum decomposition (65) of the Bethe
Hessian $\textbf{H}_{B}$. However, there might of course be values of the
critical temperature for which certain matrices in the sum are positive
definite, while others are not. Our sufficient criterion for positive
definiteness of $\textbf{H}_{B}$ (and thus convexity of $\mathcal{F}_{B}$)
requires positive definiteness of each individual matrix in that sum. We
finally arrive at the central result of this section, which is a reformulation
of Proposition 16, based on the insights that we have gained from our
analysis:

###### Theorem 23

Let $\beta$ be the inverse model temperature. The Bethe free energy
$\mathcal{F}_{B}$ is convex on the Bethe box $\mathbb{B}$ if, for all edges
$(i,j)$ in the graph with degrees $d_{i},d_{j}>2$ of their end nodes, the
following inequality is satisfied:

$\displaystyle\beta<\frac{1}{2|J_{ij}|}\operatorname{arccosh}\Big{(}1+\frac{2}{d_{i}d_{j}-d_{i}-d_{j}}\Big{)}$
(80)

In Sec. 3.2 and Sec. 3.3 we have derived conditions for the convexity of
$\mathcal{F}_{B}$ on $\mathbb{B}$, and thus provided answers to the second of
our three questions (i) - (iii) formulated in Sec. 3.1 that were the starting
point of our analysis. In the remaining part of this work, we address question
(iii): Can we explain the reliability of the Bethe approximation through its
convexity on the Bethe box $\mathbb{B}$? Or in other words, is convexity of
$\mathcal{F}_{B}$ on $\mathbb{B}$ a property that is necessary and/or
sufficient for guaranteeing a high quality of the marginal and partition
function estimates provided by the Bethe approximation (according to the
approximations made in (15) and (16), Sec. 2.3)? In Sec. 3.4 we specify this
query and discuss a few of its theoretical aspects. In Sec. 4 we present our
experimental results.

###  3.4 Theoretical Discussion

The purpose of this section is to recap a few characteristics and differences
of our main results from Sec. 3.2 (Theorem 15) and Sec. 3.3 (Theorem 23).
Furthermore, we put them into the context of similar results from the
literature (Sec. 2.4), by comparing them on a perfectly symmetric model (a
complete graph with homogenous potentials) whose simplicity allows for an
analytical specfication of its phase transition boundaries 292929More
precisely, this characterization is valid for the infinite version of a
complete graph, which is the Cayley tree. However, we argue that at least for
the case of vanishing fields (i.e., $\theta=0$) these boundaries are sharp for
the finite versions as well. [Georgii, 1988, Taga and Mase, 2006b]. Then, to
discuss theoretical aspects of question (iii), we analyze the behavior of the
Bethe free energy in dependence of the inverse temperature, by distinguishing
between its three different ’stages’: convexity, non-convexity but uniqueness
of a Bethe minimum, and the coexistence of multiple Bethe minima. This
specification refines the conceptual view on the Bethe approximation and
motivates our experimental analysis in Sec. 4.

Let us once again summarize the main results from the two previous sections,
both being sufficient conditions for the convexity of the Bethe free energy
$\mathcal{F}_{B}$ on the Bethe box $\mathbb{B}$. According to Theorem 15
(derived from diagonal dominance of the Bethe Hessian), $\mathcal{F}_{B}$ is
convex on $\mathbb{B}$ if, for all nodes $i$ in the graph, the one-dimensional
polynomial

$\displaystyle\Psi_{i}(q_{i})=-(d_{i}-1)\prod\limits_{j\in\mathcal{N}(i)}(1+\alpha_{ij}q_{i})+\sum\limits_{j\in\mathcal{N}(i)}\Big{(}(1+\alpha_{ij}q_{i}^{2})\prod\limits_{k\in\mathcal{N}(i)\setminus
j}(1+\alpha_{ik}q_{i})\Big{)}$

is strictly positive on the interval $(0,0.5]$. To apply this result in
practice, we need to verify for all polynomials $\Psi_{i}(q_{i})$ (for
$i\in\mathbf{X}$) if none of them has a root in the interval $(0,0.5]$. This
can be done efficiently by computational methods from numerical optimization
303030If an exact statement is desired, one can apply Sturm’s theorem which
calculates the total number of real roots of any one-dimensional polynomial in
an arbitrary interval [Thomas, 1941].. The degree of $\Psi_{i}(q_{i})$ is
$d_{i}+1$ and thus increases with the number of neighbors of node $i$. This
can make a numerical search for potential roots of $\Psi_{i}(q_{i})$ more
challenging (still, it remains a polynomial in one variable). The overall
computational effort increases at least linearly with the number of nodes in
the graph (since positivity must be verified for $|\mathbf{X}|$ polynomials).

Other than the above result, Theorem 23 (derived from a sum decomposition of
the Bethe Hessian) provides a condition in closed-form. It says that
$\mathcal{F}_{B}$ is convex on $\mathbb{B}$ if, for all edges $(i,j)$ in
$\mathbf{E}$ whose end nodes $i$ and $j$ have at least three neighbors, the
inverse temperature of the model satisfied the inequality

$\displaystyle\beta<\frac{1}{2|J_{ij}|}\operatorname{arccosh}\Big{(}1+\frac{2}{d_{i}d_{j}-d_{i}-d_{j}}\Big{)}.$

This condition is particularly easy to verify and additionally provides an
explicit estimate of the temperature at which the Bethe free energy might
become non-convex (by taking the greatest $\beta$ that simultanously satisfies
all inequalities – i.e., for all edges $(i,j)$ – of the above form). Note that
Theorem 15 induces such an estimate only in implicit form, given as the
greatest $\beta$ such that none of the polynomials $\Psi_{i}(q_{i})$ has a
root in $(0,0.5]$. To compute it, one can steadily increase $\beta$ and stop
when the first polynomial $\Psi_{i}(q_{i})$ has a root in $(0,0.5]$. The
computational effort for checking Theorem 23 increases linearly with the
number of edges in the graph (whereby for each edge only an elementary
calculation needs to be performed).

In summary, Theorem 23 has computational advantages over Theorem 15 and
appears to be better interpretable. However, this does not imply that it is
more powerful in the sense of inducing more accurate estimates of the critical
temperature at which $\mathcal{F}_{B}$ changes its stage from convex to non-
convex.

In the sequel, we put our results in a broader context by comparing them to
similar results from the literature. An overview of all important statements
and their connections has been given in Fig. 1 (Sec. 2.4). First, it seems
convenient to categorize them according to their purposes and properties.
Subsequently, we discuss them on complete graphs with homogenous potentials.
For further details, we refer to Mooij and Kappen [2007].

  * •

In infinite systems, the theory of Gibbs measures interprets a phase
transition as the existence of a unique Gibbs measure. The result of Dobrushin
[1968] is sufficient for uniqueness of a Gibbs measure and takes the influence
of the fields $\theta_{i}$ into account; however, its evaluation time
increases exponentially with the number of nodes in the graph. The condition
of Simon [1979] is sufficient for Dobrushin’s condition and simpler to
evaluate, but it is only based on the couplings $J_{ij}$ and thus weaker for
models that are strongly influenced by the fields $\theta_{i}$.

  * •

The relationship between Gibbs measure theory and approximate inference in
(finite) probabilistic graphical models has been established by Tatikonda and
Jordan [2002]. Namely, loopy belief propagation converges to a unique fixed
point if there exists a unique Gibbs measure on the computation tree (which is
an infinite unwrapping of the graph; in particular, Gibbs measures are well-
defined on such a tree). To the best of our knowledge, nothing is known about
the reverse direction of this statement.

  * •

The results of Ihler et al. [2005], Mooij and Kappen [2007], Watanabe and
Fukumizu [2009] are sufficient for the convergence of LBP to a unique fixed
point, and thus for the uniqueness of a Bethe minimum. Their evaluation is
based on closed-form expressions or eigenvalue problems, and can be done
efficiently. Moreover, Mooij’s criterion also takes the influence of fields
$\theta_{i}$ into account. Experimental results indicate that Mooij’s
criterion is the most powerful condition 313131By most powerful we mean that
it correctly predicts the uniqueness of an LBP fixed point for the largest
variety of models, when compared to other conditions. This usually implies
more accurate estimates of the temperature, at which multiple LBP fixed points
(or Bethe minima) arise. for predicting uniqueness of an LBP fixed point (or
Bethe minimum). Another sufficient condition for the uniqueness of an LBP
fixed point (but not convergence to it) is that of Heskes [2004]. It appears
to be less powerful than other criteria, but is of theoretical importance. Its
evaluation requires the solution of a linear program and only takes the
couplings $J_{ij}$ into account.

  * •

The results in this work (Theorem 15 and Theorem 23) are sufficient for the
convexity of the Bethe free energy $\mathcal{F}_{B}$ on a submanifold
$\mathbb{B}$ (the Bethe box) of its actual domain $\mathbb{L}$ (the local
polytope). As Proposition 1 shows, $\mathcal{F}_{B}$ can never be convex on
the entire local polytope, unless the graph contains at most one cycle.
However, convexity of $\mathcal{F}_{B}$ on $\mathbb{B}$ implies the uniqueness
of a Bethe minimum (Cor. 3, Sec. 3.1). Our results are based on the couplings
$J_{ij}$ only, because convexity of $\mathcal{F}_{B}$ is a property that is
not influenced by the fields $\theta_{i}$ (Cor. 5). Both conditions can be
evaluated efficiently.

These results fulfill different purposes, being sufficient for convergence of
LBP to a unique fixed point, uniqueness of a Gibbs measure on the computation
tree, or convexity of the Bethe free energy. Still, all of them estimate a
critical temperature $\beta^{\ast}$ – either implicitly (e.g., Theorem 15) or
explicitly (e.g., Theorem 23) – at which the Bethe free energy changes its
stage (from convex to non-convex, or from uniqueness of a minimum to multiple
minima).

Ultimately, we aim to better understand the impact of the current stage of the
Bethe free energy on its approximation accuracy (i.e., with respect to the
marginals and partition function). For many models, the quality of the Bethe
approximation significantly degrades at some critical temperature, and thus
loses its reliability. This effect can be interpreted as a finite version of a
’phase transition’ in the model 323232These are not phase transitions in the
sense of infinite models, because the partition function in finite models is
always a continuous function of the temperature (Sec. 2.4).. In our analysis,
we aim to understand if this reliabilty loss can be attributed to a specific
stage change of the Bethe free energy.

As we do not have an oracle to tell us in which stage $\mathcal{F}_{B}$ is at
a certain temperature, we must rely on the results from the literature and our
work. To get a notion about their individual power, we perform an analysis on
a particularly simple model whose behavior can be exactly characterized: a
complete graph with homogenous potentials (i.e., each node is connected to
each other with a uniform node degree $d=d_{i}$, and we have $J_{ij}=J$ for
all edges and $\theta_{i}=\theta$ for all nodes). We discuss how powerful the
above discussed results are on that simple model towards predicting a state
change of the Bethe free enegry (of course, this does not imply that they are
equally powerful for more complex models).

We first consider our condition based on diagonal dominance of the Bethe
Hessian. For this model, Theorem 15 says that $\mathcal{F}_{B}$ is convex on
$\mathbb{B}$ if the polynomial

$\displaystyle\Psi_{i}(q_{i})=-(d-1)(1+\alpha_{ij}q_{i})^{d}+d(1+\alpha_{ij}q_{i}^{2})(1+\alpha_{ij}q_{i})^{d-1}$

is positive in the interval $(0,0.5]$. Note that, due to symmetries, this
polynomial is the same for all nodes. Positivity of $\Psi_{i}(q_{i})$ is
equivalent to the positivity of

$\displaystyle-(d-1)(1+\alpha_{ij}q_{i})+d(1+\alpha_{ij}q_{i}^{2}),$

which is a quadratic polynomial and has the roots
$\frac{d-1}{2d}\pm\frac{1}{2d}\sqrt[\ ]{(d-1)^{2}-\frac{4d}{\alpha_{ij}}}$. If
the discriminant is negative, there are only complex roots and
$\Psi_{i}(q_{i})$ is positive in $(0,0.5]$. This is precisely the case if
$\alpha_{ij}<\frac{4d}{(d-1)^{2}}$ or, with $\alpha_{ij}=e^{4\beta|J|}-1$ from
(24) substituted, if

$\displaystyle\beta<\frac{1}{4|J|}\log\Big{(}\frac{(d+1)^{2}}{(d-1)^{2}}\Big{)}.$
(81)

If the discriminant is zero, $\Psi_{i}(q_{i})$ has the unique root
$\frac{d-1}{2d}$ which (for $d\geq 2$) is a number in $(0,0.5]$ (i.e.,
$\Psi_{i}(q_{i})$ is not longer positive in $(0,0.5]$). If the discriminant is
positive, we have two real roots. We claim that the lower root is always in
$(0,0.5]$ and thus the diagonal dominance criterion is violated for all values
of $\beta$ that are greater or equal to the right-hand side of (81). This
follows from the observation that the value of the lower root monotonically
decreases, if $\beta$ (and thus $\alpha_{ij}$ and the square root term)
increases. In the limit $\beta\to\infty$, the lower root approaches zero and
is consequently always in $(0,\frac{d-1}{2d})$, which is a subset of
$(0,0.5]$. We conclude that the diagonal dominance criterion is satisfied, if
and only if $\beta$ satisfies the inequality (81).

Further, we consider our second result based on the sum decomposition of the
Bethe Hessian on the symmetric model. For all edges the inequalities in
Theorem 23 are identical and $\mathcal{F}_{B}$ is convex on $\mathbb{B}$ if

$\displaystyle\beta<\frac{1}{2|J|}\operatorname{arccosh}\Big{(}1+\frac{2}{d(d-2)}\Big{)}$

is satisfied (for $d\geq 3$). We can use the trigonometric identity

$\displaystyle\operatorname{arccosh}(x)=\log\big{(}x+\sqrt[\ ]{x-1}\sqrt[\
]{x+1}\big{)}$

and perform a few simplifications to show that this is equivalent to

$\displaystyle\beta<\frac{1}{2|J|}\log\Big{(}\frac{d}{d-2}\Big{)}.$ (82)

As already mentioned, the exact behavior for the symmetric model can be
characterized. More precisely, for $\theta=0$, it has been shown that Gibbs
measures become non-unique [Georgii, 1988, Taga and Mase, 2006b] and LBP stops
converging to a unique fixed point [Mooij and Kappen, 2005], if the inverse
temperature passes the critical threshold of

$\displaystyle\frac{1}{|J|}\operatorname{arctanh}\Big{(}\frac{1}{d-1}\Big{)}.$
(83)

For $d\geq 3$, this expression is the same as the critical threshold in
inequality (82). This follows directly from the trigonometric identity

$\displaystyle\operatorname{arctanh}(x)=\frac{1}{2}\big{(}\log(1+x)-\log(1-x)\big{)}.$

Thus, our result in Theorem 23 is sharp in the sense, that it correctly
predicts the value of $\beta$ at which a modified behavior occurs 333333More
precisely, Gibbs measures become non-unique and LBP stops converging to a
unique fixed point. However, $\mathcal{F}_{B}$ might still be convex and have
a unique minimum if the inverse temperature increases beyond the predicted
thresholds. This can be observed in an antiferromagnetic model ($J<0$, Fig.
10).. Other results that are sharp on the symmetric model without fields are
those of Ihler et al. [2005], Mooij and Kappen [2007], and Watanabe and
Fukumizu [2009]. Table 2 and Fig. 9 show a comparison between the critical
thresholds that are predicted by the various criteria from the literature on
the symmetric model without fields. If $\theta\neq 0$, none of these
conditions are sharp any longer (with the results of Dobrushin and Mooij
tending to more accurate predictions, as they consider the influence of the
fields).

Table 2: Estimates of a critical inverse temperature, predicted by various results from the literature on a complete graph with homogenous potentials and vanishing fields (i.e., $J_{ij}=J$ for all edges, $\theta_{i}=0$ and $d_{i}=d$ for all nodes). The conditions in the first row are sharp on this model (i.e., they correctly predict $\beta^{\ast}$ at which the model changes its behavior). Description | Reference | $\beta^{\ast}$
---|---|---
Exact characterization |  Georgii [1988],Ihler et al. [2005] | $\frac{1}{|J|}\operatorname{arctanh}\Big{(}\frac{1}{d-1}\Big{)}$
|  Mooij and Kappen [2007], |
| Theorem 23 (this work) |
Dobrushin | Dobrushin [1968] |  $d$ odd: $\frac{1}{|J|}\operatorname{arctanh}\Big{(}\frac{1}{d}\Big{)}$
|  |  $d$ even: $\frac{1}{2|J|}\operatorname{arctanh}\Big{(}\frac{2}{d}\Big{)}$
Simon | Simon [1979] | $\frac{1}{|J|\,d}$
Diagonal Dominance | Theorem 15 (this work) | $\frac{1}{4|J|}\log\Big{(}\frac{(d+1)^{2}}{(d-1)^{2}}\Big{)}$
Heskes | Heskes [2004] | $\frac{1}{4|J|}\log\Big{(}\frac{d-1}{d-2}\Big{)}$
Figure 9: Comparison of critical values of the inverse temperature predicted
by different criteria from the literature on a complete graph with homogenous
potentials and vanishing fields ($J_{ij}=J$ for all edges, $\theta_{i}=0$ and
$d_{i}=d$ for all nodes), in dependence of the degree.

After having compared various results and their power with respect to
predicting properties of a model, we briefly discuss a few aspects of these
properties themselves; specifically, we focus on the different stages of the
Bethe free energy. As convexity of $\mathcal{F}_{B}$ on the Bethe box
$\mathbb{B}$ implies uniqueness of a Bethe minimum (but not vice versa), there
is generally a certain ’temperature gap’ between these two properties. Any
condition that predicts the convexity of $\mathcal{F}_{B}$ cannot explain this
gap, as it is caused by the presence of fields $\theta_{i}$ and increases with
stronger fields. The fields can affect the location (and number) of Bethe
minima, but not the curvature of $\mathcal{F}_{B}$. If there are no fields,
convexity of $\mathcal{F}_{B}$ and uniqueness of a Bethe minimum are supposed
to be equivalent. Although we do not have a formal proof of this statement, we
would like to formulate it as a conjecture:

###### Conjecture 24

Assume we have no external fields in our model (i.e., $\theta_{i}=0$ for all
nodes $i$). Then $\mathcal{F}_{B}$ is convex on $\mathbb{B}$, if and only if
it has a unique minimum.

Figure 10 schematically illustrates the state change behavior of the Bethe
free energy in dependence of the inverse temperature and the presence (or
absence) of fields on a perfectly symmetric model. A detailed interpretation
is given in the figure caption. Interestingly, $\mathcal{F}_{B}$ preserves its
convexity (and thus its unique minimum) in the antiferromagnetic model for any
inverse temperature. We believe that this phenomenon is due to frustrations in
the model, which are cycles in the graph with an odd number of
antiferromagnetic edges [Toulouse, 1977]. While frustrations can induce
multiple Gibbs measures and prevent loopy belief propagation from converging
[Georgii, 1988, Knoll and Pernkopf, 2017, Mooij and Kappen, 2005], they may
have ’stabilizing’ effects on the convexity (and reliability) of the Bethe
free energy. Without providing further details on this problem, we would like
to formulate another conjecture:

###### Conjecture 25

Frustrations help the Bethe free energy to preserve its convexity.

Figure 10: Stage change behavior of the Bethe free energy in dependence of the
inverse temperature $\beta$ and fields $\theta$ on a complete graph with $4$
nodes and homogenous potentials (i.e., $J_{ij}=J$ for all edges and
$\theta_{i}=\theta$ for all nodes). These plots represent a projection of
$\mathcal{F}_{B}$ onto the main diagonal of $\mathbb{B}$, i.e, points
$\bm{q}=(q_{1},q_{2},q_{3},q_{4})$ with $q_{1}=\dots=q_{4}$.
First row: For $\theta=0$ and $J>0$ (a ferromagnetic model without fields),
$\mathcal{F}_{B}$ is convex and has a unique minimum in $q_{i}=0.5$ if $\beta$
is below a critical threshold of $\beta^{\ast}\approx 0.5493$. If $\beta$
increases beyond that threshold, $\mathcal{F}_{B}$ becomes non-convex and at
the same time two new minima arise (while the ’old’ minimum disappears).
Second row: For $\theta>0$ and $J>0$ (a ferromagnetic model with fields),
$\mathcal{F}_{B}$ is convex and has a unique minimum if $\beta$ is
sufficiently small. The location of the minimum is biased by the fields and
thus not in $q_{i}=0.5$. If $\beta$ increases, $\mathcal{F}_{B}$ reaches a
point where it becomes non-convex but still has a unique minimum. If $\beta$
further increases, there arises a second Bethe minimum (with the up-to-then
minimum still existing). This shows that the fields can induce a gap between
convexity of $\mathcal{F}_{B}$ and uniqueness of a Bethe minimum.
Third row: For $\theta=0$ and $J<0$ (an antiferromagnetic model without
fields), $\mathcal{F}_{B}$ preserves its convexity (and its unique minimum) at
any temperature.

The considerations made in this section are instructive and aim to explain the
behavior of the Bethe free energy in dependence of the inverse temperature. We
primarily analyzed special kind of models to understand its different
’stages’, and discussed several criteria to predict changes between them. In
the next and final step of our analysis, we put our focus on the quality of
the Bethe approximation during these stages, regarding the accuracy of the
estimates for the marginals and partition function. We will address these
topics in Sec. 4 and present our experimental results on a broad variety of
models.

##  4 Convexity of the Bethe Free Energy: Experiments

In Sec. 3.4, we have exposed how the Bethe free energy passes through
different stages if the inverse model temperature increases, and discussed
possibilities to predict in which stage it currently is. We now focus on
analyzing the reliability of the Bethe approximation, i.e., the accuracy of
the marginal and partition function estimates that it provides.

The quality of the Bethe approximation has been analyzed, e.g., by Meshi et
al. [2009], Weller et al. [2014], Knoll and Pernkopf [2020], but rather at a
qualitative level and without specifying its different stages. We extend this
analysis by a quantitative component (i.e., concrete estimates for critical
phase transition temperatures or, in other words, reliability thresholds), and
the concept of convexity on the Bethe box $\mathbb{B}$ that we have introduced
in this work. Specifically, we analyze if the Bethe approximation is always
reliable if $\mathcal{F}_{B}$ is convex on $\mathbb{B}$; if it is still
reliable if it is non-convex but has a unique minimum; or if its current stage
has no significant influence on the quality of the Bethe approximation.

In Sec. 4.1, we show how we measure the quality of the Bethe free energy
approximation and introduce our algorithm BETHE-MIN that we use for its
minimization. In Sec. 4.2, we explain our experimental setup and present all
experimental results in detail.

###  4.1 Minimization of the Bethe Free Energy

Let us recap the essential ideas of the Bethe approximation (introduced in
Sec. 2.3) to clarify the concept of the subsequent experiments. In a finite
graphical model, we intend to find approximations to the (generally non-
accessible) partition function (4), and/or marginal distributions (5), (6).
Arising from the Gibbs variational principle, the Bethe approximation provides
estimates for these unknown quantities by making the approximations (15) and
(16), i.e., by using its minimum on the set $\mathbb{L}$ of local probability
constraints, the local polytope (20). As we have shown in Sec. 3.1, we can
reduce the dimension of this domain and instead minimize the Bethe free energy
$\mathcal{F}_{B}$ on a submanifold $\mathbb{B}$ of $\mathbb{L}$, the Bethe box
(25). In other words, to approximate the partition function and singleton
marginals 343434More precisely, $q_{i}$ (the i-th component of $\bm{q}$)
approximates the first entry of the marginal $p_{i}$, i.e., the probability
that variable $X_{i}$ takes the value $+1$. The second entry (the probability
that $X_{i}$ takes the value $-1$) can be estimated with $1-q_{i}$ (Table 1).
we use

$\displaystyle-\frac{1}{\beta}\log
Z_{B}\coloneqq\min\limits_{\tilde{\bm{q}}\,\in\,\mathbb{B}}\mathcal{F}_{B}\,\,$
$\displaystyle\approx\,\,-\frac{1}{\beta}\log Z(\beta)\,,$ (84)
$\displaystyle\\{\,q_{i}:\bm{q}=\operatorname*{argmin}\limits_{\tilde{\bm{q}}\in\mathbb{B}}\mathcal{F}_{B}\,\,|\,\,i\in\mathbf{X}\,\\}\,\,$
$\displaystyle\approx\,\,\\{\,p_{i}\,\,|\,\,i\in\mathbf{X}\,\\}.$ (85)

The approximations to the pairwise marginals
$\\{\,p_{ij}\,\,|\,\,(i,j)\in\mathbf{E}\,\\}$ can directly be computed from
the elements in the set on the left-hand side of (85), according to the
probability table in Sec. 3.1 (Table 1), where $\xi_{ij}^{\ast}(q_{i},q_{j})$
is computed from (24) as a direct function of $q_{i}$ and $q_{j}$.

To analyze the quality of the Bethe approximation, we compare the exact
partition function $Z(\beta)$, singleton marginals $p_{i}(X_{i}=\pm 1)$, and
pairwise marginals $p_{ij}(X_{i}=\pm 1,X_{j}=\pm 1)$ to the estimated
quantities $Z_{B}$, $(q_{i},1-q_{i})$, and
$(\xi_{ij}^{\ast},q_{i}-\xi_{ij}^{\ast},q_{j}-\xi_{ij}^{\ast},1+\xi_{ij}^{\ast}-q_{i}-q_{j})$,
for all nodes $i$ and all edges $(i,j)$). For that purpose, we use the
following error measures:

  * •

For the partition function, we evaluate the absolute error between the
negative logarithms of $Z(\beta)$ and $Z_{B}$, i.e.,

$\displaystyle|\\!-\log Z(\beta)+\log Z_{B}|.$ (86)

  * •

For the singleton marginals, we evaluate the average (with respect to all
nodes) $l^{1}$-error between the exact and estimated quantities, i.e.,

$\displaystyle\begin{split}&\frac{1}{|\mathbf{X}|}\sum\limits_{i\in\mathbf{X}}\big{(}|p_{i}(X_{i}=+1)-q_{i}|+|\overbrace{p_{i}(X_{i}=-1)}^{=1\,-\,p_{i}(X_{i}=+1)}-(1-q_{i})|\big{)}\\\
=\,&\frac{2}{|\mathbf{X}|}\sum\limits_{i\in\mathbf{X}}\big{(}|p_{i}(X_{i}=+1)-q_{i}|\big{)}.\end{split}$
(87)

  * •

For the pairwise marginals, we evaluate the average (with respect to all
edges) $l^{1}$-error between the exact and estimated quantities, i.e.,

$\displaystyle\begin{split}\frac{1}{|\mathbf{E}|}\sum\limits_{(i,j)\in\mathbf{E}}\big{(}&|\,p_{ij}(X_{i}=+1,X_{j}=+1)-\xi_{ij}^{\ast}\,|\\\
+&|\,p_{ij}(X_{i}=+1,X_{j}=-1)-(q_{i}-\xi_{ij}^{\ast})\,|\\\
+&|\,p_{ij}(X_{i}=-1,X_{j}=+1)-(q_{j}-\xi_{ij}^{\ast})\,|\\\
+&|\,p_{ij}(X_{i}=-1,X_{j}=-1)-(1+\xi_{ij}^{\ast}-q_{i}-q_{j})\,|\big{)}.\end{split}$
(88)

As we do generally not have access to the partition function and marginals
(that we are actually trying to approximate), we evaluate the above error
measures by performing our experiments in Sec. 4.2 on graphs that are small
enough to allow for an exact computation of these quantities. For that
purpose, we use the so-called junction tree algorithm [Lauritzen and
Spiegelhalter, 1988] that returns the correct solution, but whose
computational complexity increases exponentially with the size of the largest
clique 353535A clique is a subgraph where each node is connected to each
other. in the graph.

As discussed in previous sections, the errors induced by the Bethe
approximation strongly depend on the model parameters (in particular the graph
connectivity, potentials, and inverse temperature). We leave this problem to
the experimental analysis in Sec. 4.2 as, we first need to address another
important issue: the practical minimization of the Bethe free energy.
Following the approximations (84) and (85), we need to compute its global
minimum in the constrained set $\mathbb{B}$. This can be a challenging task in
practice, especially if we deal with large systems and numerical difficulties.
Another prohibitive aspect is the potential existence of multiple local
minima. In that case, there is no general technique that efficiently returns
the global minimum, and we must usually rely on the quality of a (possibly
sub-optimal) local minimum 363636For ferromagnetic models, Weller and Jebara
[2014a] have proposed an – at least, in theory – efficient approximation
scheme to compute the partition function up to some $\epsilon$ in polynomial
time (an FPTAS).. Besides, there is no guarantee that the global minimum, if
we manage to find it, induces a low error on the partition function or
marginals. Finally, $\mathcal{F}_{B}$ may exhibit particularly flat or steep
regimes, implying extremely small or exploding gradients and slow convergence
to a minimum (even if $\mathcal{F}_{B}$ is convex on $\mathbb{B}$).

For an efficient minimization of the Bethe free energy, various methods
373737Note that LBP can be seen as a minimization of $\mathcal{F}_{B}$ too
[Heskes, 2003]. However, its convergence rate is often not satisfying – still,