Submitted to the Annals of Statistics
arXiv: arXiv:0000.0000

arXiv:1710.11268v3 [] 11 Dec 2017

THEORETICAL AND COMPUTATIONAL GUARANTEES
OF MEAN FIELD VARIATIONAL INFERENCE
FOR COMMUNITY DETECTION
By Anderson Y. Zhang, and Harrison H. Zhou
Yale University
The mean field variational Bayes method is becoming increasingly
popular in statistics and machine learning. Its iterative Coordinate
Ascent Variational Inference algorithm has been widely applied to
large scale Bayesian inference. See Blei et al. (2017) for a recent comprehensive review. Despite the popularity of the mean field method
there exist remarkably little fundamental theoretical justifications.
To the best of our knowledge, the iterative algorithm has never been
investigated for any high dimensional and complex model. In this
paper, we study the mean field method for community detection under the Stochastic Block Model. For an iterative Batch Coordinate
Ascent Variational Inference algorithm, we show that it has a linear
convergence rate and converges to the minimax rate within log n iterations. This complements the results of Bickel et al. (2013) which
studied the global minimum of the mean field variational Bayes and
obtained asymptotic normal estimation of global model parameters.
In addition, we obtain similar optimality results for Gibbs sampling
and an iterative procedure to calculate maximum likelihood estimation, which can be of independent interest.

1. Introduction. A major challenge of large scale Bayesian inference is
the calculation of posterior distribution. For high dimensional and complex
models, the exact calculation of posterior distribution is often computationally intractable. To address this challenge, the mean field variational method
[2, 19, 30] is used to approximate posterior distributions in a wide range of
applications in many fields including natural language processing [6, 22],
computational neuroscience [14, 26], and network science [1, 8, 17]. This
method is different from Markov chain Monte Carlo (MCMC) [13, 28], another popular approximation algorithm. The variational inference approximation is deterministic for each iterative update, while MCMC is a randomized sampling algorithm, so that for large-scale data analysis, the mean
field variational Bayes usually converges faster than MCMC [7], which is
particularly attractive in the big data era.
MSC 2010 subject classifications: Primary 60G05
Keywords and phrases: mean field, variational inference, Bayesian, community detection, stochastic block model

1

2

ZHANG, ZHOU

In spite of a wide range of successful applications of the mean field variational Bayes, its fundamental theoretical properties are rarely investigated.
The existing literature [3, 8, 31, 33, 34] is mostly on low dimensional parameter estimation and on the global minimum of the variational Bayes
method. For example, in a recent inspiring paper, Wang and Blei [32] studied the frequentist consistency of the variational method for a general class
of latent variable models. They obtained consistency for low dimensional
global parameters and further showed asymptotic normality, assuming the
global minimum of the variational Bayes method can be achieved. However,
it is often computationally infeasible to attain the global minimum when the
model is high-dimensional or complex. This motivates us to investigate the
statistical properties of the mean field in high dimensional settings, and more
importantly, to understand the statistical and computational guarantees of
the iterative variational inference algorithms.
The success and the popularity of the mean field method in Bayesian
inference mainly lies in the success of its iterative algorithm: Coordinate
Ascent Variational Inference (CAVI) [7], which provides a computationally
efficient way to approximate the posterior distribution. It is important to
understand what statistical properties CAVI has and how do they compare
to the optimal statistical accuracy. In addition, we want to investigate how
fast CAVI converges for the purpose of implementation. With the ambition
of establishing a universal theory of the mean field iterative algorithm for
general models in mind, in this paper, we consider the community detection
problem [1, 4, 12, 24, 25, 35] under the Stochastic Block Model (SBM)
[4, 18, 21, 29] as our first step.
Community detection has been an active research area in recent years,
with the SBM as a popular choice of model. The Bayesian framework and
the variational inference for community detection are considered in [1, 3, 8,
11, 17, 27]. For high dimensional settings, Celisse et al. [8] and Bickel et al.
[3] are arguably the first to study the statistical properties of the mean field
for SBMs. The authors built an interesting connection between full likelihood and variational likelihood, and then studied the closeness of maximum
likelihood and maximum variational likelihood, from which they obtained
consistency and asymptotic normality for global parameter estimation. From
a personal communication with the authors of Bickel et al. [3], an implication of their results is that the variational method achieves exact community
recovery under a strong signal-to-noise (SNR) ratio. Their analysis idea is
fascinating, but it is not clear whether it is possible to extend the analysis
to other SNR conditions under which exact recovery may never be possible.
More importantly, it may not be computationally feasible to maximize the

MEAN FIELD FOR COMMUNITY DETECTION

3

variational likelihood for the SBM, as seen from Theorem 2.1.
In this paper, we consider the statistical and computational guarantees of
the iterative variational inference algorithm for community detection. The
primary goal of community detection problem is to recover the community
membership in a network. We measure the performance of the iterative variational inference algorithm by comparing its output with the ground truth.
Denote the underlying ground truth by Z ∗ . For a network of n nodes and
k communities, Z ∗ is an n × k matrix with each row a standard Euclidean
∗ }n
basis in Rk . The index of non-zero coordinate of each row {Zi,·
i=1 gives the
community assignment information for the corresponding node. We propose
an iterative algorithm called Batch Coordinate Ascent Variational Inference
(BCAVI), a slight modification of CAVI with batch updates, to make parallel and distributed computing possible. Let π (s) denote the output of the
s-th iteration, an n × k matrix with nonnegative entries. The summation
(s)
of each row {πi,· }ni=1 is equal to 1, which is interpreted as an approximate
posterior probability of assigning the corresponding node of each row into k
communities. The performance of π (s) is measured by an `1 loss `(·, ·) compared with Z ∗ .
An Informal Statement of Main Result: Let π (s) be the estimation of community membership from the iterative algorithm BCAVI after s iterations.
Under weak regularity condition, for some cn = on (1), with high probability,
we have for all s ≥ 0,
(1)

`(π (s+1) , Z ∗ ) ≤ minimax rate + cn `(π (s) , Z ∗ ).

The main contribution of this paper is Equation (1). The coefficient cn is
on (1) and is independent of s, which implies `(π (s) , Z ∗ ) decreases at a fast
linear rate. In addition, we show that BCAVI converges to the statistical
optimality [35]. It is worth mentioning that after log n iterations BCAVI attains the minimax rate, up to an error on (n−a ) for any constant a > 0. The
conditions required for the analysis of BCAVI are relatively mild. We allow
the number of communities to grow. The sizes of the communities are not
assumed to be of the same order. The separation condition on global parameters covers a wide range of settings from consistent community detection
to exact recovery.
To the best of our knowledge this provides arguably the first theoretical
justification for the iterative algorithm of the mean field variational method
in a high-dimensional and complex setting. Though we focus on the problem
of community detection in this paper, we hope the analysis would shed some

4

ZHANG, ZHOU

light on analyzing other models, which may eventually lead to a general
framework of understanding the mean field theory.
The techniques of analyzing the mean field can be extended to providing
theoretical guarantees for other iterative algorithms, including Gibbs sampling and an iterative procedure for maximum likelihood estimation, which
can be of independent interest. Results similar to Equation (1) are obtained
for both methods under the SBM.
Organization. The paper is organized as follows. In Section 2 we introduce
the mean field theory and the implementation of BCAVI algorithm for community detection. All the theoretical justifications for the mean field method
are in Section 3. Discussions on the convergence of the global minimizer and
other iterative algorithms are presented in Section 4. The proofs of theorems
are in Section 5. We include all the auxiliary lemmas and propositions and
their corresponding proofs in the supplemental material.
Notation. Throughout this paper, for any matrix X ∈ Rn×m
P , its `1 norm is
defined in analogous to that of a vector. That is, kXk1 = i,j |Xi,j |. We use
the notation Xi,· and X·,i to indicate its i-th row and column respectively.
For matrices
P X, Y of the same dimension, their inner product is defined as
hX, Y i =
i,j Xi,j Yi,j . For any set D, we use |D| for its cardinality. We
denote Ber(p) for a Bernoulli random variable with success probability p.
For two positive sequences xn and yn , xn . yn means xn ≤ cyn for some
constant c not depending on n. We adopt the notation xn  yn if xn . yn
and yn . xn . To distinguish from the probabilities p, q, we use bold p and
q to indicate distributions. The Kullback-Leibler divergence between two
distributions is defined as KL(pkq) = Eq log(p(x)/q(x)). We use ψ(·) for the
digamma function, which is defined as the logarithmic derivative of Gamma
d
function, i.e., ψ(x) = dx
[log Γ(x)]. In any Rd , we denote {ea }da=1 to be the
standard Euclidean basis with e1 = (1, 0, 0, . . .), e2 = (0, 1, 0, . . . , 0), . . . , ed =
(0, 0, 0, . . . , 1). We let 1d be a vector of length d whose entries are all 1. We use
[d] to indicate the set {1, 2, . . . , d}. Throughout this paper, the superscript
“pri” (e.g., π pri ) indicates that this is a hyperparameters of priors.
2. Mean Field Method for Community Detection. In this section,
we first give a brief introduction to the variational inference method in
Section 2.1. Then we introduce the community detection problem and the
Stochastic Block Model in Section 2.2. The Bayesian framework is presented
in Section 2.3. Its mean field approximation and CAVI updates are given in
Section 2.4 and Section 2.5 respectively. The BCAVI algorithm is introduced
in Section 2.6.

MEAN FIELD FOR COMMUNITY DETECTION

5

2.1. Mean Field Variational Inference. We first present the mean field
method in a general setting and then consider its application to the community detection problem. Let p(x|y) be an arbitrary posterior distribution
for x, given observation y. Here x can be a vector of latent variables, with
coordinates {xi }. It may be difficult to compute the posterior p(x|y) exactly. The variational Bayes ignores
Q the dependence among {xi }, by simply
taking a product measure q(x) = i qi (xi ) to approximate it. Usually each
qi (xi ) is simple and easy to compute. The best approximation is obtained
by minimizing the Kullback-Leibler divergence between q(x) and p(x|y):
(2)

q̂MF = arg min KL(qkp).
q∈Q

Despite the fact that every measure q has a simple product structure, the
global minimizer q̂MF remains computationally intractable.
To address this issue, an iterative Coordinate Ascent Variational Inference
(CAVI) is widely used to approximate the global minimum. It is a greedy
algorithm. The value of KL(qkp) decreases in each coordinate update:


Y
q̂i = min KL qi
qj p , ∀i.
(3)
qi ∈Qi

j6=i

The coordinate update has an explicit formula


(4)
q̂i (xi ) ∝ exp Eq−i [log p(xi |x−i , y)] ,

where x−i indicates
all the coordinates in x except xi , and the expectation
Q
is over q−i = j6=i qj (xj ). Equation (4) is usually easy to compute, which
makes CAVI computationally attractive, although CAVI only guarantees to
achieve a local minimum.
In summary, the mean field variational inference via CAVI can be represented in the following diagram:
approx.

approx.

p(x|y) ⇐= q̂MF (x) ⇐= q̂CAVI (x),

where q̂MF (x), the global minimum, serves mainly as an intermediate step
in the mean field methodology. What is implemented in practice to approximate global minimum is an iterative algorithm like CAVI. This motivates
us to consider directly the theoretical guarantees of the iterative algorithm
in this paper.
We refer the readers to a nice review and tutorial by Blei et al. [7] for
more detail on the variational inference and CAVI. The derivation from
Equation (3) to Equation (4) can be found in many variational inference
literatures [5, 7]. We include it in Appendix D in the supplemental material
for completeness.

6

ZHANG, ZHOU

2.2. Community Detection and Stochastic Block Model. The Stochastic
Block Model (SBM) has been a popular model for community detection.
Consider an n-node network with its adjacency matrix denoted by A.
It is an unweighted and undirected network without self-loops, with A ∈
{0, 1}n×n , A = AT and Ai,i = 0, ∀i ∈ [n]. Each edge is an independent
Bernoulli random variable with EAi,j = Pi,j , ∀i < j. In the SBM, the value of
connectivity probability Pi,j depends on the communities the two endpoints
i and j belong to. We assume Pi,j = p if both nodes come from the same
community and Pi,j = q otherwise. There are k communities in the network.
We denote z ∈ [k]n , as the assignment vector, with zi indicating the index
of community the i-th node belongs to. Thus, the connectivity probability
matrix P can be written as
Pi,j = Bzi ,zj ,
where B ∈ [0, 1]k×k with diagonal entries as p and off-diagonal entries as q.
That is, B = q1k 1Tk + (p − q)Ik . Let Z ∈ Π0 be the assignment matrix where
Π0 = {π ∈ {0, 1}n×k : kπi,· k0 = 1, ∀i ∈ [n]}.

In each row {Zi,· }ni=1 there is only one 1 with all the other coordinates as 0,
indicating the assignment of community for the corresponding node. Then
T , ∀i < j, or in a matrix form
P can be equivalently written as Pi,j = Zi,· BZj,·
Pi,j = (ZBZ T )i,j , ∀i < j.
The goal of community detection is to recover the assignment vector z,
or equivalently, the assignment matrix Z. The equivalence can be seen by
observing that there is a bijection r between z ∈ [k]n and Z ∈ Π0 which is
defined as follows,
(5)

r(z) = Z, where Zi,a = I{a = zi }, ∀i ∈ [n], a ∈ [k].

Since they are uniquely determined by each other, in our paper we may use
z directly without explicitly defining z = r−1 (Z) (or vice versa) when there
is no ambiguity.
2.3. A Bayesian Framework. Throughout the whole paper, we assume
k, the number of communities, is known. We observe the adjacency matrix
A. The global parameters p and q and the community assignment Z are
unknown. From the description of the model in Section 2.2, we can write
down the distribution of A as follows:
Y Ai,j
p(A|Z, p, q) =
Bzi ,zj (1 − Bzi ,zj )1−Ai,j ,
(6)
i<j

7

MEAN FIELD FOR COMMUNITY DETECTION

with B = q1k 1Tk + (p − q)Ik and z = r−1 (Z). We are interested in Bayesian
inference for estimating Z, with prior to be given on both p, q and Z.
We assume that {zi }ni=1 have independent categorical (a.k.a. multinomial
P
pri n
pri
with size one) priors with hyperparameters {πi,·
}i=1 , where ka=1 πi,a
=
n
1, ∀i ∈ [n]. In other words, {Zi,· }i=1 are independently distributed by
pri
P(Zi,· = eTa ) = πi,a
, ∀a = 1, 2, . . . , k,

where {ea }ka=1 are the coordinate vectors. Here we allow the priors for Zi,· to
be different for different i. If additionally πi,· = πj,· for all i 6= j is assumed,
and then this is reduced to the usual case of i.i.d. priors.
Since {Ai,j }i<j are Bernoulli, it is natural to consider a conjugate Beta
prior for p and q. Let p ∼ Beta(αppri , βppri ) and q ∼ Beta(αqpri , βqpri ). Then the
joint distribution is
(7)
p(A, Z, p, q) =

"
Y
i

×

"


#
Y
A
1−Ai,j 
π pri  Bzii,j
,zj (1 − Bzi ,zj )
i,zi

i<j

Γ(αppri + βppri ) αpri
p p −1 (1
pri
pri
Γ(αp )Γ(βp )

− p)

βppri −1

#"

Γ(αqpri + βqpri )
Γ(αqpri )Γ(βqpri )

q

αpri
q −1

(1 − q)

Our main interest is to infer Z, from the posterior distribution p(Z, p, q|A).
However, the exact calculation of p(Z, p, q|A) is computationally intractable.
2.4. Mean Field Approximation. Since the posterior distribution p(Z, p, q|A)
is computationally intractable, we apply the mean field approximation to
approximate it by a product measure,
qπ,αp ,βp ,αq ,βq (Z, p, q) = qπ (Z)qαp ,βp (p)qαq ,βq (q)
where {r−1 (Zi,· )}ni=1 areQindependent categorical variables with parameters
{πi,· }ni=1 , i.e., qπ (Z) = ni=1 qπi,· (Zi,· ) with
qπi,· (Zi,· = ea ) = πi,a , ∀i ∈ [n], a ∈ [k],

and qαp ,βp (p) and qαq ,βq (q) are Beta with parameters αp , βp , αq , βq due to
conjugacy. See Figure 1 for the graphical presentation of qπ,αp ,βp ,αq ,βq (Z, p, q).
Note that the distribution class of q is fully captured by the parameters
(π, αp , βp , αq , βq ), and then the optimization in Equation (2) is equivalent

βqpri −1

#

.

n×kn,i,jwhere
i,jhyperparameter
i,j
!
!
1 1]prior
1 1]
with
πi,a
∈[n].
[0,
,i,jwhere
n with
prior
hyperparameter
π (a.k.a.
∈k [0,
πi,a =a=1
1, ∀i
∈
It is
k
a=1
n×k
We
assume
{z
}multinomial
have
categorical
(a.k.a.
multinomi
n×k
We
{z
}its
have
categorical
with
size
one)
n
T n{Z
n
iand
iand
Tp,
estimating
with
prior
given
on
both
q
Z.
i=1
imating Z, with
given
on assume
both
p,
qbe
Z.
i=1
equivalently
to
state
that,
}
are
ind
with
hyperparameter
∈
[0,
1]
,
where
π
=
1,
∀i
priorprior
with
hyperparameter
π
∈
[0,
1]
,
where
π
=
1,
∀i
∈
[n].
It
is
equivalently
to
state
that,
{Z
}
are
independently
distributed
with
n
T
i,·
P(Z
=
e
)
=
π
,
∀u
=
1
set D,
we denote
|D|Z,
toprior
cardinality.
set
D,
we
We
denote
denote
|D|
Ber(p)
to
be
its
to
cardinality.
be
a
Bernoulli
We
denote
Ber(p)
to
be
aB
i,a
T
i,·
i,a
!
i=1
!
P(Z
=
e
)
=
π
,
∀u
=
1,
2,
.
.
.
,
k,
i=1
i,·
i,a
equivalently
to
state
that,
{Z
}
are
independ
a=1
a=1
equivalently to state that,
{Zi,·
}ai=1
are)i,a
independently
with
a=
i,·
i,·
P(Z
k πi,a , ∀u
i=1ea ) =
k πi,a , ∀u = distributed
P(Z
=
e
=
1,
2,
.
.
.
,
k,
n×k
i,·
n×k
i,·
a
n
n
n
prior
with
hyperparameter
π
∈
[0,
1]
,
where
π
n
prior
with
hyperparameter
π
∈
[0,
1]
,
where
π
=
1,
∀i
∈
[n].
It
is
i,a
i,a
We
assume
{z
}
have
categorical
(a.k.a.
multinomial
with
size
one)
a=1
T
We assume equivalently
{zi }i=1random
have to
categorical
(a.k.a.
multinomial
with
size
one)
a=1
variable
with equivalently
success
probability
random
variable
p. Forthat,
with
two =
positive
success
probability
sequ
Ti,· }
state
{Z
are
independently
state that,
{Zii,·i=1
}i=1
areto
independently
distributed
with
np.
Ti,· = positive
P(Z
eπdistributed
) ,=∀uπi,a
∀u2,
i=1
=
πsequences
∀u. k. =
1,nx2,
. .i,·
.For
,=
k, etwo
k ,2,
i,a
ai,a
!
P(Z
= ,1,
nP(Zare
kk {Z
!
=P(Z
eTa )i,·
= state
πei,aavectors.
,)∀u
=
1,
. , k,
i,·
a ) =prior
where
{e
}
are
the
coordinate
vectors.
He
kequivalently
k
to
that,
{Z
}
are
independently
dis
equivalently
to
state
that,
}
independently
distributed
with
where
{e
}
are
the
coordinate
Here
we
allow
the
paramn×k
a
n×k
i,·
i,·
a=1
where
{e
}
are
the
coordinate
vectors
a
and
y
,
x
!
y
means
x
cy
for
and
some
y
,
x
constant
!
y
means
c
not
depending
x
cy
for
on
some
n.
We
constant
c
not
depending
o
i=1
where
{e
}
are
the
coordinate
vectors.
Here
we
allow
the
prior
par
i=1
a=1
a
n hyperparameter
n n
nπ
n
n ∈ ,[n].
n isn
a
a=1
with
∈
[0,
1]
where
π
=
1,
∀i
∈
[n].
It
is
or with hyperparameter nprior
π ∈n [0,
1]
, where
π
=
1,
∀i
It
a=1
i,a
i,a
T
a=1
T
k
a=1
k
k
k
where
{e
}
are
the
coordinate
vectors.
eters
on
Z
varies.
If
we
assume
π
=
π
f
where
{e
}
are
the
coordinate
vectors.
Here
we
allow
the
prior
paramP(Z
=
e
)
=
π
,
∀u
=
1,
2,
.
.
.
,
k,
P(Z
=
e
)
=
π
,
∀u
=
1,
2,
.
.
.
,
k,
a
eters
on
Z
varies.
If
we
assume
π
=
π
for
all
i
=
̸
j,
then
this
can
be
i,·
i,·
j,·
a
where
{e
}
are
the
coordinate
vectors.
Here
a=1
i,a
are
the
vectors.
we
prior
i,a
eters
Zfrom
varies.
we
assume
πi,· =
πw
a=1
adopt{Z
the notation
xi,·nto≍where
yeters
xthat,
adopt
the
and
ynncoordinate
!Ifare
xwe
xindependently
To
≍
distinguish
yai,·
ifπHere
xon
!
yπndistributed
and
the
ythe
!
x̸=nwith
.paramTothen
distinguish
a
n
Tallow
on
Zyi,·
varies.
assume
=
for
all
i∀u
j,
this
can
a }!
T
a=1
i,·
na if {e
ni,·
n{Z
ni,·
n
nP(Z
n j,·
n If
a=1
i,·
j,·
equivalently
state
}notation
uivalently to state
that,
distributed
with
=
e
)
=
π
,
=
1,
2,
.
.
.
,
k,
P(Z
=
e
)
π
,
∀u
=
1,
2,
.
.
.
,
k,
i,·=
6
ZHANG,
6
ZHOU
ZHANG,
ZHOU
i,· }i=1 are independently
i,·
i,a
i,·
i,a
i=1
a
n
a
eters
on
Z
varies.
If
we
assume
π
=
π
eters
on
Z
varies.
If
we
assume
π
=
π
for
all
i
=
̸
j,
then
this
can
be
i,·
i,·
j
eters
Z
varies.
Ifusual
we
assume
πj,· sam
for
i,·and
j,·ip
into
the
case
with
pria
onnation
Zi,·into
varies.
Ifusual
we
assume
π
=
πon
for
all
̸=on
j,
then
this
be=
i,·
i,·
the
case
with
same
prior
all
.}πni=1
i,·degenerated
j,·i,·
probabilities p, q, wedegenerated
useeters
bold
probabilities
qp,to
q,
indicate
wedegenerated
use
bold
distributions.
nation
and
q {z
to
indicate
distributi
into
the
usual
case
with
i}
degenerated
the
usual
case
with
same
prior
on
all
{zcan
.same
i=1
k pinto
nicase
nthe
degenerated
into
usual
with
same
degenerated
into
the
usual
case
with
same
prior
on
all
{z
}
.
degenerated
into
the
usual
case
with
same
prior
where
{e
}
are
the
coordinate
vectors.
Here
we
allow
the
pri
where {ea }ka=1
are
the
coordinate
vectors.
Here
we
allow
the
prior
paramdegenerated
into
the
usual
case
with
same
prior
on
all
{z
}
.
i
T
T
a
k
i=1
i
k
Due
to
the
fact
that
p,
q
∈
(0,
1)
and
all
the
i=1
a=1
Te +)(p
Te
Due
toP(Z
the
that
p,−
qπq)I
∈
and
all
the
{A
}Bayesian
random
Due
to
fact
that
p,Bernoulli
qare
∈
(0, 1)
and
all{o
=
)(pvectors.
=
,(0,
=
1,
2,
.Due
. all
.all
,the
k,
Due
to
fact
that
∈1)
(0,
1)
and
the
{A
}are
Bernoulli
ran
{e
}p,ka=1
are
the
coordinate
vectors.
Here
we
allow
i,j
i<j
}ka=1
the
coordinate
Here
we
allow
the
prior
param=−{eπq)I
,. ∀u
=
1,
2,B
.=
.i,·
.fact
,the
k,
withP(Z
B =i,·q1=
Weare
are
with
interested
q1
1where
aa+fact
Bayesian
framework
.q∀u
We
are
interested
for
inthe
athat
framework
for
i,j
i<j
ai,a
ai,a
k 1where
k in
to
fact
that
p,
q
∈
(0,
1)
and
all
th
ka
kthat
Due
to
the
that
p,
q
∈
(0,
1)
and
the
{A
}
are
Bernoulli
random
Due
to
the
fact
p,
q
∈
(0,
1)
and
all
the
{A
Due
to
the
fact
p,
q
∈
(0,
1)
and
all
the
{A
}
are
Bernoulli
random
i,j
i<j
i,
i,j
i<j
eters
varies.
If=
we
πqis
=then
for
all
iprior
̸=
j,
then
th
eters
on Z
varies.
IfZgiven
we
assume
πZHANG,
=Mean
πZHOU
for
all
i assume
̸= variable,
j,
this
can
beconsider
2.8i,·
Field
Method.
2.
Field
Method.
variable,
it
natural
to
the
conjuga
i,·
i,·
j,·this
j,·π
it
isi,·
natural
to
consider
the
conjugate
prior:
Beta
for
it
isitZ.πnatural
to
consider
the
con
variable,
it
is
natural
to
consider
the
conjugate
Beta
prior
for
b
eters
on
Zgiven
Ifthen
we
assume
πprior:
=
πfor
for
all
i both
̸=
j,
estimating
Z,Mean
with
prior
onvariable,
estimating
both
p,
qZ
and
Z,
with
Z.
prior
on
both
p,
and
eters
on
varies.
Ifon
we
assume
πtoj,·varies.
for
all
ithe
̸=
j,
can
be
i,·
i,·prior:
j,·
i,·k variable,
i,· consider
variable,
is
natural
to
consider
the
conj
variable,
it
is
natural
to
consider
the
conjugate
variable,
it
is
natural
consider
conjugate
Beta
prior
for
both
it
is
natural
to
the
conjugate
prior:
Beta
prior
both
n .
nwe. allow
nthe
n
n
degenerated
into
the
usual
case
with
same
prior
on
all
{z
}
where
{e
}
are
the
coordinate
vectors.
Here
the
prior
paraminto
usual
case
with
same
prior
on
all
{z
}
p
and
q.
We
let
p
∼
Ber(α
,
β
)
and
q
∼
Ber
ere {ea }ka=1 degenerated
are We
theassume
coordinate
vectors.
Here
we
allow
the
prior
paramp
and
q.
We
let
p
∼
Ber(α
,
β
)
and
q
∼
i
a
{z
}
have
categorical
We
(a.k.a.
assume
multinomial
{z
}
have
with
categorical
size
one)
(a.k.a.
multinomial
with
size
one)
p
and
q.
We
let
p
∼
Ber(α
,
β
)
and
q
∼
Ber(α
,
β
)
with
hyperparameter
i
phyperparameter
pand
and
q. We
pi=1p∼same
))∼and
q∼
∼
Ber(α
,q∼
β
)p, βwith
hyperparam
pand
p
a=1
i i=1
ilet
degenerated
the
usual
case
same
allq ∼
{zqB
pβpinto
pand
degenerated
intop the
case
with
all
}βBer(α
.with
pi=1
and
let
phyperparameter
and
and
q.
We
∼
Ber(α
∼i=1
Ber(α
pp,p
pusual
q.
We
∼Ber(α
,ββ
and
q{z
, qβ
) qwith
q.and
We
let
plet
∼
Ber(α
)prior
qppon
Ber(α
)pqwith
i
p, β
p )qon
i,let
pBer(α
p )prior
qq
p ,Ber(α
qq.
qWe
i=1
!
!
pri Block
pri
pri
pri
2.1.
Stochastic
Since
2.1.
Stochastic
proposed,
Block
Stochastic
Model.
Block
Since
Model
pri Model.
k
k proposed, Stochastic Block
n×k
n×k
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
Due
to
the
fact
that
p,
q
∈
(0,
1)
and
all
the
{A
}
are
Bernou
α
,
β
,
α
β
.
Thus
the
full
likelihood
fun
eters
on
Z
varies.
If
we
assume
π
=
π
for
all
i
=
̸
j,
then
this
can
be
prior
hyperparameter
π
∈
[0,
prior
1]
with
,
where
hyperparameter
π
=
1,
π
∀i
∈
∈
[0,
[n].
1]
It
is
where
π
=
1,
∀i
∈
[n].
It
is
Due
to
the
fact
that
p,
q
∈
(0,
1)
and
all
the
{A
}
are
Bernoulli
random
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
rs on Zi,· varies.
Ifwith
we
assume
π
=
π
for
all
i
=
̸
j,
then
this
can
be
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
funct
p
p
q
q
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
α
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
p
p
q
q
i,a
i,a
i,j the
i<j{Ai,j }i<j are
parameter
i,· fact
j,· pthat
i,j
i<j
pq are
q a=1
q1) andrandom
Due
toi,·
the
fact
p,
qpq∈ Bernoulli
(0,
all
p j,·pthat
p p p q pp,
Due toi,·the
and
all
the
{A
}pi<j
q q1)
q qqpqq∈ q(0,
a=1
i,jstudied
(SBM)tohas
been
the
most
studied
(SBM)
model
has
fornbeen
community
the
most
detection.
model
Consider
for
community
detection. C
n equivalently
nwith
n"⎤
⎡
⎡
⎤
⎡
⎡⎡ with
⎡
⎤
equivalently
state
that,
{Z
}
are
independently
to
state
distributed
that,
{Z
}
are
independently
distributed
⎡
⎤
⎡
"
$
i,·
i,·
$
"
$
"
$
"
$
variable,
it
is
natural
to
consider
the
conjugate
prior:
Beta
prio
degenerated
into
the
usual
case
with
same
prior
on
all
{z
}
.
"
$
i=1
i=1
variable,
it is
natural
to
consider
the
conjugate
prior:
Beta
prior
for
both
generated into
the usual
case
with
same
prior
on
all
{z
}
.
variable,
it
is
natural
to
consider
the
conjugate
prior:
Be
"
$
"
$
variable,
it
is
natural
to
consider
the
conjugate
prior:
Beta
prior
for
both
i
i
i=1
i=1
#un#matrix
an n-node network with its adjacency
an n-node
network
denoted
with
its
A.adjacency
It1−A
is an
matrix
denoted
by
A. It1−is
# #
#
##
#by
#
#
#
#
A
#
#
#
#
A
A
i,j
A
i,j
i,j
i,j
latent
i,j
1−A
Tfact
T
A
Ai,j
⎣qπ
⎦ .)n,n
i,j ⎣
⎣
AB
⎣
p(A,
q)
=
πi,j
(1
)BB
i,j
plet
q.
We
let
p, B
∼
Ber(α
,zare
βihyperparameter
and
qi,z
∼
Ber(α
, i,j
βB
with
p) pp(A,
q.
pself
Ber(α
β
)Z,
and
qB
∼
Ber(α
, iβB⎣
hyper
p,weighted
π
BBer(α
(1
)=
Due
toi,·We
the
q=pq)
∈
(0,
1)
and
the
{A
}jwith
random
1−A
plet
and
q.
let
pand
∼
Ber(α
, are
β
).=
and
∼,=all
,i,a
β
)p,
p and
pP(Z
∼
Ber(α
,network
βπthat
)p(A,
and
qp,
∼
Ber(α
β
hyperparameter
i,j
Z,
p,
q)
=
π
B(1
(1
z)
z)
,z−
1−A
Due to the fact
thatq.
p, We
q∈
(0,
1)
and
{A
}q)We
random
Z,
q)
B
(1
−
B
za
,z
=all
ep(A,
=
,i,j
∀u
=
1,
2,
.Bernoulli
, k,
P(Z
)Ajpwith
=
π
,p(A,
1,
2,
.Bernoulli
k,
i,j
i,z
,z
i ,zwith
p
qz−
⎣
iB
j−
p
qi,z
qB
i<j
z
p.and
q−
qzi,j
pZ,
qi,z
qp(A,
i)
⎣
⎦
i∼
iZ,
,z
i<j
,zp.j)
jq
i,a
i,·
⎣
j∀u
p,
q)
=
i e
a the
Z,
q)
π
−
weighted
without
and
such
that
network
∈
without
{0,
1}
self
,⎦i π⎦
loop,
such
A
p(A,
Z,
p,
=
πloop,
(1
−
),i,z
zii,z,z
variableand undirected
Z,
p,
q)p,
=
π
B
−
B
zij,zB
z∈
zziip,
,z
j (1
i,z
i
j that
zi ,zp(A,
i ,zj{
i,zi undirected
z=
jA
i
j (1
i ,zj )zi i,zj i,zi
i<j
i
i<j
i is
i<jboth
i
i<j
T
T
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
n
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
n
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
variable,
it
is
natural
to
consider
conjugate
prior:
Beta
prior
for
α
,
β
,
α
,
β
.
Thus
the
full
likelihood
function
is
p
p
q
q
p and
p the
q i,i qconjugate
iable, it is natural
to }consider
prior:
for
i<j
p∈ [n].
q =
pwhere
p A
i) *allow
*)
= 0,p∀i
A
and
A
isthe
an=
independent
0,both
∀i
∈ [n].
Bernoulli
edge
isi<jan
independent
B
i )prior
i<j
k qA are
kq AiBeta
)areedge
*Here
)Each
*
i<j
i,i
)i )+ the
{eqa=
the Acoordinate
where
vectors.
{e
Here
}Each
allow
the
coordinate
prior
paramvectors.
we
prior
parama
a=1
a=1)we
βqqp−1
)p +
)βqpΓ(α
Γ(α
βp ) p α+
Γ(α
+*
) Γ(α
)
)
)
⎡
β
)
p +Γ(α
q Γ(α
⎡
⎤
αqp⎤
−1
β⎤
β
)
+
β
)
*
)
−1
β
−1
α
β
−1
p −1 β β
p
p
p
q
p
⎡
⎡
⎤
α
−1
"
$
α
−1
β
−1
α
−1
observation
p
"
$
p, β
and
q.weWe
pπi,ji,·∼=
Ber(α
,all
β"
and
qwith
∼βthis
Ber(α
,πPβj,·−
with
hyperparameter
pp)=
q−
p(1β
−(1
×pall
(1
p)
random
variable,
with
EA
P
∀i×
variable,
j.
EA
, ∀i
< Γ(α
nd q. We let eters
p ∼ on
Ber(α
q $∼let
Ber(α
,random
with
hyperparameter
p$
−i,·
q+then
−q)q)
.− −
(1π
p<
pi)×
q(1
q )for
i,j
i,j
i,j
"p ) Ifand
(1
p)−
pcan
p)
q)
Zi,· pvaries.
assume
eters
=qon
πβ
Zq,)#
for
varies.
̸=If
j,
weβ
then
assume
=
be
ij.×
̸=ppΓ(α
be
Γ(α
+
Γ(α
)pqcan
+
βppβ
))Γ(β
Γ(α
+
))Γ(β
) this
p)p
pp )
qq)−1
p)q+
αpα−1
αp −1
βpqΓ(α
−1
Γ(α
)Γ(β
α)pβ
β
n(n-1)/2j,· i,·
n(n-1)/2
p)Γ(β
qj,
q+
Γ(α
Γ(α
)Γ(β
)
q −1
#
#
pp−1
p
#
α
−1
β
−1
α
β
−
p
p
q
Γ(α
p
q
Γ(α
)
Γ(α
)Γ(β
p
−
×
(1
p
pconnectivity
−the
q(1−
×
(1
p) on
(1
#
#
p 1−A
qn Ppqq depends
×
(1
p)
Aall
np
A
−⎦
−⎦
(1πsame
p)
q)
i,j {z }
In SBM,
of#
connectivity
In×into
SBM,
probability
the
P
of
depends
on⎣
probability
comon−pt
i,j
1−A
i,j
i,j
αlikelihood
β#
, usual
αp,qvalue
,q)
βfunction
. ⎣Thus
the
full
likelihood
function
is
i,j
A
degenerated
into
case
with
degenerated
same
prior
on
the
all
{z
usual
}value
case
.p⎣
with
prior
.z)i,j⎦
βp , αq , βq . Thus
the full
is
Ai,j
⎣
p , the
pthe
q=
i,j
i−
i−
Γ(α
)Γ(β
)
1−A
1−A
p(A,
Z,
p,
q)
=
B
(1
B
)
Γ(α
)Γ(β
Γ(α
)Γ(β
)
i=1
i=1
p(A,
Z,
π
B
(1
B
)
i,j
Γ(α
)Γ(β
)
p
p
i,j
full
Bayesian
inference
mean
field
approximation
z
,z
i,z
,z
p
p
q
q
z
,z
⎦
i,z
z
,z
Γ(α
)Γ(β
)
Γ(α
)Γ(β
p
p
i
j
i
i
j
i
j
i
i
j
p obtaining
q) Z,
q p,
p(A,
Z,
p,
q)
=
B
(1
B
p(A,Due
Z,munities
p, q) fact
= thattwo
π q endpoints
B
(1
−
Bthe
)πp,Z,
interest
to
estimate
by
obtaining
posterio
Our
interest
to
estimate
p(Z,
i,z
zi ,zjestimate
zall
,z
zthat
,zpto.
i ,z
j=
i∈
iand
jinterest
ito
jestimate
Our
interest
Z,q|A).
by
obtaining
iOur
munities
jis
belong
two
We
endpoints
assume
pzall
iposterior
and
pis
j−
ifdistribution
belong
to.
We
assume
pp,
= pp
Z,
by
obtaining
posterior
distribution
p(Z,
to$the
p, i,z
∈i (0, 1)
Due
to
the
the
{A
fact
are
Bernoulli
q by
(0,
1)Our
and
the
{A
}is
Bernoulli
random
i,j
i,jq|A).
⎡
⎤
⎡ the
⎤
i,j }is
i<j
i,jboth
i<jtoare
irandom
i<j
iand
i<j
"
$
"
Same
as
many
other
other
problems,
calculating
Same
as
many
other
other
problems,
calculating
p(Z,
p,
q|A)
is
computationOur
interest
is
to
estimate
Z,
by
obtaining
i
i<j
i
i<j
Our
interest
is
to
estimate
Z,
by
obtaining
posterior
distribution
p(Z,
p,pq
Same
as
many
other
other
problems,
calcula
Same
as
many
other
other
problems,
calculating
p(Z,
p,
q|A)
is
computationOur
interest
is
to
estimate
Z,
by
obtaining
post
nodes
come
from
the
same
community
nodes
come
and
from
P
the
=
q
same
otherwise.
community
Usually
and
P
=
q
otherwise.
)
*
)
variable, itFig
is#
natural
to model
consider
variable,
the conjugate
it to
is#
natural
prior: Beta
to
consider
prior
for
both
conjugate
prior:
Betai,jprior p(Z,
for
)#
*the
)panel)
* both
Our
interest
is
estimate
Z,
by obtaining
posterior
distribution
p, q|A).
i,j
#
1.
Graphical
presentations
of
full
Bayesian
inference
(left
and
the
mean
A
)
*
)
ally
intractable.
)
*
)
*
A
ally
intractable.
i,j
1−A
i,j
Γ(α
β⎦β
Γ(αq + β
Γ(α
+i,j
β⎦
+
i,j
pas
p )q ),α
ally
intractable.
p∼
p )let
q+
⎣other
Same
many
other
intractable.
⎣ p(A,
β
as
many
other
problems,
calculating
p(Z,
p,
is .calculatin
computat
α
−1+
pq−1
pq−1
q)−1
let
pB∼
Ber(α
,Same
βB
and
pSame
and
qally
q.
Ber(α
We
ppB
)∼
Ber(α
hyperparameter
,βSame
βpB
)qzpiand
qqmany
Ber(α
βother
with
hyperparameter
p, q)p(right
=
π
(1
−
)βαΓ(α
other
other
A, Z, p, q) = p and
πi,zq.i We
)α1−A
p×
q−1
p+
p−1
qpp(Z,
field
approximation
for
The
edges
show
the
dependence
many
other
problems,
calculating
p,(1
isβproblems,
computationziwith
,z
i,z
,z
Γ(α
β
)+
Γ(α
βq ) calcu
zZ,
Γ(α
+
β)zppanel)
) jas
Γ(α
)∼
−problems,
×
(1
p)
j−
icommunity
jas
pq ,αβother
−
(1
q)q|A)
i ,z
j (1p−
pp)
q
i ,z
βdetection.
βq|A)
p −1 αq −1 q βp −1
q−
p −1
p −1
q −1
pΓ(α
qq α)Γ(
×Approximation.
(1 −
Γ(α
)Γ(β
)q(1
Γ(α
.
×
p)
q) of posterior
Γ(α
)isThus
)p)−
ally
intractable.
among
variables.
ally
αip , βp , αq , β
Thus
the full likelihood
α
, )Γ(β
αq(1
, βqField
.i<j
the full likelihood
function
is
pMean
pField
p
p−
qq)Γ(β
2.3.
Approximation.
The
compu
q .i<j
pp, βfunction
pintractable.
ally
intractable.
2.3.
Mean
The
computational
issue
i
ally
intractable.
Field Approximation.
) *q)) 2.3.
Γ(α
Approximation.
TheMean
issue
of posterior
Γ(α⎡p )Γ(βp )) * )2.3. Mean Field
Γ(αqdistribution
p )Γ(β
p)Γ(β
q )Γ(β
q ) The
⎤ Γ(α
⎡
⎤
* c
)
* isincomputational
one
of the
keyThe
challenges
arose
distribution is one
the$ key
challenges
arose
Bayes
world.
so- challenges
"
$
" isof
distribution
is
one
of
the
key
distribution
one
of
the
key
challenges
arose
in
Bayes
world.
The
soOur
interest
is
to
estimate
Z,
by
obtaining
posterior
distrib
Our
interest
is
to
estimate
Z,
by
obtaining
posterior
distribution
p(Z,
p,
q|A).
Γ(α
+filed
βq p+
)Field
Γ(α
+
βApproximation.
Γ(αp + βp ) #αp −1 # Aβcalled
βαqp)−1Approximation.
#
#
pΓ(α
q1−A
q ) Approximation.
2.3.
Mean
Field
The
2.3.
Mean
The
computational
issue
βcalled
−1β
αway
βThe
−1poste
mean
filed
variational
Bayes
method
is
o
−1interest
α(1
−1
p
q −1
qof
mean
variational
Bayes
method
is
one
popular
to
tackle
A
p
q⎣
q −1
2.3.
Mean
Field
com
i,j
i,j
1−A
2.3.
Mean
Field
Approximation.
The
computational
issue
of
posterior
i,j
i,j
Our
is
to
estimate
Z,
by
obtaining
posterior
distribution
⎣
⎦
⎦
p
−
q
−
.p(
×
p)
(1
q)
Our
interest
is
to
estimate
Z,
by
obtaining
posterior
distribution
p(Z,
p,
q|A).
p
−
q
−
.
×
(1
p)
(1
q)
called
mean
filed
variational
Bayes
metho
called
mean
filed
variational
Bayes
is
onecalculating
popular wayp(Z,
to tackle
p(A, Z, p, q)
πas
B
p(A,
−Z,Bp,zi q)
πas
(1other
− p,
Bmethod
zparameters
zi ,zp(Z,
i,ziover
,z
i,z
zq|A)
Same
other
problems,
p,
q|A
to =minimize
the
as
Same
many
other
problems,
calculating
computationi ,z
j (1other
j The
j )=
i many B
i ,z
j ) is
it.
philosophy
of
mean
filed
theory
is
to
appr
distribution
is
one
of
the
key
challenge
it.
The
philosophy
of
mean
filed
theory
is
to
approximate
the
posterior
disdistribution
is
one
of
the
key
challenges
arose
in
Bayes
world.
The
Γ(α
)Γ(β
)
Γ(α
)Γ(β
)
Γ(α
Γ(α
)Γ(β
p iof
q Bayes
p ) other
qone
q ) the
is
one
of
the
keyworld.
challenges
distribution
is
key
arose
in
The
soit. The
ofp(Z,
mean
filed
theory
isart
it.
philosophy
ofp(Z,
mean
filed
theory
isqphilosophy
to approximate
the
posterior
disi<j
i<jdistribution
Same
aspThe
many
other
other
problems,
calculating
p,
q|A)
ismeth
com
Same
asp )Γ(β
many
other
problems,
calculating
p,challenges
q|A)
computationally
intractable.
ally)i intractable.
tribution
p(x|y)
byissome
distribution
q(x)
with
tribution
p(x|y)
by
some
distribution
q(x)iswith
simpler
structure.
When
called
mean
filed
variational
Bayes
*
)
)
*
*
)
*
called
mean
filed
variational
Bayes
method
one
popular
way
to
ta
tribution
byβ some
distribution
q(x
mean
filed
variational
Bayes
method
tribution
p(x|y)
byβ+qsome
distribution
q(x)
with
simpler
structure.
When
(8)Our Γ(α
mean
filed
variational
Bayes
method
isp(x|y)
one
way
to
Γ(α
Γ(α
)βinvolved
)called
Γ(α
+
intractable.
are
multiple
variables
involved
inβtackle
p + βp )called
q +p h
pα
q popular
q) α
ally
intractable.
ispally
to
estimate
Z,
by
obtaining
posterior
distribution
p(Z,
p,
q|A).
i
are
variables
in
the
distribution
p(x),
r interest is to
estimate
Z, ×by interest
obtaining
posterior
distribution
p(Z,
p,
q|A).
α
−1
βmultiple
α
−1
βq −1
βphilosophy
−1
−1 po
pthere
p −1
q −1
pthere
pposterior
q −1
qthe
it.
The
of
mean
filed
theory
it.
of variables
mean
tomultiple
approximate
the
− p) philosophy
q it.
p filed
−
−there
.is the
q distribution
−posterior
(1The
×
(1involved
(1theory
q)
p) in
(1involved
q) isp(x),
are
variables
in.isat
there
aremin
multiple
posterior
MF
MF
MF it.
MF The
MF
The
philosophy
of
mean
filed
theory
to
philosophy
of
mean
filed
theory
is
to
approximate
the
posterior
dismean
field
theory
usually
assumes
the
independ
(π̂
,
α̂
,
β̂
,
α̂
,
β̂
)
=
arg
KL
q
(Z,
p,
q)
p(Z,
p,
q|A)
,
2.3.
Mean
Field
Approximation.
The
computational
Γ(α
)Γ(β
)
Γ(α
Γ(α
)Γ(β
)Γ(β
)
)
Γ(α
)Γ(β
)
mean
field
theory
usually
assumes
the
independence
of
these
variables
in
2.3.
Mean
Field
Approximation.
The
computational
issue
of
posterior
π,α
,β
,α
,β
p
p
q
p
q
p
q
q
p
p
q
q
p
p
q
q
Same as many
other
other
problems,
calculating
p(Z,
p,
q|A)
is computation+n q(x)
me as many other other problems,
calculating
p(Z,
p,
q|A)
is computation+n field
tribution
p(x|y)
by
some
distribution
q(
mean
field
theory
usually
assumes
the ind
tribution
p(x|y)
by
some
distribution
with
simpler
structure.
W
mean
usually
assumes
the
independence
oftothese
variables
in
π∈Π
1 theory
tribution
p(x|y)
by
some
distribution
q(x)
w
q(x)computational
=
q+
(x
),
a way
reduce
computat
q(x)
=
q
(x
a way
to is
reduce
The
so+
tribution
p(x|y)
by
some
distribution
q(x)
with
simpler
structure.
When
iposterior
idistribution
iq),
i=1
one
of
the
key
challenges
arose
in
Baye
n icomplexity.
distribution
is
one
of
the
key
challenges
arose
in
Bayes
world.
The
so2.3.
Mean
Field
Approximation.
The
computational
issue
of
i=1
2.3.
Mean
Field
Approximation.
The
computational
issue
of
n
α
,β
,α
,β
>0
p
p
q
to estimate
Z, by obtaining
Our
interest
posterior
is to
estimate
distribution
Z,
by
p(Z,
obtaining
p, q|A).
posterior
distribution
p(Z,
p,over
q|A).
q(x)
=
qposterior
ato way
toinvolved
reduce
com
intractable.
there
are
multiple
inp
q(x)
=field
qi (xivariables
), to
a search
way
to
reduce
computational
complexity.
The
soy intractable. Our interest isally
i (xvariables
i ),
there
are
multiple
involved
in
the
distribution
i=1
i=1
called
mean
field
theory
is
search
some
called
mean
theory
is
over
some
distribution
class
Q,
such
there
are
multiple
variables
involved
in
the
there
are
variables
involved
in
the
posterior
distribution
p(x),
n×k
called
filed
variational
Bayes
method
is
one
popu
mean
filed
variational
Bayes
method
is
one
popular
way
to
tackle
distribution
is
of
key
challenges
arose
in
Bayes
world
distribution
is
one
of problems,
the
key
challenges
arose
in
Bayes
world.
The
socalled
mean
field
theory
to
search
overin
Same as manycalled
other
other
Same
calculating
many
p(Z,
other
p,
other
q|A)
problems,
computationcalculating
p(Z,
p,
q|A)
is
computationmean
field
theory
isassumes
to
over
some
class
Q,
such
where
Π
=asmultiple
{π
∈ [0,
1]one
,mean
kπ
kthe
1}.
mean
field
theory
usually
the
1called
i,·is
mean
field
theory
usually
the
independence
ofisassumes
these
variable
that
the
Kullback-Leibler
divergence
isthe
minimized
1 =
that
the
Kullback-Leibler
divergence
issearch
minimized.
That
isdistribution
+
mean
field
theory
usually
assumes
indep
+
mean
field
theory
usually
assumes
the
independence
of
these
variables
in
n
it.
The
philosophy
of
mean
filed
theory
is
to
approximate
it.
The
philosophy
of
mean
filed
theory
is
to
approximate
the
posterior
disn
that
the
Kullback-Leibler
divergence
is
min
ally
intractable.
ally
intractable.
that
the
Kullback-Leibler
divergence
is
minimized.
That
is
+=nway
called
mean
filed
variational
Bayes
method
popular
way
calledApproximation.
mean filed
Bayes
is of
popular
toqissue
2.3. variational
Mean
Approximation.
The
computational
posterior
+
2.3. Mean Field
TheField
computational
posterior
q(x)
(xiis
), of
aone
way
to reduce
co
q(x)
=
qissue
aone
way
to
reduce
computational
complexity.
The
itackle
MF
MF
i (xi ), q̂
n method
i=1
i=1
(4)
q̂
=
min
KL(q∥p).
q(x)
=
q
(x
),
a
way
to
reduce
compu
(4)
=
min
KL(q∥p).
Here
Π
can
be
viewed
as
a
relaxation
of
Π
:
it
uses
an
`
constraint
q(x)
=
q
(x
),
a
way
to
reduce
computational
complexity.
The
soi
i
1
0
1i=1 distribution
tribution
p(x|y)
by
some
q(x)
with
simpler
i
i
tribution
p(x|y)
by
some
distribution
q(x)
with
simpler
structure.
When
MF
i=1
MF
q∈Q
q∈Q
The
philosophy
of issue
mean
filed
theory
is
to world.
approximate
the
pos
called
mean
field
theory
is
search
ove
The
philosophy
of mean
filed
theory
is
tochallenges
approximate
the
posterior
disdistribution
isit.
one
of
the
key
arose
in
Bayes
The
so(4)
q̂to
=
min Q,
KL(
called
mean
field
theory
is
to
search
over
some
distribution
class
tribution is it.
one
of Mean
the
key row
challenges
arose
in
Bayes
world.
The
so(4)
q̂
=The
min
KL(q∥p).
2.3.
Field
Approximation.
2.3.
The
computational
Field
Approximation.
of
posterior
issue
ofsearch
posterior
on each
instead
of the
`0 Mean
constraint
used
inin
Πthe
.search
The
global
minimizer
mean
field
theory
isin
to the
over
som
q∈Qsuch
0called
q∈Qcomputational
called
mean
fieldthere
theory
is
to
over
some
distribution
class
Q,
are
multiple
variables
involved
posterior
there
are
multiple
variables
involved
posterior
distribution
p(x),
that
the
Kullback-Leibler
divergence
isp(Z
md
that
the
Kullback-Leibler
divergence
isour
minimized.
That
is
Back
top(Z,
case,
the
posterior
distribution
tribution
p(x|y)
by
some
distribution
q(x)
simpler
structu
tribution
p(x|y)
by
some
distribution
q(x)
with
simpler
structure.
When
Back
to our
case,
the
posterior
distribution
p,
q|A)
iswith
the
main
intercalled
mean
filed
variational
Bayes
method
is
one
popular
way
to
tackle
led mean filed
variational
Bayes
method
is
one
popular
way
to
tackle
distribution
is
one
of
the
key
challenges
distribution
arose
is
one
in
Bayes
of
the
world.
key
challenges
The
soarose
in
Bayes
world.
The
soq
(Z)
gives
approximate
probabilities
to
classify
every
node
to
each
comMF
that
the
Kullback-Leibler
divergence
is
minim
π̂ mean field theory
thatestthe
Kullback-Leibler
divergence
is
minimized.
That
is
mean
field
theory
usually
assumes
the
independence
of
usually
assumes
the
independence
of
these
variables
in
Back
to
our
case,
the
posterior
distributi
Back
ourvariational
case,+
theway
posterior
distribution
p, q|A)way
is the
main
interestapproximate
but
it is computationally
intractable.
Conside
but
it
is computationally
Consider
called
mean
filed
variational
method
mean
isto
filed
one
popular
to
tackle
method
is distribution
onep(Z,
popular
to
tackle
+n isBayes
MFdisthere
multiple
variables
involved
in
posterior
distribut
there
multiple
variables
involved
in
posterior
distribution
p(x),
it.
The
philosophy
of
mean
filed
theory
isBayes
to
the
posterior
MF
munity.
The
optimization
inare
Equation
(8)the
can
benintractable.
shown
to
be
equivalent
to
The philosophy
of are
mean
filed
theory
tocalled
approximate
posterior
dis(4)
q̂
=
min
KL
(4)
q̂
=
min
KL(q∥p).
est
but
it
is
computationally
intractable.
C
est
but
it
is
computationally
intractable.
Consider
distribution
MF
q(x)
=
q
(x
),
a
way
to
reduce
computational
com
q(x)
=
q
(x
),
a
way
to
reduce
computational
complexity.
The
soMF
i theory
i disi theory
i The
i=1
i=1
q∈Q ′ va
it.
The
philosophy
of
mean
filed
it.
is
philosophy
to
approximate
of
mean
the
filed
posterior
is
to
approximate
posterior
dis(4)
q̂
=
min
KL(q∥p
q∈Q
(4)
q̂
=
min
KL(q∥p).
a
more
explicit
optimization
as
follows.
Recall
ψ(·)
is
the
digamma
function
′independence
′ ′,βstructure.
′ ′,α
′ ,β
′the
′ (Z)q
q
(Z,
p,
q)
=
q
′
′
′
′
′
′
′
′
q
(Z,
p,
q)
=
q
(Z)q
(p)q
(q)
mean
field
theory
usually
assumes
the
of
these
mean
field
theory
usually
assumes
the
independence
of
these
variables
in
π
,α
π
α
tribution
p(x|y)
by
some
distribution
q(x)
with
simpler
When
π ,αp ,βpstructure.
,αq ,βq
π
αp ,βp
α
bution p(x|y) by some
distribution
q(x)
with
simpler
When
p q ,β
pq q q
p ,β
q∈Q
+
mean
field
theory
to search
over
some
distribut
mean
field
theory
isp(x|y)
tocalled
search
over
someq∈Q
distribution
class
such
d distribution
tribution+
p(x|y)
by some
tribution
q(x)
simpler
by
When
simpler
structure.
When
′Q,
′ ,α
′ ,β
′ (Z, p,
q)
= qπ′ (Z
nwith
′ ,βstructure.
′ distribution
′ ,β ′ (Z, p,
′q
′ ′(p)q
′ (q)
n called
qπ′some
q) q(x)
= qis
(Z)q
with
ψ(x)
=
Γ(x)].
,α
π ′with
α
α′q ,β
p ,βp
qqq
p a
p ,αway
q in
q the
p ,βπ
p ,α
dxa[log
q(x)
=
q
(x
),
to
reduce
computational
complexity
q(x)
=
q
(x
),
way
to
reduce
computational
complexity.
The
soBack
to
our
case,
the
posterior
distribu
there
are
multiple
variables
involved
posterior
distribution
p(x),
Back
to
our
case,
the
posterior
distribution
p(Z,
p,
q|A)
is
the
main
in
re are multiple
variables
involved
in
the
posterior
distribution
p(x),
i
i
i
i
i=1
i=1
that
the
Kullback-Leibler
divergence
is
minimized.
That
i
that the
Kullback-Leibler
divergence
is
minimized.
That
is
there are multiple
variables
involved
there
in
are
the
multiple
posterior
variables
distribution
involved
p(x),
in
the
posterior
distribution
p(x),
Back
to our
case,
the
posterior
distribution
Back
to
our
case,
the
posterior
distribution
p(Z,
p,
q|A)
is
the
main
interest
but
it
is
computationally
intractable.
est
but
it
is
computationally
intractable.
Consider
distribution
called
mean
field
theory
is
to
search
over
some
distribution
clas
called
mean
field
theory
is
to
search
over
some
distribution
class
Q,
such
mean
field
theory
usually
assumes
the
independence
of
these
variables
in
an field
theory
usually
assumes
the
independence
of
these
variables
in
MF
MF
MF
MF
MF
mean
field
theory
usually
assumes
mean
the
field
independence
theory
usually
of
these
assumes
variables
the
in
independence
of
these
variables
in
2.1.
The mean computationally
fieldMF
estimator (π̂ intractable.
, α̂but
α̂computationally
, β̂q )distribution
+
+
est
intractable. Con
p , β̂it
p is,Consider
q q̂MF
+n Theorem
n est but it is +
nq̂ (4)= min KL(q∥p).
= min
KL(q∥p).
(4)
that
the
Kullback-Leibler
divergence
is
minimized.
That
is
divergence
is
minimized.
That
is
q(x)
=a Kullback-Leibler
q
(x
),
a
way
to
reduce
q(x)
=
computational
q
(x
),
complexity.
a
way
to
reduce
The
socomputational
complexity.
q(x)
=
q
(x
),
a
way
to
reduce
computational
complexity.
TheThe
x) = ni=1 qthat
way
to
reduce
computational
complexity.
The
soi
i
i
i
i
i
defined
in
Equation
(8)
is
equivalent
to
i (x
i ),the
i=1
i=1 q ′ q∈Q
q∈Q
′ ,α
′ ′,β
′ ,α′ ,β′′ (Z,
i=1
q
p,soq)so= qπ ′ (
′
′
′
′
′
′
(Z,
p,
q)
=
q
(Z)q
(p)q
π
π ,αp ,βp ,αq ,βq
π
αp ,βpp p q αq ,βq′ (q)
called
mean
field
theory
is
to
search
called
over
mean
some
field
distribution
theory
is
to
class
search
Q,
such
over
some
distribution
class
Q,
such
′
′
′
′
′
′ (Z)qα
q
(Z,
p,
q)
=
q
′ ,α
′ ,βsearch
′ ,α′ ,β
′ (Z,
′distribution
′,α
′,β
′q,β ′ (q)
q
p,
q)
=
q
(Z)q
(p)q
mean
field
theory
is
to
over
some
class
Q,
such
π
,α
,β
π
led mean field theory is called
to search
over
some
distribution
class
Q,
such
π
π
α
,β
α
p
p
q
MF
MF
p p q q
p p
q q
MF to
MFthe
MF
(4)
q̂
=
min
KL(q∥p).
Back
to
our
case,
the
posterior
distribution
p(Z,
p,
q|A)
(4)
q̂
=
min
KL(q∥p).
Back
our
case,
posterior
distribution
p(Z,
p,
q|A)
is
the
main
interthat the Kullback-Leibler
is
the
minimized.
Kullback-Leibler
That
is
divergence
is
minimized.
That
is
(π̂
, α̂pMFdivergence
, β̂pMF
, α̂that
,
β̂
)
=
arg
min
f
(π,
α
,
β
,
α
,
β
;
A),
p p
q q
q
q
that the Kullback-Leibler
divergence
is minimized.
That is
q∈Q
q∈Q
at the Kullback-Leibler divergence
is minimized.
That
is
π∈Π1

(4)

est
but it is Consider
computationally
intractable. Consider distrib
est but it is computationally
intractable.
distribution
MF
MF
α
p ,βp ,αq ,βq >0
q̂

=(4)
min KL(q∥p).

q̂

= min KL(q∥p).

MF
q∈Q distribution
q∈Qis distribution
Back
to our
case,
the KL(q∥p).
posterior
p(Z, p, q|A) is the m
Back to our
case,
theKL(q∥p).
posterior
p,qq|A)
the
main
inter(4)
q̂MF
= p(Z,
min
q̂
= min
′ ,α
′ ,β
′
′ α′′ ,β
′ ,α′ ,β ′ ,α′ ,β ′ (Z, p, q) = qπ ′ (Z)q
′ ,β
′ (p)q
′ (q)
(Z,
p, q) = qπ′ (Z)qα′p ,βp′ (p)qα′q ,β
q
π
π
α
where
q∈Q
p ,αq ,β
p p q q
p pp
qq q
q∈Q
estdistribution
but
computationally
intractable.
est but
intractable.
Consider
distribution
Back it
to is
ourcomputationally
case, the posterior
Backittoisour
p(Z,
case,
p, q|A)
the posterior
is the
main
distribution
interp(Z,Consider
p, q|A) is thedistribution
main inter1
est but
is computationally
intractable.
estp(Z,
but
Consider
itposterior
isq|A)
computationally
distribution
intractable.
Consider
distribution
T main
toβour
case,
the
distribution
p,q )]q|A)
is the main interBack to our case,
theitposterior
distribution
interf (π, αBack
= thA
− λ1p,
1Tn +
λI′isn ,′the
ππ
i′ +
[ψ(α
)′ −
ψ(β
kAk
q′; A)
q p(Z,
1
′(p)q
′ (Z,
′ (Z)q
′ (Z, p, q) =n q
′q
′ ′,β
′ ,β
p,
q)
=
q
qpπ,′β,αp ,′pα,βqp′, ,α
(Z)q
(q)
2
π
,α
,β
,α
,β
π
α′p ,βp′ (p)qα′q ,βq′ (q)
,β
π
α
α
p pp pq q
q q
q q
est
but
computationally
Consider
distribution
′ ,α
′ ,β ′ ,αit
′ ,β is
′ (Z,
′ (Z)q
′ ,β
′ intractable.
′(p)q
′ ,βα
′ ′,α
′ ′,β
′ (Z,
′ (Z)q
′ ,β ′ (p)qα′ ,β ′ (q)
but it is computationally
Consider
distribution
qπintractable.
p,
q)
=
q
q
(q)
p,
q)
=
q
n
π
α
π
,α
,β
π
α
h
i
p p q q
p p
p p q qq q X
p p
q q

qπ′ ,α′p ,βp′ ,α′q ,βq′ (Z, p, q) =

n
pri
[ψ(βq ) − ψ(αq + βq )] −
KL Categorical(πi,· )kCategorical(πi,·
)
2
′ ,α
′ ,β ′ (Z,
′
′
′
′
′
qπ′ ,αα′p ,β
p,
q)
=
q
(Z)q
(p)q
(q)
′ ,β
′ (p)q
′
′
qπ′ (Z)q
(q)
π
α ,β
αq ,βq
i=1
pp p q q αq ,βq
 p p

− KL Beta(αp , βp )kBeta(αppri , βppri ) − KL Beta(αq , βq )kBeta(αqpri , βqpri ) ,
+

9

MEAN FIELD FOR COMMUNITY DETECTION

and
t = [[ψ(αp ) − ψ(βp )] − [ψ(αq ) − ψ(βq )]] /2

(9)

λ = [[ψ(βq ) − ψ(αq + βq )] − [ψ(βp ) − ψ(αp + βp )]] /(2t).

(10)

The explicit formulation in Theorem 2.1 is helpful to understand the
global minimizer of the mean field method. However, the global minimizer
π̂ MF remains computationally infeasible as the objective function is not convex. Fortunately, there is a practically useful algorithm to approximate it.
2.5. Coordinate Ascent Variational Inference. CAVI is possibly the most
popular algorithm to approximate the global minimum of the mean field variational Bayes. It is an iterative algorithm. In Equation (8), there are latent
variables {Zi,· }ni=1 , p, q. CAVI updates them one by one. Since the distribution class of q is uniquely determined by the parameters {πi,· }ni=1 , αp , βp , αq , βq ,
equivalently we are updating those parameters iteratively. Theorem 2.2 gives
explicit formulas for the coordinate updates.
Theorem 2.2. Starts with some π, αp , βp , αq , βq , the CAVI update for
each coordinate (i.e., Equation (3) and Equation (4)) has an explicit expression as follows:
• Update on p:
αp0

=

αppri

+

k
XX

πi,a πj,a Ai,j , and

i<j a=1

βp0

=

βppri

+

k
XX
i<j a=1

πi,a πj,a (1 − Ai,j ).

• Update on q:
XX
XX
αq0 = αqpri +
πi,a πj,b Ai,j , and βq0 = βqpri +
πi,a πj,b (1 − Ai,j ).
i<j a6=b

i<j a6=b

• Update on Zi,· , ∀i = 1, 2, . . . , n:


X
pri
0
πi,a
∝ πi,a
exp 2t
πj,a (Ai,j − λ) , ∀a = 1, 2, . . . , k,
j6=i

where t and λ are defined in EquationP
(9) and Equation (10) respec0 = 1.
tively, and the normalization satisfies ka=1 πi,a

All coordinate updates in Theorem 2.2 have explicit formulas, which
makes CAVI a computationally attractive way to approximate the global
optimum q̂MF for the community detection problem.

10

ZHANG, ZHOU

2.6. Batch Coordinate Ascent Variational Inference. The Batch Coordinate Ascent Variational Inference (BCAVI) is a batch version of CAVI. The
difference lies in that CAVI updates the rows of π sequentially one by one,
0 } according to Thewhile BCAVI uses the value of π to update all rows {πi,·
orem 2.2. This makes BCAVI especially suitable for parallel and distributed
computing, a nice feature for large scale network analysis.
We define a mapping h : Π1 → Π1 as follows. For any π ∈ Π1 , we have


X
pri
(11)
[ht,λ (π)]i,a ∝ πi,a
exp 2t
πja (Ai,j − λ) ,
j6=i

with parameters t and λ. For BCAVI, we update π by π 0 = ht,λ (π) in each
batch iteration, with t, λ defined in Equations (14) and (15). See Algorithm
1 for the detailed implementation of BCAVI algorithm.
Algorithm 1: Batch Coordinate Ascent Variational Inference (BCAVI)

1

Input: Adjacency matrix A, number of communities k, hyperparameters
π pri , αppri , βppri , αqpri , βqpri , initializer π (0) , number of iterations S.
Output: Mean variational Bayes approximation π̂, α̂p , β̂p , α̂q , β̂q .
for s = 1, 2, . . . , S do
(s)
(s)
(s)
(s)
Update αp , βp , αq , βq by
(12)
αp(s) = αppri +

k X
X

(s−1)

Ai,j πi,a

(s−1)

πj,a

, βp(s) = βppri +

a=1 i<j

k X
X
(s−1) (s−1)
(1 − Ai,j )πi,a πj,a ,
a=1 i<j

(13)
αq(s) = αqpri +

XX

(s−1)

Ai,j πi,a

(s−1)

πj,b

, βq(s) = βqpri +

a6=b i<j

a6=b i<j
2

XX
(s−1) (s−1)
(1 − Ai,j )πi,a πj,b .

Define
i h
ii
1 hh
ψ(αp(s) ) − ψ(βp(s) ) − ψ(αq(s) ) − ψ(βq(s) )
2
i h
ii
1 hh
= (s) ψ(βq(s) ) − ψ(αq(s) + βq(s) ) − ψ(βp(s) ) − ψ(αp(s) + βp(s) ) ,
2t

(14)

t(s) =

(15)

λ(s)

where ψ(·) is the digamma function. Then update π (s) with
π (s) = ht(s) ,λ(s) (π (s−1) ),

3

where the mapping h(·) is defined as in Equation (11).
end
(S)
(S)
(S)
(S)
We have π̂ = π (S) , α̂p = αp , β̂p = βp , α̂q = αq , β̂q = βq .

11

MEAN FIELD FOR COMMUNITY DETECTION

Remark 2.1. The definitions of t(s) and λ(s) in Equations (14) and (15)
involve the digamma function, which costs a non-negligible computational
resources each time called. Note that we have ψ(x) ∈ (log(x− 21 ), log x) for all
x > 1/2. For the computational purpose, we propose to use the logarithmic
function instead of digamma function in Algorithm 1, i.e., Equations (14)
and (15) are replaced by
(s)

(s) (s)

(16)

t(s) =

1
αp βq
log (s) (s) ,
2
β p αq
(s)

(s)

and λ(s) =
(s)

(s)

(s)

βq (αp + βp )
1
log (s)
.
(s)
(s) (s)
2t
(αq + βq )βp

(s)

Later we show that αp , βp , αq , βq are all at least in the order of np,
which goes to infinity, and thus the error caused by using the logarithmic
function to replace the digamma function is negligible. All theoretical guarantees obtained in Section 3 for Algorithm 1 (i.e., Theorem 3.1, Theorem
3.2) still hold if we use Equation (16) to replace Equations (14) and (15).
3. Theoretical Justifications. In this section, we establish theoretical justifications for BCAVI for community detection under the Stochastic
Block Model. Though Z, p and q are all unknown, the main interest of community detection is on the recovery of the assignment matrix Z, while p
and q are nuisance parameters. As a result, our main focus is on developing
convergence rate of BCAVI for π.
3.1. Loss Function. We use `1 norm to measure the performance of recovering Z. Let Φ be the set of all the bijections from [k] to [k]. Then for
any Z, Z ∗ ∈ Π1 , the loss function is defined as
X
∗
`(Z, Z ∗ ) = inf kZ − φ ◦ Z ∗ k1 = inf
(17)
|Zi,a − Zi,φ(a)
|.
φ∈Φ

φ∈Φ

i,a

Note that the infimum over Φ addresses the issue of identifiability over the
labels. For instance, in the case of n = 4, k = 2, the assignment vector
z = (1, 1, 2, 2) and z 0 = (2, 2, 1, 1) give the same partition. In Equation (17)
two equivalent assignments give the same loss.
There are a few reasons for the choise of the `1 norm. When both Z, Z 0 ∈
Π0 , the `1 distance between Z and Z 0 is equal to the `0 norm, i.e., the Hamming distance between the corresponding assignment vectors r−1 (Z) and
r−1 (Z 0 ), which is the default metric used in community detection literature
[12, 35]. The other reason is related to the interpretation of Π1 . Since each
row of Π1 corresponds to a categorical distribution, it is natural to use the
`1 norm, the total variation distance, to measure their diffidence.

12

ZHANG, ZHOU

3.2. Ground Truth. We use the superscript asterisk (∗ ) to indicate the
ground truth. The ground truth of connectivity matrix B ∗ is
B ∗ = q ∗ 1k 1Tk + (p∗ − q ∗ )Ik ,
where p∗ is the within community connection probability and q ∗ is the between community connection probability. Throughout the paper, we assume
p∗ > q ∗ such that the network satisfies the so-called “assortative” property, with the within-community connectivity probability larger than the
between-community connectivity probability.
We further assume the network is generated by the true assignment matrix
Z ∗ in the sense that Pi,j = (Z ∗ B ∗ Z ∗T )i,j for all i 6= j. We are interested
in deriving a statistical guarantee of `(π̂ (s) , Z ∗ ). Throughout this section we
(ρ,ρ0 )
(ρ,ρ0 )
consider cases Z ∗ ∈ Π0 or Z ∗ ∈ Π0
, where Π0
is defined to be a subset
of Π0 with all the community sizes bounded between ρn/k and ρ0 n/k. That
is,
(ρ,ρ0 )

Π0

= {π ∈ Π0 : ρn/k ≤ |{i ∈ [n] : πi,a = 1}| ≤ ρ0 n/k, ∀a ∈ [k]}.

It is worth mentioning that ρ, ρ0 are not necessarily constants. We allow the
community sizes not to be of the same order in the theoretical analysis.
3.3. Theoretical Justifications for BCAVI. In Theorem 3.1, we present
theoretic guarantees of the convergence rate of BCAVI when initialized properly. Define
pri
pri
, and n̄min = min[na + nb ]/2.
w = max max πi,a
/πi,b
i∈[n] a,b∈[k]

a6=b

When w = 1, the priors for {r−1 (Zi,· )}ni=1 are i.i.d. Categorical(1/k, 1/k, . . . , 1/k)
and n̄min = n/2 when there exist only two communities. The following quantity I plays a key role in the minimax theory [35]
hp
i
p
I = −2 log
p∗ q ∗ + (1 − p∗ )(1 − q ∗ ) ,

which is the Rényi divergence of order 1/2 between two Bernoulli distributions: Ber(p∗ ) and Ber(q ∗ ). The proof of Theorem 3.1 is deferred to Section
5.3.
Theorem 3.1. Let Z ∗ ∈ Π0 . Let 0 < c0 < 1 be any constant. Assume
0 < c0 p∗ < q ∗ < p∗ = on (1),
(18) nI/[wk[n/n̄min ]2 ] → ∞, and αppri , βppri , αqpri , βqpri = on ((p∗ − q ∗ )n2 /k).

MEAN FIELD FOR COMMUNITY DETECTION

13

Under the assumption that the initializer π (0) satisfies `(π (0) , Z ∗ ) ≤ cinit n̄min
for some sufficiently small constant cinit with probability at least 1 − , there
exist some constant c > 0 and some η = on (1) such that in each iteration
for the BCAVI algorithm, we have
`(π (s) , Z ∗ )
`(π (s+1) , Z ∗ ) ≤ n exp(−(1 − η)n̄min I) + p
, ∀s ≥ 0,
nI/[wk[n/n̄min ]2 ]
1

holds uniformly with probability at least 1 − exp[−(n̄min I) 2 ] − n−c − .

Theorem 3.1 establishes a linear convergence rate for BCAVI algorithm.
The coefficient [nI/[wk[n/n̄min ]2 ]]−1/2 is independence of s, and goes to
0 when n grows. The following theorem is an immediate consequence of
Theorem 3.1.
Theorem 3.2. Under the same condition as in Theorem 3.1, for any
s ≥ s0 , [nI/k]/ log[nI/[wk[n/n̄min ]2 ]], we have
(
n exp(−(1 − o(1))ρnI/k), k ≥ 3;
(s)
∗
`(π̂ , Z ) ≤ n exp(−(1 − 2η)n̄min I) ≤
n exp(−(1 − o(1))nI/2), k = 2,
1

with probability at least 1 − exp[−(n̄min I) 2 ] − n−c − .
Theorem 3.2 shows that BCAVI provably attains the statistical optimality from the minimax lower bound in Theorem 3.3 after at most s0 iterations. When the network is sparse, i.e., p∗ and q ∗ are at most in an order
of (log n)/n, the quantity s0 can be shown to be o(log n), and then BCAVI
converges to be minimax rate within log n iterations. When the network is
dense, i.e., p∗ and q ∗ are far bigger than (log n)/n, log n iterations are not
enough to attain the minimax rate. However, `(π (s) , Z ∗ ) = o(n−a ) for any
a > 0 when s ≥ log n, and thus all the nodes can be correctly clustered
with high probability by clustering each note to a community with the highest assignment probability. Therefore, it is enough to pick the number of
iterations to be log n in implementing BCAVI.
Under the assumption nI/(k log k) → ∞, we have
(
n exp(−(1 − o(1))ρnI/k), k ≥ 3;
inf sup E`(π̂, Z ∗ ) ≥
.
π̂
n exp(−(1 − o(1))nI/2), k = 2,
(ρ,ρ0 )
Z ∗ ∈Π

Theorem 3.3.

0

14

ZHANG, ZHOU

Theorem 3.3 gives the minimax lower bound for community detection
problems with respect to the `(·, ·) loss. In Theorem 3.2, under the addi(ρ,ρ0 )
tional assumption that Z ∗ ∈ Π0
, it immediately reveals that BCAVI
converges to the minimax rate after s0 iterations. As a consequence, BCAVI
is not only computationally efficient, but also achieves statistical optimality.
The minimax lower bound in Theorem 3.3 is almost identical to the minimaxity established in [35]. The only difference is that [35] consider a `0 loss
function. The proof of Theorem 3.3 is just a routine extension of that in
[35]. Therefore, we omit the proof.
To help understand Theorem 3.1, we add a remark on conditions on model
parameters and priors, and a remark on initialization.
Remark 1 (Conditions on model parameters and priors). The community
sizes are not necessarily of the same order in Theorem 3.1. If we further
pri
assume ρ, ρ0 are constants, and the prior πi,a
 1/k, ∀i ∈ [n], a ∈ [k] (for
example, uniform prior), and then the first condition in Equation (18) is
equivalent to
nI/k 3 → ∞,
noting that n/n̄min  k and w  1. This condition is necessary for consistent
community detection [35] when k is finite. The assumptions in Equation
(18) is slightly stronger than the assumption in [23], which is essentially
nI ≥ Ck 2 log k for a sufficient large constant C.
Under the assumption nI/k 3 → ∞, since we have I  (p∗ − q ∗ )2 /p∗ , it
can be shown that p∗ , q ∗ are far bigger than n−1 , and then the second part
of Equation (18) can also be easily satisfied. For instance, we can simply set
αppri , βppri , αqpri , βqpri all equals to 1, i.e., consider non-informative priors.
Remark 2 (Initialization). The requirement on the initializers for BCAVI
in Theorem 3.1 is relatively weak. When k is a constant and the community sizes are of the same order, the condition needed is `(π (0) , Z ∗ ) ≤ cn
for some small constant c. Many existing methodologies in community detection literature can be used. One popular choice is spectral clustering.
Established in [9, 12, 21], the spectral clustering has a mis-clustering error bound as O(k 2 /I). From Equation (18), the error is o(n̄min ), and then
the condition that Theorem 3.1 requires for initialization is satisfied. The
semidefinite programming (SDP), another popular method for community
detection, also enjoys satisfactory theoretical guarantees [10, 16], and is suitable as an initializer.
4. Discussion.

MEAN FIELD FOR COMMUNITY DETECTION

15

4.1. Statistical Guarantee of Global Minimizer. Though it is often challenging to obtain the global minimizer of the mean field method, it is still interesting to understand the statistical property of the global minimizer π̂ MF .
Assume that both p∗ and q ∗ are known, the optimization problem stated in
Theorem 2.1 can be further simplified. The posterior
distribution becomes
Q
p(Z|A). We use a product measure qπ (Z) = i qi (πi,· ) for approximation,
and then π̂ MF = arg minπ∈Π1 KL[qπ (Z)kp(Z|A)]. Theorem 4.1 reveals that
π̂ MF is rate-optimal, not surprisingly given the theoretical results obtained
for BCAVI, an approximation of π̂ MF .
Theorem 4.1. Assume p∗ and q ∗ are known. Under the assumption
ρnI/[wk 2 [n/n̄min ]2 ] → ∞, there exist some constant c > 0 and η = on (1)
such that
`(π̂ MF , Z ∗ ) ≤ n exp(−(1 − η)n̄min I)
1

with probability at least 1 − exp[−(n̄min I) 2 ] − n−c .
4.2. Gibbs Sampling. In Section 3.3 we analyze an iterative algorithm,
BCAVI, and establish its linear convergence towards statistical optimality.
The framework and methodology we establish is not limited to BCAVI, but
can be extended to other iterative algorithms, including Gibbs sampling.
As a popular Markov chain Monte Carlo (MCMC) algorithm, Gibbs sampling has been widely used in practice to approximate the posterior distribution. There is a strong tie between Gibbs sampling and the mean field
variational inference: both implement coordinate updates using conditional
distributions. Using the general notation introduced in Section 2.1, to approximate p(x|y), Gibbs sampling obtains the update on xi by a random generation from the conditional distribution p(xi |x−i ,y), while the variational

inference updates in a deterministic way with exp Eq−i log p(xi |x−i , y) .
We present a batched version of Gibbs sampling for community detection.
It involves iterative updates with
• Generate p(s) by sampling from p(p|q (s−1) , Z (s−1) , A);
• Generate q (s) by sampling from p(q|p(s−1) , Z (s−1) , A);
(s)
(s−1) (s−1)
• Generate Zi,· independently by sampling from p(Zi,· |Z−i,· , p(s) , q (s) , A),
for i ∈ [n].
We include the detailed implementation as Algorithm 2 in the supplemental
material (Section A.1). The similarity between Algorithm 1 and Algorithm 2
makes it possible for us to analyze the output of Gibbs sampling in a similar
way as we did for the variational inference.

16

ZHANG, ZHOU

Theorem 4.2. Assume the initializer Z (0) satisfies `(Z (0) , Z ∗ ) ≤ cinit n̄min
for some sufficiently small constant cinit with probability at least 1−. Under
the same condition as in Theorem 3.1, there exist some constant c > 0 and
some η, η 0 = on (1) that go to 0 slowly, such that for all s ≥ 0 of the batched
Gibbs sampling (Algorithm 2), we have
h
i
EZ (s+1) `(Z (s+1) , Z ∗ ) A, Z (0) ≤ n exp(−(1 − η)n̄min I) + csn `(Z (0) , Z ∗ ) + (s + 1)nbn
1

−c
holds
probability
1 − exp[−(n̄min I)
p2 )] − n − , where bn =
 with

at least
02
2
02
2
exp −η n̄min + exp −η n I and cn = 1/ nI/[wk[n/n̄min ]2 ]. Consequently, for s = [nI/k]/ log[nI/[wk[n/n̄min ]2 ]], we have
h
i
EZ (s+1) `(Z (s+1) , Z ∗ ) A, Z (0) ≤ exp(−(1 − 2η)n̄min I),
1

with probability at least 1 − exp[−(n̄min I) 2 ] − n−c − .
Theorem 4.2 establishes theoretical justification for batched Gibbs sampling for community detection. Despite that we have the same cn and similar
convergence as Theorem 3.1, some extra efforts are needed due to the existence of randomness in each iterative update. The additional term of bn is
necessary to handle the extreme events due to random generation. Note that
(s + 1)nbn is dominated by n exp(−(1 − η)n̄min I) as long as s ≤ en . Thus,
when s ≤ en , we have similar “linear convergence” results as in Theorem
3.1.
4.3. An Iterative Algorithm for Maximum Likelihood Estimation. Maximum likelihood estimator (MLE) usually yields statistical optimality. However, the maximization of the likelihood p(A|Z, p, q) over Z, p, q is computationally infeasible. Inspired by the procedures proposed in Algorithm 1 and
Algorithm 2, we may approach max p(A|Z, p, q) by alternating maximization. We use a batched coordinate maximization:
• Maximize p(A|p, q (s−1) , Z (s−1) ) over p to obtain p(s) ;
• Maximize p(A|p(s−1) , q, Z (s−1) ) over q to obtain q (s) ;
(s−1)
(s)
• Maximize p(A|p(s−1) , q (s−1) , Zi,· , Z−i,· ) over Zi,· to obtain Zi,· , for
each i ∈ [n].
We include its detailed implementation in Algorithm 3 in the supplemental
material (Section A.2). We have the following theoretical guarantee of this
iterative algorithm to approximate the MLE.

17

MEAN FIELD FOR COMMUNITY DETECTION

Theorem 4.3. Assume the initializer Z (0) satisfies `(Z (0) , Z ∗ ) ≤ cinit n̄min
for some sufficiently small constant cinit with probability at least 1−. Under
the same condition as in Theorem 3.1, there exist some constant c > 0 and
some η = on (1), such that in each iteration of the BCAVI algorithm,
`(Z (s) , Z ∗ )
`(Z (s+1) , Z ∗ ) ≤ n exp(−(1 − η)n̄min I) + p
, ∀s ≥ 0,
nI/[wk[n/n̄min ]2 ]
1

holds with probability at least 1 − exp[−(n̄min I) 2 ] − n−c − .

Algorithm 3 is essentially the same with the procedure proposed in [12].
However, [12] can only analyze the performance of one single iteration from
Z (0) (i.e., `(Z (1) , Z ∗ )), and it requires extra data splitting steps. Theorem
4.3 provides a stronger and cleaner result compared with that of [12].
5. Proofs of Main Theorems. In this section, we give proofs of the
theorems in Section 2 and Section 3. We first present the proof of Theorem
2.1 in Section 5.1. Then we give the proof Theorem 2.2 in Section 5.2. The
proof of Theorem 3.1 is given in Section 5.3.
5.1. Proof of Theorem 2.1. From Equation (8), by some algebra (see
Equation (53) in Appendix D for detailed derivation) we have
(19)
(π̂ MF , α̂pMF , β̂pMF , α̂qMF , β̂qMF ) =

arg min
π∈Π1
αp ,βp ,αq ,βq >0

Eq [log p(A|Z, p, q)] − KL(q(Z, p, q)kp(Z, p, q)),

where we use q instead of qπ,αp ,βp ,αq ,βq for simplicity. From the conditional
distribution in Equation (6), the log-likelihood function can be simplified as


XX
Bab
log p(A|Z, p, q) =
Zia Zjb Ai,j log
+ log(1 − Bab ) .
1 − Bab
a,b i<j

Due to the independence of Z and p, q under q, we have



XX


Eq [log p(A|Z, p, q)] = Eq(p,q) Eq(Z)
Zi,a Zj,b Ai,j log
a,b i<j



XX
= Eq(p,q) 
πi,a πj,b Ai,j log
a,b i<j

Bab
1 − Bab



+ log(1 − Bab ) 



Bab
+ log(1 − Bab )  .
1 − Bab

18

ZHANG, ZHOU

Since Ba,a = p, ∀a ∈ [k] and Ba,b = q, ∀a 6= b, we have
(20)





XX
1−p 
p(1 − q)
+ log
Eq [log p(A|Z, p, q)] = Eq(p,q) 
πi,a πj,a Ai,j log
q(1 − p)
1−q
a i<j


i
h
XX
q
+ Eq(p,q) 
+ log(1 − q)  .
πi,a πj,b Ai,j log
1−q
a,b i<j

By properties of Beta distribution, we obtain
Eq(p,q) log

p(1 − q)
= Eq(p) [log p − log(1 − p)] − Eq(q) [log q − log(1 − q)]
q(1 − p)
= [ψ(αp ) − ψ(βp )] − [ψ(αq ) − ψ(βq )] ,

and
Eq(p,q) log

1−q
= Eq(q) log(1 − q) − Eq(p) log(1 − p)
1−p
= [ψ(βq ) − ψ(αq + βq )] − [ψ(βp ) − ψ(αp + βp )] .

This leads to
(21)







X
X
XX
1−p 
p(1 − q)
+ log
= 2t 
Eq(p,q) 
πi,a πj,a (Ai,j − λ)
πi,a πj,a Ai,j log
q(1
−
p)
1
−
q
a
a
i<j

i<j

= thA − λ1n 1Tn + λIn , ππ T i.

Similarly we can obtain
(22)


h
XX
Eq(p,q) 
πi,a πj,b Ai,j log


a,b i<j

q
1−q


i
+ log(1 − q) 

X
X

XX
q
= Eq(q) log
Ai,j
πi,a πj,b + Eq(q) log(1 − q)
πi,a πj,b
1−q
i<j

a,b

1
n
= [ψ(αq ) − ψ(βq )] kAk1 + [ψ(βq ) − ψ(αq + βq )] ,
2
2

i<j a,b

19

MEAN FIELD FOR COMMUNITY DETECTION

where we use the fact that kπi,· k1 = 1, ∀i ∈ [n]. Now consider the KullbackLeibler divergence between q(Z, p, q) and p(Z, p, q). Due to the independence
of p, q and {Zi,· }ni=1 in both distributions, we have
(23)
KL(q(Z, p, q)kp(Z, p, q)) = KL(q(Z)kp(Z)) + KL(q(p)kp(p)) + KL(q(q)kp(q))
n
h
i
X
pri
=
KL Categorical(πi,· )kCategorical(πi,·
)
i=1





+ KL Beta(αp , βp )kBeta(αppri , βppri ) + KL Beta(αq , βq )kBeta(αqpri , βqpri ) .

By Equations (19) - (23), we conclude with the desired result.
5.2. Proof of Theorem 2.2. Note that


" k
#
X
X
Bzi ,zj =
Zi,a Zj,a p + 
Zi,a Zj,b  q.
a=1

a6=b

We rewrite the joint distribution p(p, q, z, A) in Equation (7) as follows,
(24)
p(p, q, Z, A)


" n
#
Pk
Pk
Y pri
Y
Y



Z
Z
Z
Z
i,a j,b
i,a j,a


q Ai,j (1 − q)1−Ai,j a6=b
=
πi,zi 
pAi,j (1 − p)1−Ai,j a=1
i=1

×

"

i<j

Γ(αppri + βppri )

Γ(αppri )Γ(βppri )

p

αpri
p −1

(1 − p)

βppri −1

#"

i<j

Γ(αqpri + βqpri )

Γ(αqpri )Γ(βqpri )

q

αpri
q −1

(1 − q)

βqpri −1

#

.

Updates on p and q. From Equation (24), p has conditional probability as

"
#
P
pri
pri
Y
k

pri
pri
Γ(α
+
β
)
Z
Z
p
p
i,a
j,a

pαp −1 (1 − p)βp −1 .
p(p|q, Z, A) ∝ 
pAi,j (1 − p)1−Ai,j a=1
pri
pri
Γ(α
)Γ(β
)
p
p
i<j

Then the CAVI update in Equation (4) leads to


q̂(p) ∝ exp Eq(q,Z) log p(p|q, Z, A)

"
#
k
pri
pri
XX
 A

pri
pri
Γ(α
+
β
)
p
p
∝ exp Eq(Z)
Zi,a Zj,a log p i,j (1 − p)1−Ai,j 
pαp −1 (1 − p)βp −1
pri
pri
Γ(α
)Γ(β
)
p
p
i<j a=1

"
#
k
pri
pri
XX
 A

pri
pri
Γ(α
+
β
)
p
p
= exp 
πi,a πj,a log p i,j (1 − p)1−Ai,j 
pαp −1 (1 − p)βp −1 .
pri
pri
Γ(α
)Γ(β
)
p
p
i<j a=1

20

ZHANG, ZHOU

It can be written as
"
#
h P Pk
i Γ(αpri + β pri ) pri
P
Pk
pri
p
p
q̂(p) ∝ p i<j a=1 πi,a πj,a Ai,j (1 − p) i<j a=1 πi,a πj,a (1−Ai,j )
pαp −1 (1 − p)βp −1 .
Γ(αppri )Γ(βppri )
The distribution of p is still Beta p ∼ Beta(αp0 , βp0 ), with
αp0

=

αppri

+

k
XX

πi,a πj,a Ai,j , and

βp0

i<j a=1

=

βppri

+

k
XX
i<j a=1

πi,a πj,a (1 − Ai,j ).

Similar analysis on q yields updates on αq0 and βq0 . Hence, its proof is omitted.
Updates on {Zi,· }ni=1 .
Zi,· is

From Equation (24), the conditional distribution on

pri
p(Zi,· |Z−i,· , p, q, A) ∝ πi,z
i



Y Ai,j
 Bzi ,zj (1 − Bzi ,zj )1−Ai,j  .
j6=i

Consequently, up to a constant not depending on i, we have
log P(Zi,a = 1|Z−i,· , p, q, A)


X
pri
= log πi,a + log 
Zj,a Ai,j log


 XX


q
p
+ log(1 − p) +
Zj,b Ai,j log
+ log(1 − q) 
1−p
1−q
j6=i b6=a
j6=i






X
X
1
−
q
q
p(1
−
q)
pri
= log πi,a
− log
+
Ai,j log
+ log(1 − q)  .
+ log 
Zj,a Ai,j log
q(1 − p)
1−p
1−q
j6=i

j6=i

Then the CAVI update from Equation (4) leads to

0
πi,a
= q̂Zi,· (Zi,a = 1)


∝ exp Eq(p,q,z−i ) log P(Zi,a = 1|Z−i,· , p, q, A)
h
i
= exp Eq(p) Eq(q) Eq(Z−i,· ) log P(Zi, = 1|Z−i,· , p, q, A)




X
1
−
q
p(1
−
q)
pri
,
(25)
− log
∝ πi,a
exp Eq(p) Eq(q)
πj,a Ai,j log
q(1 − p)
1−p
j6=i

where we use the property that p, q, Z are all independent of each other
under q. Recall that p ∼ Beta(αp , βp ) and q ∼ Beta(αq , βq ). It can be shown
that
p
Eq(p) log
= ψ(αp ) − ψ(βp ), and Eq(p) log(1 − p) = ψ(βp ) − ψ(αp + βp ),
1−p

MEAN FIELD FOR COMMUNITY DETECTION

21

where ψ(·) is digamma function. Similar results hold for Eq(q) log(q/(1 − q))
and Eq(q) log(1 − q). Plug in these expectations to Equation (25), we have


X
pri
0
πi,a
∝ πi,a
exp 2t
πj,a (Ai,j − λ) .
j6=i

5.3. Proof of Theorem 3.1. Theorem 3.1 gives a theoretical justification
for all iterations in the BCAVI algorithm. Due to the limit of pages, in this
section we assume `(π (0) , Z ∗ ) = o(n̄min ). The proof of the case `(π (0) , Z ∗ )
in a constant order of n̄min is essentially the same with slight modification,
and we defer it to Section B.1 in the supplemental material.
To prove the theorem, it is sufficient if we are able to show the loss `(·, Z ∗ )
decreases in a desired way for one BCAVI iteration, when the community
assignment is in an appropriate neighborhood of the truth. Let γ = o(1) be
any sequence that goes to zero when n grows. Define t∗ and λ∗ as the true
counterparts of t and λ, by
t∗ =

1
p∗ (1 − q ∗ )
1
1 − q∗
∗
log ∗
,
and
λ
=
log
.
2
q (1 − p∗ )
2t∗
1 − p∗

The proof of Theorem 3.1 involves three parts as follows.

Part One: One Iteration. Consider any π ∈ Π1 such that kπ − Z ∗ k1 ≤
γ n̄min . Let η 0 be any sequence such that η 0 = o(1). Consider any t and λ
with |t−t∗ | ≤ η 0 (p∗ −q ∗ )/p∗ and |λ−λ∗ | ≤ η 0 (p∗ −q ∗ ). We define F to be the
event, that after applying the mapping ht,λ (·), there exists some η = o(1)
such that
kπ − Z ∗ k1
kht,λ (π) − Z ∗ k1 ≤ n exp(−(1 − η)n̄min I) + p
,
nI/[wk[n/n̄min ]2 ]

holds uniformly over all the eligible π, t and λ. We have
1

P(F) ≥ 1 − exp[−(n̄min I) 2 )] − n−r ,
for some constant r > 0. We defer its proof to the later part of this section.
Part Two: Consistency of Model Parameters. Consider any π ∈ Π1
such that kπ − Z ∗ k1 ≤ γ n̄min . Define
(26)
αp = αppri +

k X
X
a=1 i<j

Ai,j πi,a πj,a ,

βp = βppri +

k X
X
(1 − Ai,j )πi,a πj,a ,
a=1 i<j

22

ZHANG, ZHOU

and
(27)
αq = αqpri +

XX

Ai,j πi,a πj,b ,

a6=b i<j

βq = βqpri +

XX
(1 − Ai,j )πi,a πj,b ,
a6=b i<j

and consequently,
(28)
(29)

1
[[ψ(αp ) − ψ(βp )] − [ψ(αq ) − ψ(βq )]]
2
1
λ=
[[ψ(βq ) − ψ(αq + βq )] − [ψ(βp ) − ψ(αp + βp )]] .
2t
t=

From Lemma C.1, we have a concentration of t, λ towards t∗ , λ∗ . That is,
there exists some η 0 = o(1), such that with probability at least 1 − e3 5−n ,
the following inequalities hold
|t − t∗ | ≤ η 0 (p∗ − q ∗ )/p∗ , and |λ − λ∗ | ≤ η 0 (p∗ − q ∗ ),
uniformly over all the eligible π.
Part Three: Multiple Iterations. Consider any π ∈ Π1 such that kπ − Z ∗ k1 ≤
γ n̄min . Define αp , βp , αq , βq , t, λ as Equations (26) - (29). A combination of
results from Part One and Part Two immediately implies that
(30)

kπ − Z ∗ k1
,
kht,λ (π) − Z ∗ k1 ≤ n exp(−(1 − η)n̄min I) + p
nI/[wk[n/n̄min ]2 ]

1

holds uniformly over all the eligible π with probability at least 1−exp[−(n̄min I) 2 )]−
n−r . This is sufficient to show Theorem 3.1.
The only thing left to be proved, the most critical part towards the proof
of Theorem 3.1, is the claim we made in Part One. We are going to prove
the claim as follow.
Proof Sketch of Part One. The error associated with the [ht,λ (π)]i,· is a
function of π and Ai,· . It can be decomposed into a summation of two terms,
one only involves the ground truth Z ∗ and the other involves the deviation
π − Z ∗ . That is,
∗
[ht,λ (π)]i,· − Zi,·

1

≤ fi,1 (Z ∗ , Ai,· ) + fi,2 (π − Z ∗ , Ai,· ).

23

MEAN FIELD FOR COMMUNITY DETECTION

Consequently,
(31)

∗

kht,λ (π) − Z k1 ≤

n
X

|i=1

∗

fi,1 (Z , Ai,· ) +
{z

involves

Z∗

}

n
X

|i=1

fi,2 (π − Z ∗ , Ai,· ) .
{z

involves π−Z ∗

}

With a proper choice of f·,1 and f·,2 , the first term on the RHS of Equation
(31) leads to the minimax rate n exp(−(1 − η)n̄min I). Up to a constant not
dependent on π, Z ∗ or A, the second term can be written as
n
X
i=1

fi,2 (π − Z ∗ , Ai,· ) .

X
a

∗ T
∗
(π·,a − Z·,a
) (A − EA)(A − EA)T (π·,a − Z·,a
).

In this way it is all about the random
EA and
P there exist∗ sharp
P matrix A −
∗ 2 ≤
bounds on kA − EAkop . Note that a π·,a − Z·,a
a π·,a − Z·,a 1 ≤
∗
kπ − Z k1 . The second term ends up being upper bounded by kπ − π ∗ k1
multiplied by a coefficient factor.
Proof of Part One. Denote z = r−1 (Z ∗ ). By the definition of ht,λ (·) in
Equation (11), we have
h P
i
P
pri
2 a6=zi πi,a
exp 2t j6=i πj,a (Ai,j − λ)
∗
i
h P
≤ P
[ht,λ (π)]i,· − Zi,·
1
pri
π
π
(A
−
λ)
exp
2t
j,a
i,j
a i,a
j6=i


X
X
≤ 2w
1 ∧ exp 2t
(πj,a − πj,zi )(Ai,j − λ) .
a6=zi

j6=i

Define f (x) = 1 ∧ exp(−x). It can be shown P
that for any x0 < 0 and any
integer m ≥ 1 we have f (x) ≤ exp(x0 ) + m−1
l=0 exp(lx0 /m)I{x ≥ (l +
1)x0 /m}, which can be seen as a stepwise approximation of theP
continuous
function f (x). By taking x0 = −(na +nzi )I/2 and letting x = 2t j6=i (πj,a −
πj,zi )(Ai,j − λ), we have
"




m−1
X
X
(n
+
n
)I
l(na + nzi )I
a
zi
∗
[ht,λ (π)]i,· − Zi,· 1 ≤ 2w
exp −
+ 2w
exp −
2
2m
a6=zi
l=0
#
X  X
(l + 1)(na + nzi )I
×
I 2t
(πj,a − πj,zi )(Ai,j − λ) ≥ −
.
2m
a6=zi

j6=i

We choose some m → ∞ slowly such that
(32)

m = o(n̄min I) and m = o([wnI/[k[n/n̄min ]2 ]1/4 ).

24

ZHANG, ZHOU

Thus, we have
∗

kht,λ (π) − Z k1 ≤ 2wnk exp(−n̄min I) + 2w

m−1
k X
XX
l=0 a=1 b6=a

"



l(na + nb )I
exp −
2m



#
X X
(l + 1)(na + nb )I
×
I
(πj,a − πj,b )(Ai,j − λ) ≥ −
4mt

(33)

i:zi =b

j6=i

where we use the fact that mina6=b (na + bb )/2 ≥ n̄min .
The key to the rest of the analysis is to understand
Equation (33) through
P
the decomposition of the critical quantity j6=i (πj,a − πj,b )(Ai,j − λ). We
will show for any pair of a, b ∈ [k] such that a 6= b, and any i ∈ [n] such
that zi = b, it is equal to a summation of two terms: one only involves the
ground truth Z ∗ , and the other involves the deviation π − Z ∗ . The former
remains steady along iterations and contributes to the minimax rate, while
the latter needs to be connected with the error kπ − Z ∗ k1 .
∗ + Z∗ −
Let θa,b be a vector of length n such that [θa,b ]j = πj,a − Zj,a
j,b
πj,b , ∀j ∈ [n]. Then we have
(34)
X
X
X
∗
∗
∗
∗
− πj,b )(Ai,j − λ)
)(Ai,j − λ) +
(πj,a − Zj,a
+ Zj,b
(πj,a − πj,b )(Ai,j − λ) =
(Zj,a
− Zj,b
j6=i

j6=i

j6=i

=

X
j6=i

=

X
j6=i

|

∗
∗
)(Ai,j
(Zj,a
− Zj,b
∗
(Zj,a

−

∗
Zj,b
)(Ai,j

X
− λ) +
(Ai,j − λ)[θa,b ]j
j6=i

− λ) + (Ai,· − EAi,· )θa,b +

{z

}

involves Z ∗

|

With the help of Equation (34), Equation (33) can be written as
kht,λ (π) − Z ∗ k1
≤ 2wnk exp(−n̄min I) + 2w

m−1
k X
XX
l=0 a=1 b6=a

"



l(na + nb )I
exp −
2m

X
(EAi,j − λ)[θa,b ]j .
j6=i

{z

involves π−Z ∗

#
X X
X
(l
+
3/2)(n
+
n
)I
a
b
∗
∗
×
I
(Zj,a
− Zj,b
)(Ai,j − λ) ≥ −
−
(EAi,j − λ)[θa,b ]j
4mt
i:zi =b
j6=i
j6=i
"" m−1

#
#
k X
X
X
X 
l(na + nb )I
n̄min I
+ 2w
exp −
×
I (Ai,· − EAi,· )θa,b ≥
.
2m
4mt
a=1 b6=a

l=0

i:zi =b

}

MEAN FIELD FOR COMMUNITY DETECTION

Equations (18) and (32) imply
we have

Pm−1
l=0

exp [−l(na + nb )I/(2m)] ≤ 2. Thus,

kht,λ (π) − Z ∗ k1 ≤ 2wnk exp(−n̄min I) + 2wLsum
+ 4wLsum
,
| {z1 }
| {z2 }
involves Z ∗

where

25

Lsum
1

,

m−1
k X
XX
l=0 a=1 b6=a



l(na + nb )I
exp −
2m

 X

involves π−Z ∗

L1,i (a, b, l),

i:zi =b

P
∗
∗
with
L
(a,
b,
l)
,
I[
1,i
j6=i (Zj,a −Zj,b )(Ai,j −λ) ≥ −(l+3/2)(na +nb )I/(4mt)−
P
j6=i (EAi,j − λ)[θa,b ]j ], and
Lsum
2

,

k X X
X

a=1 b6=a i:zi =b



I (Ai,· − EAi,· )θa,b


n̄min I
≥
.
4mt

In this way we turn kht,λ (π) − Z ∗ k1 into calculations on Lsum
and Lsum
1
2 ,
where the former only involves the ground truth Z ∗ and the latter only
involves the deviation π − Z ∗ .
and Lsum
as follows. Their proofs
We can obtain upper bounds on Lsum
1
2
are deferred to the end of this section.
00
• For Lsum
1 , there exists a sequence η = o(1) such that with probability
1
at least 1 − exp[−2(n̄min I) 2 ], we have


(35)
Lsum
≤ nmk exp −(1 − 2η 00 )n̄min I .
1

• For Lsum
2 , there exist constants c and r such that with probability at
least 1 − n−r − exp(−5np∗ ), we have
(36)

Lsum
≤
2

cknp∗ kπ − Z ∗ k1 cn2 kp∗ exp(−5np∗ )
+
.
(n̄min I/(mt∗ ))2
n̄min I/(mt∗ )

Thus, we have


kht,λ (π) − Z ∗ k1 ≤ 2wnk exp(−n̄min I) + 2wnmk exp −(1 − 2η 00 )n̄min I
+

4cwknp∗ kπ − Z ∗ k1 4cwkn2 p∗ exp(−5np∗ )
,
+
(n̄min I/(mt∗ ))2
n̄min I/(mt∗ )
1

with probability at least 1−exp[−2(n̄min I) 2 ]−n−r −exp(−5np∗ ). By Propositions C.2 and C.3, we have p∗ t∗2  I. Then due to Equation (32), we have
"
#


wknp∗
1
n 2 k
2
 wm
=o p
,
(n̄min I/(mt∗ ))2
n̄min nI
nI/[wk[n/n̄min ]2 ]

26

ZHANG, ZHOU

and

√ ∗
n
wkn2 p∗ exp(−5np∗ )
np
n exp(−5np∗ ) ≤ n exp(−5n̄min I).
 wmk √
n̄min I/(mt∗ )
n̄
nI
min
1

Thus, with probability at least 1 − exp[−(n̄min I) 2 ] − n−r , there exists some
η = o(1), such that
kπ − Z ∗ k1
kht,λ (π) − Z ∗ k1 ≤ n exp(−(1 − η)n̄min I) + p
.
nI/[wk[n/n̄min ]2 ]

The proof for Part One is complete. The very last thing remained to be
obtained is upper bounds on Lsum
and Lsum
1
2 , i.e., Equations (35) and (36).
Recall the definition of θa,b . We have some properties on θa,b which will be
useful in the analysis for Lsum
and Lsum
1
2 : kθa,b k∞ ≤ 2 and
(37)

∗
kθa,b k1 ≤ π·,a − Z·,a

1

∗
+ π·,b − Z·,b

1

≤ kπ − Z ∗ k1 ≤ γ n̄min ,

1

≤ 2k kπ − Z ∗ k1 .

and
(38)

k X
X
a=1 b6=a

kθa,b k1 ≤ 2k

1. Bounds on Lsum
1 .

X
a

∗
π·,a − Z·,a

By applying Markov inequality, we have

EL1,i (a, b, l)


∗
X
X
t (l + 3/2)(na + nb )I
∗
∗
− t∗
(EAi,j − λ)[θa,b ]j 
= P t∗
(Zj,a
− Zj,b
)(Ai,j − λ) ≥ −
4mt
j6=i
j6=i


 ∗

X
t (l + 3/2)(na + nb )I
∗
∗
≤ exp
+ t∗ (EAi,j − λ1Tn )θa,b E exp t∗
(Zj,a
− Zj,b
)(Ai,j − λ) .
4mt
j6=i

With the help of Proposition C.1, we have


X
∗
∗
)(Ai,j − λ)
E exp t∗
(Zj,a
− Zj,b
j6=i

= exp(−t∗ (λ − λ∗ )(na − nb )) exp(−t∗ λ∗ (na − nb ))

Y
j6=i

∗
∗
E exp(t∗ (Zj,a
− Zj,b
)Ai,j )

b

 na −n
tX
2
b
 tX −tY  na +n
Ee
∗
∗
−tλ
2
= exp(−t (λ − λ )(na − nb )) e
Ee
Ee
Ee−tY


(na + nb )I
= exp(−t∗ (λ − λ∗ )(na − nb )) exp −
.
2

27

MEAN FIELD FOR COMMUNITY DETECTION

Hence
(39)
ELsum
1
m−1
k X
XX

"









X

t∗ (l

+ 3/2)(na + nb )I
+ t∗
(EAi,j − λ)[θa,b ]j 
4mt
l=0 a=1 b6=a
j6=i

#
(na + nb )I
× exp(−t∗ (λ − λ∗ )(na − nb )) exp −
2


m−1
k X
t∗ (l+3/2)
l
XX
X
(1 + m − 2mt )(na + nb )I
≤
exp −
(EAi,j − λ)[θa,b ]j  .
− t∗ (λ − λ∗ )(na − nb ) + t∗
2

=

exp −

l(na + nb )I
exp 
2m

l=0 a=1 b6=a

j6=i

We are going to show −(1−η 00 )n̄min I upper bounds terms in the exponent of
RHS of Equation (39) by some η 00 = o(1). We first present some properties
of λ∗ , t∗ and I that will be helpful:
(40)
(41)
(42)

I  (p∗ − q ∗ )2 /p∗ ,
λ∗ ∈ (q ∗ , p∗ ),

and t∗  (p∗ − q ∗ )/p∗ .

Here Equations (40) and (41) are proved by Propositions C.2 and C.3 respectively. Equation (42) is due to t∗  log(1 + (p∗ − q ∗ )/q ∗ )  (p∗ − q ∗ )/p∗
under the assumption that p∗ , q ∗ = o(1), p∗  q ∗ .
The first term in the exponent of Equation (39) is upper bounded by
−(1 − 7/(8m))n̄min I by the assumption t∗ /t = 1 + o(1). Since |t∗ (λ − λ∗ )| ≤
η 0 t∗ (p∗ − q ∗ ), by Equations (40) and (42) the second term is upper bounded
by η 0 n̄min I up to a constant factor. For the last term in the exponent of
Equation (39), since |λ − λ∗ | ≤ η 0 (p∗ − q ∗ ) we have
t∗

X
X
X
(λ∗ − λ)[θa,b ]i
(EAi,j − λ)[θa,b ]i ≤ t∗
(EAi,j − λ∗ )[θa,b ]i + t∗
j6=i

j6=i

j6=i

0

∗

∗

∗

≤ (1 + η )t (p − q ) kθa,b k1
≤ (1 + η 0 )t∗ (p∗ − q ∗ )γ n̄min
. γ n̄min I,

where we use Equations (37) and (40) - (42).
As a consequence, there exists a sequence η 00 = o(1) that goes to zero
slower than m−1 , γ, η 0 , such that the summation of three terms in the exponent of the RHS of Equation (39) is upper bounded by −(1 − η 00 )n̄min I.

28

ZHANG, ZHOU

Thus, Equation (39) can be written as


ELsum
≤ nmk exp −(1 − η 00 )n̄min I .
1

1

Since η 00 goes to 0 slower than m−1 , we have η 00 ≥ m−1 ≥ (n̄min I) 4 by
Equation (32). Then by applying Markov inequality, we have
h
i





1
P Lsum
≥ nmk exp −(1 − 2η 00 )n̄min I ≤ exp −η 00 n̄min I ≤ exp −2(n̄min I) 2 .
1
1

That is, with probability at least 1 − exp[−2(n̄min I) 2 ], Equation (35) holds.

2. Bounds on Lsum
2 . Depending on whether the network is dense or sparse,
we consider two scenarios.
(1) Dense Scenario: q ∗ ≥ (log n)/n. In this scenario, we have a sharp
bound on kA − EAkop . First we observe that
X
X
T
[(Ai,· − EAi,· )θa,b ]2 = θa,b
[(Ai,· − EAi,· )T (Ai,· − EAi,· )]θa,b
i:zi =b

i:zi =b

≤
=

T
θa,b

X

i
T
θa,b [(A

[(Ai,· − EAi,· )T (Ai,· − EAi,· )]θa,b
− EA)T (A − EA)]θa,b .

By applying Markov inequality, we have
Lsum
≤
2

k X T
X
θa,b [(A − EA)T (A − EA)]θa,b
a=1 b6=a

(n̄min I/(4mt))2

.

Since kθa,b k∞ ≤ 2, we have kθa,b k2 ≤ 2 kθa,b k1 . Lemma C.3 shows kA −
√
EAkop ≤ c1 np holds with probability at least 1 − n−r for some constants
c1 , r > 0. Together with Equation (38), we have
k X
X
a=1 b6=a

T
θa,b
[(A

T

− EA) (A − EA)]θa,b ≤
≤

k X
X
a=1 b6=a

k X
X
a=1 b6=a

kA − EAk2op kθa,b k2
2c1 np kθa,b k1

≤ 4c1 knp kπ − Z ∗ k1 .
Thus, with probability at least 1 − n−r ,
Lsum
≤
2

4c1 knp kπ − Z ∗ k1
.
(n̄min I/(4mt))2

29

MEAN FIELD FOR COMMUNITY DETECTION

(2) Sparse Scenario: q ∗ < (log n)/n. When the network is sparse, the
previous upper bound on kA − EAkop no longer holds. Instead, removing
nodes with large degrees is required
to yield provably sharp bound on kA −
P
EAkop . Define S = {i ∈ [n], j Ai,j ≥ 20np∗ }. We define Ã, P̃ such that
Ãi,j = Ai,j I{i, j ∈
/ S} and P̃i,j = (EAi,j )I{i, j ∈
/ S}. Then we have the
decomposition as

X 
n̄min I
L2 (a, b) ,
I (Ai,· − EAi,· )θa,b ≥
4mt
i:zi =b


X
n̄min I
≤
I (Ãi,· − P̃i,· )θa,b ≥
8mt
i:zi =b


X
X
n̄min I 
+
I  (Ai,j − EAi,j )[θa,b ]i,j I{i ∈ S or j ∈ S} ≥
8mt
i:zi =b

j6=i

, L2,1 (a, b) + L2,2 (a, b).
Pk P
Define Lsum
2,1 ,
a=1
b6=a L2,1 (a, b). We have
Lsum
2,1

≤

k X T
X
θa,b [(Ã − P̃ )T (Ã − P̃ )]θa,b

(n̄min I/(8mt))2

a=1 b6=a

≤

k X
X
2kÃ − P̃ k2op kθa,b k
a=1 b6=a

(n̄min I/(8mt))2

1

.

√
Lemma C.4 shows kÃ− P̃ kop ≤ c2 np holds with probability at least 1−n−1
for some constant c2 > 0. Then we have
Lsum
2,1 ≤

4c2 knp kπ − Z ∗ k1
.
(n̄min I/(8mt))2

P
Lemma C.5 shows i,j |Ai,j − EAi,j |I{i ∈ S} ≤ 20n2 p∗ exp(−5np∗ ) holds
with probability at least 1 − exp(−5np∗ ). Then by applying Markov inequality, we have


k
X
X

Lsum
L2,2 (a, b)
2,2 ,
a=1

≤
≤
≤

k
X

b6=a

n
X

a=1 i,j=1

k
X
4
a=1

|Ai,j − EAi,j ||[θa,b ]i,j |I{i ∈ S or j ∈ S}
n̄min I/(8mt)

P

i,j

|Ai,j − EAi,j |I{i ∈ S}
n̄min I/(8mt)

80n2 kp∗ exp(−5np∗ )
.
n̄min I/(8mt)

30

ZHANG, ZHOU

As a consequence, we have
sum
Lsum
≤ Lsum
2
2,1 + L2,2 ≤

4c2 knp∗ kπ − Z ∗ k1 80n2 kp∗ exp(−5np∗ )
+
,
(n̄min I/(8mt))2
n̄min I/(8mt)

with probability at least 1 − n−1 − exp(−5np∗ ). By the bounds on Lsum
and
1
∗ = 1 + o(1), we obtain Equation (36).
Lsum
,
and
due
to
t/t
2
SUPPLEMENTARY MATERIAL
Supplement A: Supplement to “Theoretical and Computational
Guarantees of Mean Field Variational Inference for Community
Detection”
(url to be specified). In the supplement [36], we provide the detailed implementations of the batched Gibbs sampling and an iterative algorithm for
MLE in Algorithm 2 and Algorithm 3 respectively. We include proof of Theorem 4.1, Theorem 4.2 and Theorem 4.3. We also include all the auxiliary
propositions and lemmas in the supplement.
References.
[1] Edoardo M Airoldi, David M Blei, Stephen E Fienberg, and Eric P Xing. Mixed
membership stochastic blockmodels. Journal of Machine Learning Research, 9(Sep):
1981–2014, 2008.
[2] Matthew James Beal. Variational algorithms for approximate Bayesian inference.
University of London, 2003.
[3] Peter Bickel, David Choi, Xiangyu Chang, and Hai Zhang. Asymptotic normality
of maximum likelihood and its variational approximation for stochastic blockmodels.
The Annals of Statistics, 41(4):1922–1943, 2013.
[4] Peter J Bickel and Aiyou Chen. A nonparametric view of network models and
Newman-Girvan and other modularities. Proceedings of the National Academy of
Sciences, 106(50):21068–21073, 2009.
[5] Christopher M Bishop. Pattern recognition and machine learning. springer, 2006.
[6] David M Blei, Andrew Y Ng, and Michael I Jordan. Latent Dirichlet allocation.
Journal of machine Learning research, 3(Jan):993–1022, 2003.
[7] David M Blei, Alp Kucukelbir, and Jon D McAuliffe. Variational inference: A review
for statisticians. Journal of the American Statistical Association, (just-accepted),
2017.
[8] Alain Celisse, Jean-Jacques Daudin, and Laurent Pierre. Consistency of maximumlikelihood and variational estimators in the stochastic block model. Electronic Journal
of Statistics, 6:1847–1899, 2012.
[9] Peter Chin, Anup Rao, and Van Vu. Stochastic block model and community detection
in sparse graphs: A spectral algorithm with optimal rate of recovery. In COLT, pages
391–423, 2015.
[10] Yingjie Fei and Yudong Chen. Exponential error rates of SDP for block models:
Beyond Grothendieck’s inequality. arXiv preprint arXiv:1705.08391, 2017.
[11] Chao Gao, Aad W van der Vaart, and Harrison H Zhou. A general framework for
bayes structured linear models. arXiv preprint arXiv:1506.02174, 2015.

MEAN FIELD FOR COMMUNITY DETECTION

31

[12] Chao Gao, Zongming Ma, Anderson Y Zhang, and Harrison H Zhou. Achieving optimal misclassification proportion in stochastic block model. The Journal of Machine
Learning Research, 18(60):1–45, 2017.
[13] Alan E Gelfand and Adrian FM Smith. Sampling-based approaches to calculating
marginal densities. Journal of the American statistical association, 85(410):398–409,
1990.
[14] Agnieszka Grabska-Barwińska, Simon Barthelmé, Jeff Beck, Zachary F Mainen,
Alexandre Pouget, and Peter E Latham. A probabilistic approach to demixing odors.
Nature neuroscience, 20(1):98–106, 2017.
[15] Alexandre Grothendieck. Résumé de la théorie métrique des produits tensoriels
topologiques. Resenhas do Instituto de Matemática e Estatı́stica da Universidade
de São Paulo, 2(4):401–481, 1996.
[16] Olivier Guédon and Roman Vershynin. Community detection in sparse networks via
Grothendiecks inequality. Probability Theory and Related Fields, 165(3-4):1025–1049,
2016.
[17] Jake M Hofman and Chris H Wiggins. Bayesian approach to network modularity.
Physical review letters, 100(25):258701, 2008.
[18] Paul W Holland, Kathryn Blackmond Laskey, and Samuel Leinhardt. Stochastic
blockmodels: First steps. Social networks, 5(2):109–137, 1983.
[19] Michael I Jordan, Zoubin Ghahramani, Tommi S Jaakkola, and Lawrence K Saul.
An introduction to variational methods for graphical models. Machine learning, 37
(2):183–233, 1999.
[20] Beatrice Laurent and Pascal Massart. Adaptive estimation of a quadratic functional
by model selection. Annals of Statistics, pages 1302–1338, 2000.
[21] Jing Lei and Alessandro Rinaldo. Consistency of spectral clustering in stochastic
block models. The Annals of Statistics, 43(1):215–237, 2015.
[22] Percy Liang, Slav Petrov, Michael I Jordan, and Dan Klein. The infinite PCFG using
hierarchical dirichlet processes. In EMNLP-CoNLL, pages 688–697, 2007.
[23] Yu Lu and Harrison H Zhou. Statistical and computational guarantees of Lloyd’s
algorithm and its variants. arXiv preprint arXiv:1612.02099, 2016.
[24] Elchanan Mossel, Joe Neeman, and Allan Sly. Stochastic block models and reconstruction. arXiv preprint arXiv:1202.1499, 2012.
[25] Mark EJ Newman. Modularity and community structure in networks. Proceedings
of the national academy of sciences, 103(23):8577–8582, 2006.
[26] William D Penny, Nelson J Trujillo-Barreto, and Karl J Friston. Bayesian fMRI time
series analysis with spatial priors. NeuroImage, 24(2):350–362, 2005.
[27] Zahra S Razaee, Arash A Amini, and Jingyi Jessica Li. Matched bipartite block
model with covariates. arXiv preprint arXiv:1703.04943, 2017.
[28] Christian P Robert. Monte carlo methods. Wiley Online Library, 2004.
[29] Karl Rohe, Sourav Chatterjee, and Bin Yu. Spectral clustering and the highdimensional stochastic blockmodel. The Annals of Statistics, 39(4):1878–1915, 2011.
[30] Martin J Wainwright and Michael I Jordan. Graphical models, exponential families,
and variational inference. Foundations and Trends R in Machine Learning, 1(1–2):
1–305, 2008.
[31] Bo Wang and DM Titterington. Convergence properties of a general algorithm for
calculating variational bayesian estimates for a normal mixture model. Bayesian
Analysis, 1(3):625–650, 2006.
[32] Yixin Wang and David M Blei. Frequentist consistency of variational bayes. arXiv
preprint arXiv:1705.03439, 2017.
[33] Ted Westling and Tyler H McCormick. Establishing consistency and improving un-

32

ZHANG, ZHOU

certainty estimates of variational inference through M-estimation. arXiv preprint
arXiv:1510.08151, 2015.
[34] Chong You, John T Ormerod, and Samuel Müller. On variational bayes estimation
and variational information criteria for linear regression models. Australian & New
Zealand Journal of Statistics, 56(1):73–87, 2014.
[35] Anderson Y Zhang and Harrison H Zhou. Minimax rates of community detection in
stochastic block models. The Annals of Statistics, 44(5):2252–2280, 2016.
[36] Anderson Y Zhang and Harrison H Zhou. Supplement to “theoretical and computational guarantees of mean field variational inference for community detection”. 2017.

SUPPLEMENT TO “THEORETICAL AND COMPUTATIONAL
GUARANTEES OF MEAN FIELD VARIATIONAL INFERENCE
FOR COMMUNITY DETECTION”
BY Anderson Y. Zhang and Harrison H. Zhou
Yale University
APPENDIX A: ADDITIONAL ALGORITHMS
In this section, we provide the detailed implementations of the batched
Gibbs sampling and an iterative algorithm of MLE for community detection.
A.1. Batched Gibbs Sampling.
Algorithm 2: Batched Gibbs Sampling

1

Input: Adjacency matrix A, number of communities k, hyperparameters
π pri , αppri , βppri , αqpri , βqpri , some initializers Z (0) , number of iterations S.
Output: Gibbs sampling Ẑ, p̂, q̂.
for s = 1, 2, . . . , S do
(s)
(s)
(s)
(s)
Update αp , βp , αq , βq by
αp(s) = αppri +

k X
X

(s−1)

, βp(s) = βppri +

k X
X
(s−1) (s−1)
(1 − Ai,j )Zi,a Zj,a ,

(s−1)

, βq(s) = βqpri +

XX
(s−1) (s−1)
(1 − Ai,j )Zi,a Zj,b .

(s−1)

Zj,a

(s−1)

Zj,b

Ai,j Zi,a

a=1 i<j

a=1 i<j

αq(s) = αqpri +

XX

Ai,j Zi,a

a6=b i<j

a6=b i<j
(s)

2

(s)

(s)

(s)

Then generate p(s) ∼ Beta(αp , βp ) and q (s) ∼ Beta(αq , βq ) independently.
Define
t(s) =

p(s) (1 − q (s) )
1
log
,
2
(1 − p(s) )q (s)

and λ(s) =

1
1 − q (s)
log
.
(s)
2t
1 − p(s)

Then update π (s) with
π (s) = ht(s) ,λ(s) (Z (s−1) ),
where ht,λ (·) is defined as in Equation (11). Independently generate each row of
Z (s) from distributions
(s)

(s)

P(Zi,· = ea ) = πi,a , ∀a ∈ [k], ∀i ∈ [n].

3

end
We have ẑ = z (S) , p̂ = p(S) and q̂ = q (S) .

2

ZHANG, ZHOU

A.2. An Iterative Algorithm for Maximum Likelihood Estimation. We first define a mapping h0 : Π0 → Π0 as follows


X
[h0λ (Z)]i,a = I a = arg max
Zi,b (Ai,j − λ) .
(43)
b

j6=i

Here if the maximizer is not unique, we simply pick the smallest index.
Algorithm 3: An Iterative Algorithm for MLE

1

Input: Adjacency matrix A, number of communities k, some initializers z (0) ,
number of iterations S.
Output: Estimation Ẑ, p̂, q̂.
for s = 1, 2, . . . , S do
Update p(s) , q (s) by
Pk
p

(s)

P

a=1

= Pk

a=1

(s−1)
Zj,a
(s−1) (s−1)
Ai,j )Zi,a Zj,a
(s−1)

i<j

P

Ai,j Zi,a

i<j (1

−

P

Ai,j Zi,a

and
P
q

(s)

= P

a6=b

a6=b

2

(s−1)

i<j

P

i<j (1

(s−1)

Zj,b

(s−1)

− Ai,j )Zi,a

(s−1)

.

Zj,b

Define
t(s) =

p(s) (1 − q (s) )
1
log
,
2
(1 − p(s) )q (s)

and λ(s) =

1
1 − q (s)
log
.
(s)
2t
1 − p(s)

Then update π (s) with
Z (s) = h0λ(s) (Z (s−1) ),
where h0λ (·) is defined as in Equation (43).
3

end
We have ẑ = z (S) , p̂ = p(S) and q̂ = q (S) .

APPENDIX B: PROOFS OF OTHER THEOREMS
In this section, we first validate Theorem 3.1 when `(π (0) , π ∗ ) is in a
constant order of n̄min , which complements the proof presented in Section
5.3. The we give proofs of theorems stated in Section 4, including Theorem
4.1, Theorem 4.2 and Theorem 4.3.
B.1. Proof of Theorem 3.1 for the case `(π (0) , π ∗ ) in a constant
order of n̄min . For any π such that `(π, π ∗ ) ≤ cinit n̄min , we are going to

3

MEAN FIELD FOR COMMUNITY DETECTION

show when cinit is sufficiently small
(44)

`(π, Z ∗ )
`(ht,λ (π), Z ∗ ) ≤ n exp(−n̄min I/25) + p
,
2 nI/[wk[n/n̄min ]2 ]

with probability at least 1−exp(−n̄min I/10)−n−r for some constant r > 0. If
it holds, for any π (0) such that `(π (0) , Z ∗ ) = cn̄min for some
pconstant c ≤ cinit ,
(0)
∗
the term n exp(−n̄min I/25) is dominated by `(π , Z )/ nI/[wk[n/n̄min ]2 ]
which implies
`(π (1) , Z ∗ ) ≤ n exp(−(1 − η)/n̄min I) + p

`(π (0) , Z ∗ )
nI/[wk[n/n̄min ]2 ]

.

It also implies `(π (1) , Z ∗ ) = o(n̄min ), which means after the first iteration,
the results in Section 5.3 can be directly applied and the proof is complete.
The proof of Equation (44) mainly follows the proof of Part One in Section
5.3. We have


X
X
∗
[ht,λ (π)]i,· − Zi,·
≤ 2w
1 ∧ exp 2t
(πj,a − πj,zi )(Ai,j − λ) .
1
a6=zi

j6=i

Note that the inequality 1 ∧ exp(−x) ≤ f (x0 ) + I{x ≥ x0 } holds for any
x0 ≥ 0. By taking x0 = (na + nzi )I/4, we have





X
X
(n
+
n
)I
(n
+
n
)I
zi
a
zi
∗
exp − a
 ,
≤ 2w
+ I  (πj,a − πj,zi )(Ai,j − λ) ≥ −
[ht,λ (π)]i,· − Zi,·
1
4
8t
a6=zi

j6=i

and consequently,

kht,λ (π) − Z ∗ k1 ≤ 2wnk exp(−n̄min I/2)
#
k X X X
X
(na + nzi )I
+ 2w
.
I
(πj,a − πj,b )(Ai,j − λ) ≥ −
8t
a=1 b6=a i:zi =b

j6=i

Define θa,b the same way as in Section 5.3, and by the same argument, we
have

k X X 
X
n̄min I
∗
kht,λ (π) − Z k1 ≤ 2wnk exp(−n̄min I/2) + 2w
I (Ai,· − EAi,· )θa,b ≥
8t
a=1 b6=a i:zi =b


k X X X
X
(na + nb )I X
∗
∗
+ 2w
I
(Zj,a − Zj,b )(Ai,j − λ) ≥ −
−
(EAi,j − λ)[θa,b ]j .
4t
a=1 b6=a i:zi =b

j6=i

j6=i

4

ZHANG, ZHOU

From Lemma C.1, when cinit is sufficiently small, with probability at least
1 − e3 5−n we have


|λ − λ∗ |
|t − t∗ |
(45)
,
≤ 24c0 cinit .
max
(p∗ − q ∗ )/p∗ (p∗ − q ∗ )
Proposition C.3 shows that λ∗ ∈ (q ∗ + c(p∗ − q ∗ ), q ∗ + (1 − c)(p∗ − q ∗ )) for
some positive constant 0 < c < 1/2. Therefore, when cinit is sufficiently
small, we have λ ∈ (q ∗ , p∗ ). Thus,
X
j6=i

(EAi,j − λ)[θa,b ]j ≤ (p∗ − q ∗ ) kθa,b k1 ≤ (p∗ − q ∗ ) kπ − Z ∗ k1 ≤ cinit (p∗ − q ∗ )n̄min ,

where we use Equation (37). By Equations (40) - (42), it is smaller than
(na + nzi )/(8t) when cinit is sufficiently small. As a consequence, we have
∗

kht,λ (π) − Z k1 ≤ 2wnk exp(−n̄min I/2) + 2w
+ 2w

k X
X

a=1 b6=a i:zi =b



I (Ai,· − EAi,· )θa,b

n̄min I
≥
8t


X X
(na + nb )I
∗
∗
)(Ai,j − λ) ≥ −
(Zj,a
− Zj,b
I
.
8t

a=1 b6=a i:zi =b

Define Lsum
1
sum
and L2 =

k X X
X

Pk

=
Pk

j6=i

a=1

P

P

b6=a

P

P

i:zi =b I

hP

b6=a

i:zi =b I [(Ai,·

∗
j6=i (Zj,a

−

∗ )(A
Zj,b
i,j

i
− λ) ≥ −(na + nb )I/(8t)

− EAi,· )θa,b ≥ n̄min I/(8t)]. Our analysis on them is quite similar to that in Section 5.3. By Markov inequality,


k X X
X
X
∗
∗
)(Ai,j − λ) ≥ −t∗ (na + nb )I/(8t)
(Zj,a
− Zj,b
P t∗
ELsum
=
1
a=1

a=1 b6=a i:zi =b

≤
≤

k X X
X

a=1 b6=a i:zi =b

k X
X



X

a=1 b6=a i:zi =b

j6=i

exp





t∗ (n




X
a + nb )I
∗
∗
− t∗ (λ − λ∗ )(na − nb ) E exp t∗
(Zj,a
− Zj,b
)(Ai,j − λ∗ )
8t
j6=i


t∗ (na + nb )I
(na + nb )I
∗
∗
− t (λ − λ )(na − nb ) −
.
exp
8t
2

By Equations (40) - (42) and (45), when cinit is small enough, t∗ /t ≤ 2 and
t∗ |λ − λ∗ | ≤ I/6. Thus
ELsum
≤ nk exp(−n̄min I/12).
1

5

MEAN FIELD FOR COMMUNITY DETECTION

Hence, with probability at least 1 − exp(−n̄min I/24),
Lsum
≤ nk exp(−n̄min I/24).
1

For Lsum
we use the same argument as in Section 5.3 and obtain
2
Lsum
≤
2

4c2 knp∗ kπ − Z ∗ k1 80n2 kp∗ exp(−5np∗ )
+
,
(n̄min I/(8t))2
n̄min I/(8t)

with probability at least 1−n−r −exp(−5np∗ ) for some constants r, c1 , c2 > 0.
Recall that
kht,λ (π) − Z ∗ k1 ≤ 2wnk exp(−n̄min I/2) + 2wLsum
+ 2wLsum
1
2 .
Using the same argument as in Section 5.3, we conclude with
1
kπ − Z ∗ k1 ,
kht,λ (π) − Z ∗ k1 ≤ n exp(−n̄min I/25) + p
2 nI/[wk[n/n̄min ]2 ]

with probability at least 1 − exp(−n̄min I/10) − n−r .
B.2. Proof of Theorem 4.1. Define t∗ =
1
2t∗

∗

1
2

∗

∗

(1−q )
∗ =
log qp∗ (1−p
∗ ) and λ

1−q
log 1−p
∗ . By the same simplification we derive in Theorem 2.1, we have

π̂ MF = arg max f 0 (π; A),
π∈Π1

where
f 0 (π; A) = hA + λ∗ In − λ∗ 1n 1Tn , ππ T i −

n
1 X
pri
KL(Categorical(πi,· )kCategorical(πi,·
)).
t∗
i=1

Recall the definition of ht,λ (·) as in Equation (11). A key observation is
that π̂ MF = ht∗ ,λ∗ (π̂ MF ), otherwise if there exists some i ∈ [n] such that
MF . This indicates the implementation of CAVI
[ht∗ ,λ∗ (π̂ MF )]i,· not equal to π̂i,·
update on the i-th row of π will make change, leading to the decrease of
f 0 (·; A). This contradicts with the fact that π̂ MF is the global minimizer.
The fixed-point property of π̂ MF is the key to our analysis. It involves
three steps.
• Step One. For any π such that `(π, Z ∗ ) = o(n̄min ), by the same analysis
as in the proof of Theorem 3.1, we are able to show that there exist
constant r > 0 and sequence η = o(1) such that
kht∗ ,λ∗ (π) − Z ∗ k1 ≤ n exp(−(1 − η)n̄min I) + p
1

kπ − Z ∗ k1

nI/[wk[n/n̄min ]2 ]

with probability at least 1 − exp[−(n̄min I) 2 ] − n−r .

,

6

ZHANG, ZHOU

• Step Two. Lemma C.6 presents some loose upper bound for `(π̂ MF , Z ∗ ).
That is, under the assumption ρnI/[wk 2 [n/n̄min ]2 ] → ∞, with probability at least 1 − e3 5−n , we have
`(π̂ MF , Z ∗ ) ≤ o(n̄min ).
• Step Three. Using the property that ht∗ ,λ∗ (π̂ MF ) = π̂ MF , we have
π̂ MF − Z ∗

1

≤ n exp(−(1 − η)n̄min I) + p

π̂ MF − Z ∗

1

nI/[wk[n/n̄min ]2 ]
1

holds with probability at least 1 − exp[−(n̄min I) 2 ] − n−r . Then we
obtain the desired result by simple algebra.
B.3. Proof of Theorem 4.2. By law of total expectation, we have
(46)
EZ (s+1)

h

Z

(s+1)

−Z

∗

1

A, Z

(0)

i



h
i
(0)
(s+1)
∗
(s+1)
(0)
= Eπ(s+1) EZ (s+1) Z
A, Z
−Z
π
, A, Z
1
h
i
= Eπ(s+1) π (s+1) − Z ∗
A, Z (0) ,
1

where the first equation is due to that the conditional expectation of Z (s+1)
is π (s+1) . We are going to build the connection between π (s) and π (s+1) . In
Algorithm 2, there are intermediate steps between π (s) and π (s+1) as follows:
π (s)

Z (s)

(p(s+1) , q (s+1) ) → (t(s+1) , λ(s+1) ) → π (s+1) ,

where we use the plain right arrow (→) to indicate deterministic generation
and the curved right arrow ( ) to indicate random generation. Despite a
slight abuse of notation, we define π (0) = Z (0) .
Analogous to the proof of Theorem 3.1 in Section 5.3, we assume `(Z (0) , Z ∗ ) =
o(n̄min ). The proof for the case `(Z (0) , Z ∗ ) in the same order of n̄min is similar
and thus is omitted.
Let γ = o(1) be any sequence goes to 0 when n grows. We define a series
of events as follows:
• global event F: We define F exactly the same way as we define in the
proof of Theorem 3.1 in Section 5.3 with respect to sequences γ and
1
η 0 , and we have P(F) ≥ 1 − exp[−(n̄min I) 2 )] − n−r for some constant
r > 0. We have η 0 = o(1) whose value will be determined later.

7

MEAN FIELD FOR COMMUNITY DETECTION

• global event G: Consider any Z ∈ Π1 such that kZ − Z ∗ k1 ≤ γ n̄min .
Define
αp = αppri +

k X
X

Ai,j Zi,a Zj,a , βp = βppri +

a=1 i<j

αq = αqpri +

XX

k X
X
(1 − Ai,j )Zi,a Zj,a ,
a=1 i<j

Ai,j Zi,a Zj,b , βq = βqpri +

a6=b i<j

XX
a6=b i<j

(1 − Ai,j )Zi,a Zj,b .

Define G be the event that


αq
αp
∗
∗
−p ,
−q
≤ η 00 (p∗ − q ∗ )
max
αp + βp
αq + β q
holds uniformly over all the eligible Z for some sequence η 00 = o(1).
Then by the same analysis as in Lemma C.1, we have P(G) ≥ 1−e3 5−n .
(s)
(s)
• local events {H1 }Ss=1 : We define H1 = { π (s) − Z ∗ 1 ≥ γ n̄min /2}.
(s)

(s)

• local events {H2 }Ss=1 : We define H2 = { Z (s) − Z ∗
the conditional probability, we have
(s)

≤P

≥ γ n̄min }. For

(s)

P(H2 = 1|H1 = 0)
" n
X h (s)
∗
≤P
Zi,· − Zi,·
"

1

i=1
n h
X

(s)
Zi,·

i=1

Since π (s) − Z ∗
we have
(s)
P(H2

=

−

∗
Zi,·

1

−

∗
− Zi,·

−

(s)
πi,·

∗
Zi,·

−

1

1

i

i

≥ γ n̄min − π (s) − Z ∗
≥ γ n̄min /2

(s)
H1

#

1

(s)
H1

=0

≤ γ n̄min /2 given H1 = 0 by Bernstein inequality,
"

(γ n̄min )2 /8
= 0) ≤ exp − (s)
π − Z ∗ 1 + γ n̄min /6


≤ exp −3(γ n̄min )2 /16 .
(s)

#

=0

(s)

1

(s)
1|H1

(s)

1

(s)
πi,·

#

• local events {H3 }Ss=1 : We define H3 = {|t(s) −t∗ | ≥ η 0 (p∗ −q ∗ )/p∗ , or |λ(s) −
(s)
λ∗ | ≥ η 0 (p∗ − q ∗ )}. If the global event G holds and the local event H2
does not hold, we have
(
)
(s+1)
(s+1)
αp
αq
∗
∗
max
− p , (s+1)
−q
≤ η 00 (p∗ − q ∗ ).
(s+1)
(s+1)
(s+1)
αp
+ βp
αq
+ βq

8

ZHANG, ZHOU

P
P
(s) (s)
(s+1)
(s+1)
Note that αp
+ βp
= αppri + βppri + ka=1 i<j Zi,a Zj,a ≥ n2 /k.
Using the tail bound of Beta distribution (Lemma C.7) we are able to
show
"
#
(s+1)
α
p
(s)
P p(s+1) − (s+1)
≥ η 00 (p∗ − q ∗ ) H2 = 0, G = 1
(s+1)
αp
+ βp


∗
∗ 2
002 2 (p − q )
≤ exp −η n
2p∗
 002 2

≤ exp −η n I/2 ,
where the last inequality is due to Proposition C.2. This leads to
h
i


(s)
P p(s+1) − p∗ ≥ 2η 00 (p∗ − q ∗ ) H2 = 0, G = 1 ≤ exp −η 002 n2 I/2 .

And similar result holds for q (s+1) . Then by the same analysis as in the
proof of Lemma C.1, max{|p(s+1) − p∗ |, |q (s+1) − q ∗ |} ≤ 2η 00 (p∗ − q ∗ )
leads to
(
)
∗
|t(s+1) − t∗ | |λ(s+1)−λ |
max
,
≤ 16c0 η 00 .
(p∗ − q ∗ )/p∗ p∗ − q ∗
By taking η 0 = 16c0 η 00 , we obtain
(s+1)

P(H3



(s)
= 1|H2 = 0, G = 1) ≤ 2 exp −η 002 n2 I/2 .

Note that events F and G are about the adjacency matrix A. The events
(s+1)
and H3
are for π (x) , Z (s) and (p(s+1) , q (s+1) ) respectively. With
all the above events defined, we can continue our analysis for Equation (46).
(s)
(s)
(s+1) C
Under the event F ∩ G ∩ (H1 ∪ H2 ∪ H3
) we have

(s)
(s)
H1 , H2

π (s+1) − Z ∗

(47)

1

≤ n exp(−(1 − η)n̄min I) + cn π (s) − Z ∗

1

,

where cn = [nI/[wk[n/n̄min ]2 ]]−1/2 . As a consequence, under the event F ∩
Q
(v)
(v)
(v+1) C
) , we have
G ∩ ( sv=0 H1 ∪ H2 ∪ H3
π (s+1) − Z ∗

1

≤ n exp(−(1 − 2η)n̄min I) + csn π (0) − Z ∗

1

.

Therefore, we have
(48)
Eπ(s+1)

h

i
(0)
H1 = 0, F = 1, G = 1 ≤ n exp(−(1 − 2η)n̄min I)
1
#
" s
Y (v)
(0)
(v)
(v+1)
H1 = 0, F = 1, G = 1 .
+ nP
H1 ∪ H 2 ∪ H 3

π (s+1) − Z ∗

+ csn π (0) − Z ∗

1

v=1

MEAN FIELD FOR COMMUNITY DETECTION

9

Due to the small value of cn , if π (s) − Z ∗ 1 ≤ γ n̄min , Equation (47) immediately implies π (s+1) − Z ∗ 1 ≤ γ n̄min . This implies that under the event
F ∪ G we have
(s+1)

H1

(s)

(s)

(s+1)

⊂ H1 ∪ H2 ∪ H3

, ∀s ≥ 0,

and consequently,
s
Y

v=0

(v)

(v)

(v+1)

H1 ∪ H 2 ∪ H 3

(0)

⊂ H1

s
Y

v=0

(v)

(v+1)

H2 ∪ H3

, ∀s ≥ 1.

Thus,
(49)

P

"

s
Y

v=0

≤P
≤

"

(v)
H1

s
Y

v=0

s
X

∪

(v)
H2

(v)
H2

∪

(v)

∪

(v+1)
H3

(v+1)
H3
(v)

(0)
H1

(0)
H1

= 0, F = 1, G = 1
#

= 0, F = 1, G = 1

P(H2 = 1|H1 = 0) +

v=0

n
X

(v+1)

P(H3

v=0

(v)

= 1|H2 = 0, G = 1)





≤ (s + 1) exp −3(γ n̄min )2 /16 + 2 exp −η 002 n2 I/2 .
(0)



#

1

Note that P(H1 = 0, F = 1, G = 1) ≥ 1−exp[−(n̄min I) 2 )]−n−r −e3 5−n −.
Recall we define π (0) = Z (0) . By Equations (46), (48) and (49), we have
h
i
EZ (s+1) Z (s+1) − Z ∗
A, Z (0) ≤ n exp(−(1 − 2η)n̄min I) + csn Z (0) − Z ∗ + (s + 1)nbn ,
1

1

1

−r
3 −n
2
with probability at least
1 − exp[−(n̄

 002 2min I)
 )] − n − e 5 − , where bn =
2
exp −3(γ n̄min ) /16 + 2 exp −η n I/2 .

B.4. Proof of Theorem 4.3. Note the similarity between Algorithm
3 and Algorithm 1. We can prove Theorem 4.3 with almost the identical
argument used in the proof of Theorem 3.1, thus omitted.
APPENDIX C: STATEMENTS AND PROOFS OF AUXILIARY
LEMMAS AND PROPOSITIONS
We include all the auxiliary propositions and lemmas in this section.

10

ZHANG, ZHOU

C.1. Statements and Proofs of Lemmas and Propositions for
Theorem 3.1.
Lemma C.1. Let cinit be some sufficiently small constant. Consider any
π ∈ Π1 such that kπ − Z ∗ k1 ≤ cinit n/k. Let αp , βp , αq , βq , t, λ be the outputs
after one step CAVI iteration from π described in Algorithm 1. That is, they
are defined as Equations (26) - (29). Define
P P
P Pk
i<j
a6=b πi,a πj,b Ai,j
i<j
a=1 πi,a πj,a Ai,j
, and q̂ = P P
p̂ = P Pk
.
i<j
a6=b πi,a πj,b
a=1 πi,a πj,a
i<j

Under the same assumption as in Theorem 3.1, there exists some sequence
 = o(1) such that with probability at least 1 − e3 5−n , the following inequality
holds


kπ − Z ∗ k1
|p̂ − p∗ | |q̂ − q ∗ |
|t − t∗ |
|λ − λ∗ |
max
,
,
,
≤

+
24c
,
0
p∗ − q ∗ p∗ − q ∗ (p∗ − q ∗ )/p∗ p∗ − q ∗
n/k

uniformly over all the eligible π. In addition if we further assume cinit goes
to 0, the LHS of the above inequality will be simply upper bounded by .
Proof. We are going to obtain tight bounds on |p̂ − p∗ | and |q̂ − q ∗ | first.
Note that we have the “variance-bias” decomposition as in
P Pk
P P
| i<j ka=1 πi,a πj,a (Ai,j − EAi,j )|
i<j
a=1 πi,a πj,a EAi,j
∗
+
− p∗ .
|p̂ − p | ≤
P Pk
P Pk
a=1 πi,a πj,a
a=1 πi,a πj,a
i<j
i<j

We have concentration inequality holds for the numerator in the first term
by Lemma C.2. That is, with probability at least 1 − e3 5−n , we have
k
XX
i<j a=1

πi,a πj,a (Ai,j − EAi,j ) =

p
1
hA − EA, ππ T i ≤ 3n np∗
2

holds uniformly over all π ∈ Π1 . For the denominator, we have
k

k

n2 X X
1X
n2
≥
πi,a πj,a =
kπ·,a k21 ≥
,
2
2
2k
i<j a=1

a=1

P
since ka=1 kπ·,a k1 = n. Thus, we are able to obtain an upper bound on the
first term as
r
P P
| i<j ka=1 πi,a πj,a (Ai,j − EAi,j )|
k 2 p∗
.
≤6
P Pk
n
π
π
i,a
j,a
i<j
a=1

11

MEAN FIELD FOR COMMUNITY DETECTION

For the second term, since EAi,j = p∗
we have
P

Pk

a=1 πi,a πj,a EAi,j
Pk
i<j
a=1 πi,a πj,a

i<j

P

Pk

∗
∗
∗
a=1 Zi,a Zj,a + q (1 −

− p∗ = (p∗ − q ∗ )

hP

P

∗
∗
a=1 Zi,a Zj,a ),

k
a=1 πi,a πj,a

i<j

P

hππ T , 11T −
= (p∗ − q ∗ ) P Pk

= (p∗ − q ∗ )

Pk

i hP

Pk
i<j
a=1 πi,a πj,a
Z ∗ Z ∗T i

i<j
a=1 πi,a πj,a
hππ T − Z ∗ Z ∗T , 11T −

P

i<j

k
a=1 1

Pk

Z ∗ Z ∗T i

a=1 πi,a πj,a

∗ Z∗
− Zi,a
j,a

i

,

where in the last inequality we use the orthogonality between Z ∗ Z ∗T and
11T − Z ∗ Z ∗T . For its numerator, we have
hππ T − Z ∗ Z ∗T , 11T − Z ∗ Z ∗T i ≤ ππ T − Z ∗ Z ∗T

1

∗

≤ kπ − Z k1 (kπk1 + kZ ∗ k1 )

≤ kπ − Z ∗ k1 (2 kZ ∗ k1 + kπ − Z ∗ k1 )
≤ 3n kπ − Z ∗ k1 .

This leads to
P Pk
3n kπ − Z ∗ k1 (p∗ − q ∗ )
i<j
a=1 πi,a πj,a EAi,j
≤ 3kn−1 (p∗ − q ∗ ) kπ − Z ∗ k1 .
− p∗ ≤
P Pk
2 /k
n
π
π
a=1 i,a j,a
i<j
Thus,

|p̂ − p∗ | ≤ 6

r

k 2 p∗
+ 3kn−1 (p∗ − q ∗ ) kπ − Z ∗ k1 ≤
n

Similar result holds for |q̂ − q ∗ |. Denote η0 =

q

"s

#
3 kπ − Z ∗ k1
k 2 p∗
+
(p∗ − q ∗ ).
n(p∗ − q ∗ )2
n/k

k 2 p∗
n(p∗ −q ∗ )2

+

3kπ−Z ∗ k1
,
n/k

thus

max{|p̂ − p∗ |, |q̂ − q ∗ |} ≤ η0 (p∗ − q ∗ ).
By the assumption of nI in Equation (18) and Proposition C.2, we have
n(p∗ − q ∗ )2 /(k 2 p∗ )  nI/k 2 → ∞. Therefore, the first term in η0 goes to 0.
The second term in η0 is at most 3cinit which implies η0 ≤ 4cinit .
By the fact that the digamma function satisfies ψ(x) ∈ (log(x−1/2), log x), ∀x ≥

12

ZHANG, ZHOU

1/2, we have
αp − 1/2
βp
i  hP P
i
 h pri
P P
k
αp − 1/2 + i<j ka=1 πi,a πj,a Ai,j
π
π
i<j
a=1 i,a j,a
h
i  hP P
i 
= log 
P
P
k
1 + βppri − i<j ka=1 πi,a πj,a Ai,j
π
π
i<j
a=1 i,a j,a
hP P
i


k
p̂ + (αppri − 1/2)
i<j
a=1 πi,a πj,a

h
i .
= log
Pk
pri  P
1 − p̂ + βp
i<j
a=1 πi,a πj,a

ψ(αp ) − ψ(βp ) ≥ log

P P
Recall that we have shown i<j ka=1 πi,a πj,a lies in the interval of (n2 /(2k), n2 /2).
By Equation (18), there exists a sequence η 0 = o(1) such that αp , βp ≤
η 0 (p∗ − q ∗ )n2 /k. Then we have
ψ(αp ) − ψ(βp ) ≥ log

p∗ − |p∗ − p̂| − η 0 (p∗ − q ∗ )
.
1 − p∗ + |p∗ − p̂| + η 0 (p∗ − q ∗ )

Similar analysis leads to
ψ(αq ) − ψ(βq ) ≤ log

q ∗ + |q ∗ − q̂| + η 0 (p∗ − q ∗ )
.
1 − q ∗ − |q ∗ − q̂| − η 0 (p∗ − q ∗ )

Together we have

 ∗
p − |p∗ − p̂| − η 0 (p∗ − q ∗ ) 1 − q ∗ − |q ∗ − q̂| − η 0 (p∗ − q ∗ )
∗
− t∗
t − t ≥ log
1 − p∗ + |p∗ − p̂| + η 0 (p∗ − q ∗ ) q ∗ + |q ∗ − q̂| + η 0 (p∗ − q ∗ )
"
#

|p∗ − p̂| + η 0 (p∗ − q ∗ ) 4 p∗ (1 − q ∗ )
≥ log 1 −
− t∗
q∗
q ∗ (1 − p∗ )


∗
∗
0 p −q
= 4 log 1 − (η0 + η )
.
q∗
Recall that we assume c0 p∗ < q ∗ < p∗ . Thus (η0 + η 0 )(p∗ − q ∗ )/p∗ ≤ 5cinit c0 .
When cinit is sufficiently small, we have (η0 + η 0 )(p∗ − q ∗ )/p∗ ≤ 1/2. Then
using the fact −x ≥ log(1 − x) ≥ −2x, ∀x ∈ (0, 1/2). We have
t − t∗ ≥ −8(η0 + η 0 )(p∗ − q ∗ )/q ∗ .
Analogously we can obtain the same upper bound on t̂ − t∗ , and then
|t − t∗ | ≤ 8c0 (η0 + η 0 )

p∗ − q ∗
.
p∗

MEAN FIELD FOR COMMUNITY DETECTION

13

Identical analysis can be applied towards bounds on |λ̂ − λ∗ |. Note that
h
i 

Pk
pri  P
1
−
p̂
+
β
π
π
p
i,a
j,a
i<j
a=1
βp
i ,
= log 
log
 hP Pk
αp + βp
1 + (αppri + βppri )
π π
i<j

a=1

i,a j,a

similarly for αq , βq . Omitting the immediate steps, we end up with

|λ − λ∗ | = | [ψ(βq ) − ψ(αq + βq )] − [ψ(βp ) − ψ(αp + βp )] − λ∗ | ≤ 8(η0 + η 0 )(p∗ − q ∗ ).
The proof is complete after we unify and rephrase all the aforementioned
results.
Lemma C.2. Let A ∈ [0, 1]n×n such that A = AT and Ai,i = 0, ∀i ∈ [n].
Assume {Ai,j }i<j are independent
random variable, and there exists p ≤ 1
P
2
−1
such that 9n ≤ n(n−1) i<j Var(Ai,j ) ≤ p, and then we have
√
sup hA − EA, ππ T i ≤ 6n np,

π∈Π1

with probability at least 1 − e3 5−n .
Proof. This result is a direct consequence of Grothendieck inequality
[15] (see also Theorem 3.1 of [16] for a rephrased statement) on the matrix
A−EA. The Lemma 4.1 of [16] proves that with probability at least 1−e3 5−n ,
X
√
sup
(Ai,j − EAi,j )si tj ≤ 3n np.
s,t∈{−1,1}n

i,j

Then by applying Grothendieck inequality we obtain
X
√
sup
(Ai,j − EAi,j )XiT Xj ≤ 3cn np,
kXi k2 ≤1,∀i∈[n]

i,j

where c is a positive constant smaller than 2. This concludes with
√
sup hA − EA, ππ T i ≤ 6n np,
π∈Π1

Proposition C.1. Assume 0 < q < p < 1. Let X ∼ Ber(q) and Y ∼
p(1−q)
1−q
1
Ber(p). Recall the definition λ = log 1−p
/ log p(1−q)
q(1−p) , t = 2 log q(1−p) and
p
√
I = −2 log[ pq + (1 − p)(1 − q)]. Then the following two equations hold
(50)

e

tλ

=



EetX
Ee−tY

 12

, and EetX Ee−tY = exp(−I).

14

ZHANG, ZHOU

Proof. The proof is straightforward and all by calculation. Note that
E exp(tX) = pet + 1 − p and E exp(tY ) = qet + 1 − q. We can easily obtain
p
√
EetX Ee−tY = (pet + 1 − p)(qe−t + 1 − q) = ( pq + (1 − p)(1 − q))2 = exp(−I).

We can justify the first part of Equation (50) in a similar way.

Lemma C.3. [Theorem 5.2 of [21]] Let A ∈ {0, 1}n×n be a symmetric binary matrix with Ai,i = 0, ∀i ∈ [n], and {Ai,j }i<j are independent Bernoulli
random variable. If p , maxi,j EAi,j ≥ log n/n. Then there exist constants
c, r > 0 such that
√
kA − EAkop ≤ c np,
with probability at least 1 − n−r .
The following lemma on the operator norm of sparse networks is from [9].
In the original statement of Lemma 12 in [9], “with probability 1 − o(1)” is
stated. However, its proof in [9] gives explicit form of the probability that
the statement holds, which is at least 1 − n−1 .
Lemma C.4. [Lemma 12 of [9]] Suppose M is random symmetric matrix
with zero on the diagonal whose entries above the diagonal are independent
with the following distribution
(
1 − pi,j , w.p. pi,j ;
Mi,j =
−pi,j , w.p. 1 − pi,j .
Let p , maxi,j pi,j and M̃ be the matrix obtained from M by zeroing out all
the rows and columns having more than 20np positive entries. Then there
exists some constant c > 0 such that
√
kM̃ kop ≤ c np,
holds with probability at least 1 − n−1 .
Lemma C.5. Let A ∈ {0, 1}n×n be a symmetric binary matrix with Ai,i =
0, ∀i ∈ [n], and {Ai,j }i<j are independent
Bernoulli random variable.
Let
P
P
p ≥ maxi,j EAi,j . Define S = {i ∈ [n], j Ai,j ≥ 20np} and Zi = j |Ai,j −
EAi,j |I{i ∈ S}. Then with probability at least 1 − exp(−5np), we have
X
Zi ≤ 20n2 p exp(−5np).
i

MEAN FIELD FOR COMMUNITY DETECTION

15

P
Proof. Note that E j |Ai,j − EAi,j | ≤ 2np(1 − p) ≤ 2np. For any s ≥
20np, we have


X
X
P(Zi > s) ≤ P 
|Ai,j − EAi,j | − E
|Ai,j − EAi,j | > s − 2np
j

"

≤ exp −

1
2 (s

− 2np)2
np + 13 (s − 2np)

#

j

≤ exp(−s/2),

by implementing Bernstein inequality. Applying Bernstein inequality again
we have


X
P(Zi > 0) = P 
Ai,j ≥ 20np
j



X
X
≤ P
Ai,j − E
Ai,j ≥ 18np
j



≤ exp −

j

(18np)2 /2

np + 18np/3
≤ exp(−21np/2).



Thus, we are able to bound EZi with
Z 20np
Z ∞
EZi ≤
P(Zi > 0) ds +
P(Zi > s) ds
0
20np
Z ∞
≤ 20np exp(−21np/2) +
exp(−s/2)
20np

≤ 20np exp(−10np).

By Markov inequality, we have


"
#
X
X
2
2
P
|Ai,j − EAi,j |I{i ∈ S} ≥ 20n p exp(−5np) = P
Zi ≥ 20n p exp(−5np)
i,j

i

≤

nEZ1
2
20n p exp(−5np)

≤ exp(−5np).

16

ZHANG, ZHOU

Proposition
C.2. Under the iassumption that 0 < q < p = o(1). For
h√
p
I = −2 log pq + (1 − p)(1 − q) we have
√
√
I = (1 + o(1))( p − q)2 .

Consequently, (p − q)2 /(4p) ≤ I ≤ (p − q)2 /p.
Proof. It is a partial result of Lemma B.1 in [35].
p(1−q)
1−q
Proposition C.3. Define λ = log 1−p
. For any p, q > 0 such
/ log q(1−p)
that p, q = o(1) and p  q, there exists a constant 0 < c < 1/2 such that

λ−q
∈ (c, 1 − c).
p−q
Proof. First we are going to establish the lower bound. Let x = p − q,
and then we can rewrite λ as
λ=

1
1+

log(1+x/q)
log(1+x/(1−q−x))

.

Case I: x ≥ q/10. Define s = (p − q)/q. Since p  q we have s ≥ 1/10 and
also upper bounded by some constant. We have


λ−q
1 1
1
=
− 1
log(1+s)
p−q
s q1+
log(1+sq/(1−(s+1)q))


1 (1 − q) log(1 + sq/(1 − (s + 1)q)) − q log(1 + s)
=
s
q log(1 + sq/(1 − (s + 1)q)) + q log(1 + s)
sq
1 (1 − q) 1−(s+1)q − q log(1 + s)
≥
s
2q log(1 + s)
1 1−q
≥
,
8 log(1 + s)
which is lower bounded by some constant c > 0.
Case II: x < q/10. By Taylor theorem, there exist constants 0 ≤ 1 , 2 ≤
1/10 such that


 
x
x 1 − 1 x 2
log 1 +
,
= −
q
q
2
q



2
x
x
1 − 2
x
and log 1 +
=
−
.
1−q−x
1−q−x
2
1−q−x

MEAN FIELD FOR COMMUNITY DETECTION

17

Thus, we have
log(1 + xq )
log(1 +

x
1−q−x )

=


q(1 − q)2 − 2q(1 − q) +
q 2 (1 −

1−1
2
2 (1 − q)
2 2
q) − 3−
2 q x



x + c1 x2 + c2 x3

,

where c1 = (1 − 1 )(1 − q) + q and c2 = −(1 − 1 )/2. Thus,
"
#
2 2
q 2 (1 − q) − 3−
q x
1
λ−q
2


=
−q
3−2 2
1
2
2
3
p−q
x q(1 − q) − 2q(1 − q) + 1−
2 (1 − q) + 2 q x + c1 x + c2 x
1

2 2
1
2
2
2 q(1 − q) + 2 q (1 − q) − 2 (1 − q) q + c1 qx + c2 qx
=
3−2 2
1
2
2
3
q(1 − q) − 2q(1 − q) + 1−
2 (1 − q) + 2 q x + c1 x + c2 x
Note that |c1 |, |c2 | ≤ 1. We have

1
q(1 − q)
λ−q
≥ 4
≥ 1/8.
p−q
2q(1 − q)

By using exactly the same discussion, we can show (p − λ)/(p − q) > c.
Thus, we proved the desired bound stated in the proposition.
C.2. Statements and Proofs of Lemmas and Propositions for
Theorem 4.1.
Lemma C.6. Let Z ∗ ∈ Π0 . Assume p∗ , q ∗ = o(1) and p∗  q ∗ . Define
t∗ , λ∗ and π̂ MF the same way as in Theorem 4.1. If nI/[k log kw] → ∞, we
have with probability at least 1 − e3 5−n ,
√
Z ∗ Z ∗T − π̂ MF (π̂ MF )T 1 . n2 / nI.
(ρ,ρ0 )

If we further assume Z ∗ ∈ Π0
probability at least 1 − e3 5−n ,

with arbitrary ρ, ρ0 , and then we have with

`(π̂ MF , Z ∗ ) . ρ−1 n

p
k 2 /(nI).

Proof. Form Lemma C.2, with probability at least 1 − e3 5−n , we have
uniformly for all π ∈ Π1
p
(51)
|hA − EA, ππ T i| ≤ 6n np∗ .

In the remaining part of the proof, we always assume
holds.
P the above event
pri
Denote f 0 (π) = hA+λ∗ In −λ∗ 1n 1Tn , ππ T i−(t∗ )−1 ni=1 KL(πi,· kπi,·
) for any

18

ZHANG, ZHOU

pri
pri
π ∈ Π1 . Here we adopt the notation KL(πi,· kπi,·
) short for KL(Categorical(πi,· )kCategorical(πi,·
)),
and we do it in the same way in the rest part of the proof. Thus,
p
hEA + λ∗ In − λ∗ 1n 1Tn , π̂ MF (π̂ MF )T i ≥ hA + λ∗ In − λ∗ 1n 1Tn , π̂ MF (π̂ MF )T i − 6n np∗
n
X
p
pri
0 MF
∗ −1
MF
∗
= f (π̂ ) − 6n np + (t )
KL(π̂i,·
kπi,·
)
i=1

n
X
p
pri
∗ −1
MF
0
∗
∗
KL(π̂i,·
kπi,·
)
≥ f (Z ) − 6n np + (t )
i=1
T
∗ ∗T
1n 1n , Z Z i

p
≥ hEA + λ In − λ
− 12n np∗
n
n
X
X
pri
pri
∗ −1
MF
∗ −1
∗
+ (t )
KL(π̂i,· kπi,· ) − (t )
KL(Zi,·
kπi,·
),
∗

∗

i=1

i=1

where we use Equation (51) twice in the first and last inequality. Note that
for any π ∈ Π1 , we have
X
X
pri
pri
|KL(πi,· kπi,·
)| ≤ |
πi,j log πi,j | + |
πi,j log πi,j
| ≤ log k + log w,
j

j

P
where the second inequality is due to 0 ≥ j πi,j log πi,j = KL(πi,· kk −1 1k )−
log k ≥ − log k, where k −1 1k can be explicitly written as a length-k vector
(1/k, 1/k, . . . , 1/k). Then we have
n
X
i=1

pri
MF
KL(π̂i,·
kπi,·
)

−

n
X
i=1

pri
∗
KL(Zi,·
kπi,·
) ≤ 2n log kw.

Thus,
hEA + λ∗ In − λ∗ 1n 1Tn , Z ∗ Z ∗T − π̂ MF (π̂ MF )T i ≤ 12n
By Proposition C.4, we have
∗

hEA + λ In − λ

∗

1n 1Tn , Z ∗ Z ∗T

− π̂

MF

(π̂

p
np∗ + 2(t∗ )−1 n log kw.




λ∗ − q ∗
λ∗ − q ∗
α+ ∗
γ ,
) i ≥ 2(p − q ) 1 − ∗
p − q∗
p − q∗

MF T

∗

∗

where α = hZ ∗ Z ∗T − π̂ MF (π̂ MF )T , Z ∗ Z ∗T − In i/2 and γ = hπ̂ MF (π̂ MF )T −
Z ∗ Z ∗T , 1n 1Tn − Z ∗ Z ∗T i/2. By Proposition C.3, there exists a constant c > 0
such that
(52)

hEA + λ∗ In − λ∗ 1n 1Tn , Z ∗ Z ∗T − π̂ MF (π̂ MF )T i ≥ 2c(p∗ − q ∗ )(α + γ).

MEAN FIELD FOR COMMUNITY DETECTION

19

Note that the following inequality holds
2(α + γ) = Z ∗ Z ∗T − π̂ MF (π̂ MF )T
≥ Z ∗ Z ∗T − π̂ MF (π̂ MF )T

1
1

− hZ ∗ Z ∗T − π̂ MF (π̂ MF )T , In i/2

− n/2.

These together lead to
Z ∗ Z ∗T − π̂ MF (π̂ MF )T

1

≤

i
h
p
1
∗ + 2(t∗ )−1 n log kw + c(p∗ − q ∗ )n/2 .
np
12n
c(p∗ − q ∗ )

Note that t∗  (p∗ − q ∗ )/p∗ when p∗  q ∗ . Together by Proposition C.2,
as long as nI/[k log kw] → ∞, the last two terms in the RHS of the above
formula is dominated by the first term. Thus,
Z ∗ Z ∗T − π̂ MF (π̂ MF )T
(ρ,ρ0 )

If we further assume Z ∗ ∈ Π0
to

1

n2
.√ .
nI

, Proposition C.5 and Equation (52) lead

hEA + λ∗ In − λ∗ 1n 1Tn , Z ∗ Z ∗T − π̂ MF (π̂ MF )T i ≥

ρcn(p∗ − q ∗ )
`(π̂ MF , Z ∗ ).
8k

So we have
p
8k
(12n np∗ + 2(t∗ )−1 n log kw)
∗
∗
ρcn(p − q )
s
192k
np∗
≤
.
∗
ρc
(p − q ∗ )2

`(π̂ MF , Z ∗ ) ≤

Before we state the remaining lemmas and propositions used in the Proof
of Lemma C.6, we first introduce two definitions. For any π, π 0 ∈ [0, 1]n×k ,
0 0
0 0
0 0
define α(π; π 0 ) = hπ π T −ππ T , π π T −In i/2 and γ(π; π 0 ) = hππ T −π π T , 1n 1Tn −
0 0T
π π i/2.
Proposition C.4. Define P = Z ∗ BZ ∗T − pIn , with B = q1k 1Tk + (p −
q)Ik . We have the equation



λ−q
λ−q
α(π; Z ∗ ) +
γ(π; Z ∗ ) .
hP + λIn − λ1n 1Tn , Z ∗ Z ∗T − ππ T i = 2(p − q) 1 −
p−q
p−q

20

ZHANG, ZHOU

Proof. Note that Z ∗ BZ ∗T − pIn = (p − q)Z ∗ Z ∗T + q1n 1Tn . We have
λ−q
λ−p
1n 1Tn +
In , Z ∗ Z ∗T − ππ T i
p−q
p−q
− In , Z ∗ Z ∗T − ππ T i

hP + λIn − λ1n 1Tn , Z ∗ Z ∗T − ππ T i = (p − q)hZ ∗ Z ∗T −
= (p − q)hZ ∗ Z ∗T

+ (λ − q)hIn − 1n 1Tn , Z ∗ Z ∗T − ππ T i

= (p − λ)hZ ∗ Z ∗T − In , Z ∗ Z ∗T − ππ T i

+ (λ − q)hZ ∗ Z ∗T − 1n 1Tn , Z ∗ Z ∗T − ππ T i

= 2(p − q)α(π; Z ∗ ) + 2(λ − q)γ(π; Z ∗ ).
Consequently, we obtain the desired bound.
(ρ,ρ0 )

If Z ∗ ∈ Π0

Proposition C.5.

, π ∈ Π1 , we have

α(π; Z ∗ ) + γ(π; Z ∗ ) ≥

ρn
`(π, Z ∗ ).
16k

Proof. We use α, γ instead of α(π; Z ∗ ), γ(π; Z ∗ ) for simplicity. Without
∗
∗ = 1}
loss of generality
Z ∗ ). Define Cu = {i : Zi,u
P we assume kπ − Z k1 = `(π, P
and Lu,v = i∈Cu πi,v . We have the equality v Lu,v = |Cu | and also
α=
and γ =

1X
2

u

1X
2



|Cu |2 −
X

X X

i,j∈Cu w

X

u6=v i∈Cu ,j∈Cv w



"
#
X
X
1
1X X
πi,w πj,w  =
|Cu |2 −
L2u,w =
Lu,w Lu,w0
2 u
2 u
0
w

πi,w πj,w =

1 XX
2

w6=w

Lu,w Lv,w .

u6=v w

We define [k] into two disjoint subsets S1 and S2 where
o
n
3
S1 = u ∈ [k] : ∀v 6= u, Lu,v ≤ |Cu | ,
4
n
o
3
and S2 = i ∈ [k] : ∃v 6= u, Lu,v > |Cu | .
4

P
Define L = v6=u Lu,v . For any u ∈ S1 , if Lu,u ≥ |Cu |/4, we have |Cu |2 −
P 2
P 2 u
1
2
w Lu,w ≥
w Lu,w ≥ Lu,u Lu ≥ |Cu |Lu /4. If Lu,u < 4 |Cu | we have |Cu | −
3
2
8 |Cu | ≥ |Cu |Lu /4 as well. This leads to
"
#
X
1 X
1 X
α≥
|Cu |2 −
L2u,w ≥
|Cu |Lu .
2
8
w
u∈S1

u∈S1

MEAN FIELD FOR COMMUNITY DETECTION

21

For any u ∈ S2 there exists a v 6= u such that Lu,v > 43 |Cu |. We must have
Lu,u + Lv,v ≥ Lu,v + Lv,u otherwise kπ − Z ∗ k1 = `(π, Z ∗ ) does not hold since
we can switch the u-th and v-th columns of π to make
− Z ∗ k1 smaller.
P kπP
Consequently, we have Lv,v ≥ Lu /2. So we have
u0 6=u
w Lu,w Lu0 ,w ≥
Lu,v Lv,v ≥ 3|Cu |Lu /8. Then we have
γ≥

1 X XX
3 X
Lu,w Lu0 ,w ≥
|Cu |Lu .
2
8
0
w
u∈S2 u 6=u

u∈S2

Thus,
α+γ ≥

ρn X
ρn
ρn
1 X
|Cu |Lu ≥
Lu ≥
kπ − Z ∗ k1 =
`(π, Z ∗ ).
16 u
16k u
16k
16k

C.3. Statements and Proofs of Lemmas and Propositions for
Theorem 4.2.
Lemma C.7. Let X ∼ Beta(α, β) where α = n2 p and β = n2 (1 − p) with
p = o(1). Let η = o(1). Then we have
P(|X − p| ≥ ηp) ≤ exp(−η 2 n2 p/2).
Proof. Note X has the same distribution as Y /(Y + Z) where Y and
Z are independent χ2 random variables with Y ∼ χ2 (2α) and Z ∼ χ2 (2β).
Then by using tail bound of χ2 distribution (i.e., Proposition C.6)
P(|X − p| ≥ ηp) ≤ P(|Y − 2n2 p| ≥ 2ηn2 p) + P(|Y + Z − 2n2 | ≥ ηn2 )
≤ 2 exp(−η 2 n2 p/4) + 2 exp(−η 2 n2 /16)
≤ exp(−η 2 n2 p/2).

Proposition C.6. Let X ∼ χ2 (k) we have


P |X − k| ≥ kt ≤ 2 exp(−kt2 /8), ∀t ∈ (0, 1).
Proof. See Lemma 1 of [20].

22

ZHANG, ZHOU

APPENDIX D: GENERAL DERIVATIONS OF CAVI FOR
VARIATIONAL INFERENCE
In this section, we provide the derivation from Equation (3) to Equation
(4). First we have
(53)

h
q(x) i
KL(q(x)kp(x|y)) = Eq(x) log
p(x|y)
= Eq(x) [log q(x)] − Eq(x) [log p(x|y)]

= Eq(x) [log q(x)] − Eq(x) [log p(x, y)] + log p(y)

= −(Eq(x) [log p(x, y)] − Eq(x) [log q(x)]) + log p(y)


= − Eq(x) [log p(y|x)] − KL(q(x)kp(x)) + log p(y).

Thus, to minimize KL(q(x)kp(x|y)) w.r.t. q(x) is equivalent to maximize
Eq(x) [log p(y|x)] − KL(q(x)kp(x)).
n
Recall we have independence under both p and
Q q for {xi }i=1 . For simplicity, denote x−i to be {xj }j6=i and q−i to be j6=i qj . We have the decomposition

bi (qi ) , Eq(x) [log p(x, y)] − Eq(x) [log q(x)]




= Eqi Eq−i [log p(xi , x−i , y)] − Eqi Eq−i [log q(xi , x−i )]


= Eqi Eq−i [log p(xi |x−i , y)] − Eqi [log qi (xi )] + const
log qi (xi )

 + const,
= −Eqi log −1
c exp Eq−i [log p(xi |x−i , y)]


P
where the constant includes all terms not depending on xi and c = xi exp Eq−i [log p(xi |x−i , y)]
which is also independent of xi . It is obvious that to solve Equation (3) is
equivalent to
q̂i = arg max bi (qi )
qi



= arg min KL qi kc−1 exp Eq−i [log p(xi |x−i , y)] .
qi



Immediately we have q̂i (xi ) = c−1 exp Eq−i [log p(xi |x−i , y)] . Or we may
write it as


q̂i (xi ) ∝ exp Eq−i [log p(xi |x−i , y)] .
Department of Statistics
Yale University
New Haven, CT 06511
E-mail: ye.zhang@yale.edu
E-mail: huibin.zhou@yale.edu
URL: http://www.stat.yale.edu/˜hz68/

