Theoretical and Computational Guarantees of Mean Field Variational Inference for Community Detection

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 is used to approximate posterior distributions in a wide range of applications in many fields including natural language processing , computational neuroscience , and network science . This method is different from Markov chain Monte Carlo (MCMC) , 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 , which is particularly attractive in the big data era.

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 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 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) , 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 under the Stochastic Block Model (SBM) 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 . For high dimensional settings, Celisse et al. and Bickel et al. 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. , 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 variational likelihood for the SBM, as seen from Theorem 2.1.

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 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.

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

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.

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

Despite the fact that every measure $\mathbf{q}$ has a simple product structure, the global minimizer $\mathbf{\hat{q}}^{\text{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 $\text{KL}(\mathbf{q}\|\mathbf{p})$ decreases in each coordinate update:

The coordinate update has an explicit formula

where $x_{-i}$ indicates all the coordinates in $x$ except $x_{i}$, and the expectation is over $\mathbf{q}_{-i}=\prod_{j\neq i}\mathbf{q}_{j}(x_{j})$. 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:

where $\mathbf{\hat{q}}^{\text{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. 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 . We include it in Appendix D in the supplemental material for completeness.

2 Community Detection and Stochastic Block Model

The Stochastic Block Model (SBM) has been a popular model for community detection.

where $B\in^{k\times k}$ with diagonal entries as $p$ and off-diagonal entries as $q$. That is, $B=q1_{k}1_{k}^{T}+(p-q)I_{k}$. Let $Z\in\Pi_{0}$ be the assignment matrix where

In each row $\{Z_{i,\cdot}\}_{i=1}^{n}$ there is only one 1 with all the other coordinates as 0, indicating the assignment of community for the corresponding node. Then $P$ can be equivalently written as $P_{i,j}=Z_{i,\cdot}BZ_{j,\cdot}^{T},\forall i<j$, or in a matrix form

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\in[k]^{n}$ and $Z\in\Pi_{0}$ which is defined as follows,

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.

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:

with $B=q1_{k}1_{k}^{T}+(p-q)I_{k}$ 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 $\{z_{i}\}_{i=1}^{n}$ have independent categorical (a.k.a. multinomial with size one) priors with hyperparameters $\{\pi^{\text{pri}}_{i,\cdot}\}_{i=1}^{n}$, where $\sum_{a=1}^{k}\pi^{\text{pri}}_{i,a}=1,\forall i\in[n]$. In other words, $\{Z_{i,\cdot}\}_{i=1}^{n}$ are independently distributed by

where $\{e_{a}\}_{a=1}^{k}$ are the coordinate vectors. Here we allow the priors for $Z_{i,\cdot}$ to be different for different $i$. If additionally $\pi_{i,\cdot}=\pi_{j,\cdot}$ for all $i\neq j$ is assumed, and then this is reduced to the usual case of i.i.d. priors.

Since $\{A_{i,j}\}_{i<j}$ are Bernoulli, it is natural to consider a conjugate Beta prior for $p$ and $q$. Let $p\sim\text{Beta}(\alpha^{\text{pri}}_{p},\beta^{\text{pri}}_{p})$ and $q\sim\text{Beta}(\alpha^{\text{pri}}_{q},\beta^{\text{pri}}_{q})$. Then the joint distribution is

Our main interest is to infer $Z$, from the posterior distribution $\mathbf{p}(Z,p,q|A)$. However, the exact calculation of $\mathbf{p}(Z,p,q|A)$ is computationally intractable.

4 Mean Field Approximation

Since the posterior distribution $\mathbf{p}(Z,p,q|A)$ is computationally intractable, we apply the mean field approximation to approximate it by a product measure,

where $\{r^{-1}(Z_{i,\cdot})\}_{i=1}^{n}$ are independent categorical variables with parameters $\{\pi_{i,\cdot}\}_{i=1}^{n}$, i.e., $\mathbf{q}_{\pi}(Z)=\prod_{i=1}^{n}\mathbf{q}_{\pi_{i,\cdot}}(Z_{i,\cdot})$ with

and $\mathbf{q}_{\alpha_{p},\beta_{p}}(p)$ and $\mathbf{q}_{\alpha_{q},\beta_{q}}(q)$ are Beta with parameters $\alpha_{p},\beta_{p},\alpha_{q},\beta_{q}$ due to conjugacy. See Figure 1 for the graphical presentation of $\mathbf{q}_{\pi,\alpha_{p},\beta_{p},\alpha_{q},\beta_{q}}(Z,p,q)$.

Note that the distribution class of $\mathbf{q}$ is fully captured by the parameters $(\pi,\alpha_{p},\beta_{p},\alpha_{q},\beta_{q})$, and then the optimization in Equation (2) is equivalent to minimize over the parameters as

The mean field estimator $(\hat{\pi}^{\text{MF}},\hat{\alpha}_{p}^{\text{MF}},\hat{\beta}_{p}^{\text{MF}},\hat{\alpha}_{q}^{\text{MF}},\hat{\beta}_{q}^{\text{MF}})$ defined in Equation (8) is equivalent to

The explicit formulation in Theorem 2.1 is helpful to understand the global minimizer of the mean field method. However, the global minimizer $\hat{\pi}^{\text{MF}}$ remains computationally infeasible as the objective function is not convex. Fortunately, there is a practically useful algorithm to approximate it.

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 $\{Z_{i,\cdot}\}_{i=1}^{n},p,q$. CAVI updates them one by one. Since the distribution class of $\mathbf{q}$ is uniquely determined by the parameters $\{\pi_{i,\cdot}\}_{i=1}^{n},\alpha_{p},\beta_{p},\alpha_{q},\beta_{q}$, equivalently we are updating those parameters iteratively. Theorem 2.2 gives explicit formulas for the coordinate updates.

Starts with some $\pi,\alpha_{p},\beta_{p},\alpha_{q},\beta_{q}$, the CAVI update for each coordinate (i.e., Equation (3) and Equation (4)) has an explicit expression as follows:

Update on $Z_{i,\cdot},\forall i=1,2,\ldots,n$:

where $t$ and $\lambda$ are defined in Equation (9) and Equation (10) respectively, and the normalization satisfies $\sum_{a=1}^{k}\pi^{\prime}_{i,a}=1$.

All coordinate updates in Theorem 2.2 have explicit formulas, which makes CAVI a computationally attractive way to approximate the global optimum $\mathbf{\hat{q}}^{\text{MF}}$ for the community detection problem.

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 $\pi$ sequentially one by one, while BCAVI uses the value of $\pi$ to update all rows $\{\pi^{\prime}_{i,\cdot}\}$ according to Theorem 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:\Pi_{1}\rightarrow\Pi_{1}$ as follows. For any $\pi\in\Pi_{1}$, we have

with parameters $t$ and $\lambda$. For BCAVI, we update $\pi$ by $\pi^{\prime}=h_{t,\lambda}(\pi)$ in each batch iteration, with $t,\lambda$ defined in Equations (14) and (15). See Algorithm 1 for the detailed implementation of BCAVI algorithm.

The definitions of $t^{(s)}$ and $\lambda^{(s)}$ in Equations (14) and (15) involve the digamma function, which costs a non-negligible computational resources each time called. Note that we have $\psi(x)\in(\log(x-\frac{1}{2}),\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

Later we show that $\alpha^{(s)}_{p},\beta^{(s)}_{p},\alpha^{(s)}_{q},\beta^{(s)}_{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).

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 $\pi$.

Note that the infimum over $\Phi$ 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^{\prime}=(2,2,1,1)$ give the same partition. In Equation (17) two equivalent assignments give the same loss.

2 Ground Truth

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.

It is worth mentioning that $\rho,\rho^{\prime}$ are not necessarily constants. We allow the community sizes not to be of the same order in the theoretical analysis.

3 Theoretical Justifications for BCAVI

In Theorem 3.1, we present theoretic guarantees of the convergence rate of BCAVI when initialized properly. Define

When $w=1$, the priors for $\{r^{-1}(Z_{i,\cdot})\}_{i=1}^{n}$ are i.i.d. $\text{Categorical}(1/k,1/k,\ldots,1/k)$ and $\bar{n}_{\text{min}}=n/2$ when there exist only two communities. The following quantity $I$ plays a key role in the minimax theory

which is the Rényi divergence of order $1/2$ between two Bernoulli distributions: $\text{Ber}(p^{*})$ and $\text{Ber}(q^{*})$. The proof of Theorem 3.1 is deferred to Section 5.3.

Let $Z^{*}\in\Pi_{0}$. Let $0<c_{0}<1$ be any constant. Assume $0<c_{0}p^{*}<q^{*}<p^{*}=o_{n}(1)$,

holds uniformly with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-c}-\epsilon$.

Theorem 3.1 establishes a linear convergence rate for BCAVI algorithm. The coefficient $[nI/[wk[n/\bar{n}_{\text{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.

Under the same condition as in Theorem 3.1, for any $s\geq s_{0}\triangleq[nI/k]/\log[nI/[wk[n/\bar{n}_{\text{min}}]^{2}]]$, we have

with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-c}-\epsilon$.

Under the assumption $nI/(k\log k)\rightarrow\infty$, we have

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 assume $\rho,\rho^{\prime}$ are constants, and the prior $\pi^{\text{pri}}_{i,a}\asymp 1/k,\forall i\in[n],a\in[k]$ (for example, uniform prior), and then the first condition in Equation (18) is equivalent to

noting that $n/\bar{n}_{\text{min}}\asymp k$ and $w\asymp 1$. This condition is necessary for consistent community detection when $k$ is finite. The assumptions in Equation (18) is slightly stronger than the assumption in , which is essentially $nI\geq Ck^{2}\log k$ for a sufficient large constant $C$.

Under the assumption $nI/k^{3}\rightarrow\infty$, since we have $I\asymp(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 $\alpha^{\text{pri}}_{p},\beta^{\text{pri}}_{p},\alpha^{\text{pri}}_{q},\beta^{\text{pri}}_{q}$ all equals to 1, i.e., consider non-informative priors.

Discussion

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 $\hat{\pi}^{\text{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 $\mathbf{p}(Z|A)$. We use a product measure $\mathbf{q}_{\pi}(Z)=\prod_{i}\mathbf{q}_{i}(\pi_{i,\cdot})$ for approximation, and then $\hat{\pi}^{\text{MF}}=\argmin_{\pi\in\Pi_{1}}\text{KL}[\mathbf{q}_{\pi}(Z)\|\mathbf{p}(Z|A)]$. Theorem 4.1 reveals that $\hat{\pi}^{\text{MF}}$ is rate-optimal, not surprisingly given the theoretical results obtained for BCAVI, an approximation of $\hat{\pi}^{\text{MF}}$.

Assume $p^{*}$ and $q^{*}$ are known. Under the assumption $\rho nI/[wk^{2}[n/\bar{n}_{\text{min}}]^{2}]\rightarrow\infty$, there exist some constant $c>0$ and $\eta=o_{n}(1)$ such that

with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-c}$.

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.

We present a batched version of Gibbs sampling for community detection. It involves iterative updates with

Generate $p^{(s)}$ by sampling from $\mathbf{p}(p|q^{(s-1)},Z^{(s-1)},A)$;

Generate $q^{(s)}$ by sampling from $\mathbf{p}(q|p^{(s-1)},Z^{(s-1)},A)$;

Generate $Z^{(s)}_{i,\cdot}$ independently by sampling from $\mathbf{p}(Z^{(s-1)}_{i,\cdot}|Z^{(s-1)}_{-i,\cdot},p^{(s)},q^{(s)},A)$, for $i\in[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.

holds with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}})]-n^{-c}-\epsilon$, where $b_{n}=\exp\left[-\eta^{\prime 2}\bar{n}_{\text{min}}^{2}\right]+\exp\left[-\eta^{\prime 2}n^{2}I\right]$ and $c_{n}=1/\sqrt{nI/[wk[n/\bar{n}_{\text{min}}]^{2}]}$. Consequently, for $s=[nI/k]/\log[nI/[wk[n/\bar{n}_{\text{min}}]^{2}]]$, we have

with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-c}-\epsilon$.

Theorem 4.2 establishes theoretical justification for batched Gibbs sampling for community detection. Despite that we have the same $c_{n}$ 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 $b_{n}$ is necessary to handle the extreme events due to random generation. Note that $(s+1)nb_{n}$ is dominated by $n\exp(-(1-\eta)\bar{n}_{\text{min}}I)$ as long as $s\leq e^{n}$. Thus, when $s\leq e^{n}$, we have similar “linear convergence” results as in Theorem 3.1.

3 An Iterative Algorithm for Maximum Likelihood Estimation

Maximum likelihood estimator (MLE) usually yields statistical optimality. However, the maximization of the likelihood $\mathbf{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\mathbf{p}(A|Z,p,q)$ by alternating maximization. We use a batched coordinate maximization:

Maximize $\mathbf{p}(A|p,q^{(s-1)},Z^{(s-1)})$ over $p$ to obtain $p^{(s)}$;

Maximize $\mathbf{p}(A|p^{(s-1)},q,Z^{(s-1)})$ over $q$ to obtain $q^{(s)}$;

Maximize $\mathbf{p}(A|p^{(s-1)},q^{(s-1)},Z_{i,\cdot},Z^{(s-1)}_{-i,\cdot})$ over $Z_{i,\cdot}$ to obtain $Z_{i,\cdot}^{(s)}$, for each $i\in[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.

holds with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-c}-\epsilon$.

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.

From Equation (8), by some algebra (see Equation (53) in Appendix D for detailed derivation) we have

where we use $\mathbf{q}$ instead of $\mathbf{q}_{\pi,\alpha_{p},\beta_{p},\alpha_{q},\beta_{q}}$ for simplicity. From the conditional distribution in Equation (6), the log-likelihood function can be simplified as

Due to the independence of $Z$ and $p,q$ under $\mathbf{q}$, we have

Since $B_{a,a}=p,\forall a\in[k]$ and $B_{a,b}=q,\forall a\neq b$, we have

By properties of Beta distribution, we obtain

where we use the fact that $\left\|{\pi_{i,\cdot}}\right\|_{1}=1,\forall i\in[n]$. Now consider the Kullback-–Leibler divergence between $\mathbf{q}(Z,p,q)$ and $\mathbf{p}(Z,p,q)$. Due to the independence of $p,q$ and $\{Z_{i,\cdot}\}_{i=1}^{n}$ in both distributions, we have

By Equations (19) - (23), we conclude with the desired result.

2 Proof of Theorem 2.2

We rewrite the joint distribution $\mathbf{p}(p,q,z,A)$ in Equation (7) as follows,

From Equation (24), $p$ has conditional probability as

Then the CAVI update in Equation (4) leads to

The distribution of $p$ is still Beta $p\sim\text{Beta}(\alpha^{\prime}_{p},\beta^{\prime}_{p})$, with

Similar analysis on $q$ yields updates on $\alpha^{\prime}_{q}$ and $\beta^{\prime}_{q}$. Hence, its proof is omitted.

Updates on {Zi,⋅}i=1n\{Z_{i,\cdot}\}_{i=1}^{n}

From Equation (24), the conditional distribution on $Z_{i,\cdot}$ is

Consequently, up to a constant not depending on $i$, we have

Then the CAVI update from Equation (4) leads to

where we use the property that $p,q,Z$ are all independent of each other under $\mathbf{q}$. Recall that $p\sim\text{Beta}(\alpha_{p},\beta_{p})$ and $q\sim\text{Beta}(\alpha_{q},\beta_{q})$. It can be shown that

3 Proof of Theorem 3.1

The proof of Theorem 3.1 involves three parts as follows.

Part One: One Iteration. Consider any $\pi\in\Pi_{1}$ such that $\left\|{\pi-Z^{*}}\right\|_{1}\leq\gamma\bar{n}_{\text{min}}$. Let $\eta^{\prime}$ be any sequence such that $\eta^{\prime}=o(1)$. Consider any $t$ and $\lambda$ with $|t-t^{*}|\leq\eta^{\prime}(p^{*}-q^{*})/p^{*}$ and $|\lambda-\lambda^{*}|\leq\eta^{\prime}(p^{*}-q^{*})$. We define $\mathcal{F}$ to be the event, that after applying the mapping $h_{t,\lambda}(\cdot)$, there exists some $\eta=o(1)$ such that

holds uniformly over all the eligible $\pi,t$ and $\lambda$. We have

for some constant $r>0$. We defer its proof to the later part of this section.

Part Two: Consistency of Model Parameters. Consider any $\pi\in\Pi_{1}$ such that $\left\|{\pi-Z^{*}}\right\|_{1}\leq\gamma\bar{n}_{\text{min}}$. Define

From Lemma C.1, we have a concentration of $t,\lambda$ towards $t^{*},\lambda^{*}$. That is, there exists some $\eta^{\prime}=o(1)$, such that with probability at least $1-e^{3}5^{-n}$, the following inequalities hold

Part Three: Multiple Iterations. Consider any $\pi\in\Pi_{1}$ such that $\left\|{\pi-Z^{*}}\right\|_{1}\leq\gamma\bar{n}_{\text{min}}$. Define $\alpha_{p},\beta_{p},\alpha_{q},\beta_{q},t,\lambda$ as Equations (26) - (29). A combination of results from Part One and Part Two immediately implies that

holds uniformly over all the eligible $\pi$ with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{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 $[h_{t,\lambda}(\pi)]_{i,\cdot}$ is a function of $\pi$ and $A_{i,\cdot}$. It can be decomposed into a summation of two terms, one only involves the ground truth $Z^{*}$ and the other involves the deviation $\pi-Z^{*}$. That is,

With a proper choice of $f_{\cdot,1}$ and $f_{\cdot,2}$, the first term on the RHS of Equation (31) leads to the minimax rate $n\exp(-(1-\eta)\bar{n}_{\text{min}}I)$. Up to a constant not dependent on $\pi,Z^{*}$ or $A$, the second term can be written as

Proof of Part One. Denote $z=r^{-1}(Z^{*})$. By the definition of $h_{t,\lambda}(\cdot)$ in Equation (11), we have

We choose some $m\rightarrow\infty$ slowly such that

where we use the fact that $\min_{a\neq b}(n_{a}+b_{b})/2\geq\bar{n}_{\text{min}}$.

The key to the rest of the analysis is to understand Equation (33) through the decomposition of the critical quantity $\sum_{j\neq i}(\pi_{j,a}-\pi_{j,b})(A_{i,j}-\lambda)$. We will show for any pair of $a,b\in[k]$ such that $a\neq b$, and any $i\in[n]$ such that $z_{i}=b$, it is equal to a summation of two terms: one only involves the ground truth $Z^{*}$, and the other involves the deviation $\pi-Z^{*}$. The former remains steady along iterations and contributes to the minimax rate, while the latter needs to be connected with the error $\left\|{\pi-Z^{*}}\right\|_{1}$.

Let $\theta_{a,b}$ be a vector of length $n$ such that $[\theta_{a,b}]_{j}=\pi_{j,a}-Z^{*}_{j,a}+Z^{*}_{j,b}-\pi_{j,b},\forall j\in[n]$. Then we have

With the help of Equation (34), Equation (33) can be written as

Equations (18) and (32) imply $\sum_{l=0}^{m-1}\exp\left[-l(n_{a}+n_{b})I/(2m)\right]\leq 2$. Thus, we have

In this way we turn $\left\|{h_{t,\lambda}(\pi)-Z^{*}}\right\|_{1}$ into calculations on $L_{1}^{\text{sum}}$ and $L_{2}^{\text{sum}}$, where the former only involves the ground truth $Z^{*}$ and the latter only involves the deviation $\pi-Z^{*}$.

We can obtain upper bounds on $L_{1}^{\text{sum}}$ and $L_{2}^{\text{sum}}$ as follows. Their proofs are deferred to the end of this section.

For $L_{1}^{\text{sum}}$, there exists a sequence $\eta^{\prime\prime}=o(1)$ such that with probability at least $1-\exp[-2(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]$, we have

For $L_{2}^{\text{sum}}$, there exist constants $c$ and $r$ such that with probability at least $1-n^{-r}-\exp(-5np^{*})$, we have

with probability at least $1-\exp[-2(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-r}-\exp(-5np^{*})$. By Propositions C.2 and C.3, we have $p^{*}t^{*2}\asymp I$. Then due to Equation (32), we have

Thus, with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-r}$, there exists some $\eta=o(1)$, such that

The proof for Part One is complete. The very last thing remained to be obtained is upper bounds on $L_{1}^{\text{sum}}$ and $L_{2}^{\text{sum}}$, i.e., Equations (35) and (36). Recall the definition of $\theta_{a,b}$. We have some properties on $\theta_{a,b}$ which will be useful in the analysis for $L_{1}^{\text{sum}}$ and $L_{2}^{\text{sum}}$: $\left\|{\theta_{a,b}}\right\|_{\infty}\leq 2$ and

1. Bounds on $L_{1}^{\text{sum}}$. By applying Markov inequality, we have

With the help of Proposition C.1, we have

We are going to show $-(1-\eta^{\prime\prime})\bar{n}_{\text{min}}I$ upper bounds terms in the exponent of RHS of Equation (39) by some $\eta^{\prime\prime}=o(1)$. We first present some properties of $\lambda^{*},t^{*}$ and $I$ that will be helpful:

Here Equations (40) and (41) are proved by Propositions C.2 and C.3 respectively. Equation (42) is due to $t^{*}\asymp\log(1+(p^{*}-q^{*})/q^{*})\asymp(p^{*}-q^{*})/p^{*}$ under the assumption that $p^{*},q^{*}=o(1)$, $p^{*}\asymp q^{*}$.

The first term in the exponent of Equation (39) is upper bounded by $-(1-7/(8m))\bar{n}_{\text{min}}I$ by the assumption $t^{*}/t=1+o(1)$. Since $|t^{*}(\lambda-\lambda^{*})|\leq\eta^{\prime}t^{*}(p^{*}-q^{*})$, by Equations (40) and (42) the second term is upper bounded by $\eta^{\prime}\bar{n}_{\text{min}}I$ up to a constant factor. For the last term in the exponent of Equation (39), since $|\lambda-\lambda^{*}|\leq\eta^{\prime}(p^{*}-q^{*})$ we have

where we use Equations (37) and (40) - (42).

As a consequence, there exists a sequence $\eta^{\prime\prime}=o(1)$ that goes to zero slower than $m^{-1},\gamma,\eta^{\prime}$, such that the summation of three terms in the exponent of the RHS of Equation (39) is upper bounded by $-(1-\eta^{\prime\prime})\bar{n}_{\text{min}}I$. Thus, Equation (39) can be written as

Since $\eta^{\prime\prime}$ goes to 0 slower than $m^{-1}$, we have $\eta^{\prime\prime}\geq m^{-1}\geq(\bar{n}_{\text{min}}I)^{\frac{1}{4}}$ by Equation (32). Then by applying Markov inequality, we have

That is, with probability at least $1-\exp[-2(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]$, Equation (35) holds.

2. Bounds on $L_{2}^{\text{sum}}$. Depending on whether the network is dense or sparse, we consider two scenarios.

Thus, with probability at least $1-n^{-r}$,

Define $L_{2,1}^{\text{sum}}\triangleq\sum_{a=1}^{k}\sum_{b\neq a}L_{2,1}(a,b)$. We have

with probability at least $1-n^{-1}-\exp(-5np^{*})$. By the bounds on $L_{1}^{\text{sum}}$ and $L_{2}^{\text{sum}}$, and due to $t/t^{*}=1+o(1)$, we obtain Equation (36).

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 , 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

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.2 An Iterative Algorithm for Maximum Likelihood Estimation

We first define a mapping $h^{\prime}:\Pi_{0}\rightarrow\Pi_{0}$ as follows

Here if the maximizer is not unique, we simply pick the smallest index.

B Proofs of Other Theorems

The proof of Equation (44) mainly follows the proof of Part One in Section 5.3. We have

Define $\theta_{a,b}$ the same way as in Section 5.3, and by the same argument, we have

From Lemma C.1, when $c_{\text{init}}$ is sufficiently small, with probability at least $1-e^{3}5^{-n}$ we have

Proposition C.3 shows that $\lambda^{*}\in(q^{*}+c(p^{*}-q^{*}),q^{*}+(1-c)(p^{*}-q^{*}))$ for some positive constant $0<c<1/2$. Therefore, when $c_{\text{init}}$ is sufficiently small, we have $\lambda\in(q^{*},p^{*})$. Thus,

where we use Equation (37). By Equations (40) - (42), it is smaller than $(n_{a}+n_{z_{i}})/(8t)$ when $c_{\text{init}}$ is sufficiently small. As a consequence, we have

By Equations (40) - (42) and (45), when $c_{\text{init}}$ is small enough, $t^{*}/t\leq 2$ and $t^{*}|\lambda-\lambda^{*}|\leq I/6$. Thus

Hence, with probability at least $1-\exp(-\bar{n}_{\text{min}}I/24)$,

For $L_{2}^{\text{sum}}$ we use the same argument as in Section 5.3 and obtain

with probability at least $1-n^{-r}-\exp(-5np^{*})$ for some constants $r,c_{1},c_{2}>0$. Recall that

Using the same argument as in Section 5.3, we conclude with

with probability at least $1-\exp(-\bar{n}_{\text{min}}I/10)-n^{-r}$.

B.2 Proof of Theorem 4.1

Define $t^{*}=\frac{1}{2}\log\frac{p^{*}(1-q^{*})}{q^{*}(1-p^{*})}$ and $\lambda^{*}=\frac{1}{2t^{*}}\log\frac{1-q^{*}}{1-p^{*}}$. By the same simplification we derive in Theorem 2.1, we have

Recall the definition of $h_{t,\lambda}(\cdot)$ as in Equation (11). A key observation is that $\hat{\pi}^{\text{MF}}=h_{t^{*},\lambda^{*}}(\hat{\pi}^{\text{MF}})$, otherwise if there exists some $i\in[n]$ such that $[h_{t^{*},\lambda^{*}}(\hat{\pi}^{\text{MF}})]_{i,\cdot}$ not equal to $\hat{\pi}^{\text{MF}}_{i,\cdot}$. This indicates the implementation of CAVI update on the $i$-th row of $\pi$ will make change, leading to the decrease of $f^{\prime}(\cdot;A)$. This contradicts with the fact that $\hat{\pi}^{\text{MF}}$ is the global minimizer.

The fixed-point property of $\hat{\pi}^{\text{MF}}$ is the key to our analysis. It involves three steps.

with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-r}$.

Step Three. Using the property that $h_{t^{*},\lambda^{*}}(\hat{\pi}^{\text{MF}})=\hat{\pi}^{\text{MF}}$, we have

holds with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}}]-n^{-r}$. Then we obtain the desired result by simple algebra.

B.3 Proof of Theorem 4.2

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

where we use the plain right arrow ($\rightarrow$) to indicate deterministic generation and the curved right arrow ($\leadsto$) to indicate random generation. Despite a slight abuse of notation, we define $\pi^{(0)}=Z^{(0)}$.

Let $\gamma=o(1)$ be any sequence goes to 0 when $n$ grows. We define a series of events as follows:

global event $\mathcal{G}$: Consider any $Z\in\Pi_{1}$ such that $\left\|{Z-Z^{*}}\right\|_{1}\leq\gamma\bar{n}_{\text{min}}$. Define

local events $\{\mathcal{H}_{1}^{(s)}\}_{s=1}^{S}$: We define $\mathcal{H}_{1}^{(s)}=\{\left\|{\pi^{(s)}-Z^{*}}\right\|_{1}\geq\gamma\bar{n}_{\text{min}}/2\}$.

local events $\{\mathcal{H}_{2}^{(s)}\}_{s=1}^{S}$: We define $\mathcal{H}_{2}^{(s)}=\{\left\|{Z^{(s)}-Z^{*}}\right\|_{1}\geq\gamma\bar{n}_{\text{min}}\}$. For the conditional probability, we have

Since $\left\|{\pi^{(s)}-Z^{*}}\right\|_{1}\leq\gamma\bar{n}_{\text{min}}/2$ given $\mathcal{H}_{1}^{(s)}=0$ by Bernstein inequality, we have

local events $\{\mathcal{H}_{3}^{(s)}\}_{s=1}^{S}$: We define $\mathcal{H}_{3}^{(s)}=\{|t^{(s)}-t^{*}|\geq\eta^{\prime}(p^{*}-q^{*})/p^{*},\text{ or }|\lambda^{(s)}-\lambda^{*}|\geq\eta^{\prime}(p^{*}-q^{*})\}$. If the global event $\mathcal{G}$ holds and the local event $\mathcal{H}_{2}^{(s)}$ does not hold, we have

Note that $\alpha_{p}^{(s+1)}+\beta_{p}^{(s+1)}=\alpha^{\text{pri}}_{p}+\beta^{\text{pri}}_{p}+\sum_{a=1}^{k}\sum_{i<j}Z_{i,a}^{(s)}Z_{j,a}^{(s)}\geq n^{2}/k$. Using the tail bound of Beta distribution (Lemma C.7) we are able to show

where the last inequality is due to Proposition C.2. This leads to

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^{*}|\}\leq 2\eta^{\prime\prime}(p^{*}-q^{*})$ leads to

By taking $\eta^{\prime}=16c_{0}\eta^{\prime\prime}$, we obtain

Note that events $\mathcal{F}$ and $\mathcal{G}$ are about the adjacency matrix $A$. The events $\mathcal{H}_{1}^{(s)},\mathcal{H}_{2}^{(s)}$ and $\mathcal{H}_{3}^{(s+1)}$ are for $\pi^{(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). Under the event $\mathcal{F}\cap\mathcal{G}\cap(\mathcal{H}_{1}^{(s)}\cup\mathcal{H}_{2}^{(s)}\cup\mathcal{H}_{3}^{(s+1)})^{C}$ we have

where $c_{n}=[nI/[wk[n/\bar{n}_{\text{min}}]^{2}]]^{-1/2}$. As a consequence, under the event $\mathcal{F}\cap\mathcal{G}\cap(\prod_{v=0}^{s}\mathcal{H}_{1}^{(v)}\cup\mathcal{H}_{2}^{(v)}\cup\mathcal{H}_{3}^{(v+1)})^{C}$, we have

Due to the small value of $c_{n}$, if $\left\|{\pi^{(s)}-Z^{*}}\right\|_{1}\leq\gamma\bar{n}_{\text{min}}$, Equation (47) immediately implies $\left\|{\pi^{(s+1)}-Z^{*}}\right\|_{1}\leq\gamma\bar{n}_{\text{min}}$. This implies that under the event $\mathcal{F}\cup\mathcal{G}$ we have

with probability at least $1-\exp[-(\bar{n}_{\text{min}}I)^{\frac{1}{2}})]-n^{-r}-e^{3}5^{-n}-\epsilon$, where $b_{n}=\exp\left[-3(\gamma\bar{n}_{\text{min}})^{2}/16\right]+2\exp\left[-\eta^{\prime\prime 2}n^{2}I/2\right]$.

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.

C Statements and Proofs of Auxiliary Lemmas and Propositions

We include all the auxiliary propositions and lemmas in this section.

Let $c_{\text{init}}$ be some sufficiently small constant. Consider any $\pi\in\Pi_{1}$ such that $\left\|{\pi-Z^{*}}\right\|_{1}\leq c_{\text{init}}n/k$. Let $\alpha_{p},\beta_{p},\alpha_{q},\beta_{q},t,\lambda$ be the outputs after one step CAVI iteration from $\pi$ described in Algorithm 1. That is, they are defined as Equations (26) - (29). Define

Under the same assumption as in Theorem 3.1, there exists some sequence $\epsilon=o(1)$ such that with probability at least $1-e^{3}5^{-n}$, the following inequality holds

uniformly over all the eligible $\pi$. In addition if we further assume $c_{\text{init}}$ goes to 0, the LHS of the above inequality will be simply upper bounded by $\epsilon$.

We are going to obtain tight bounds on $|\hat{p}-p^{*}|$ and $|\hat{q}-q^{*}|$ first. Note that we have the “variance-bias” decomposition as in

We have concentration inequality holds for the numerator in the first term by Lemma C.2. That is, with probability at least $1-e^{3}5^{-n}$, we have

holds uniformly over all $\pi\in\Pi_{1}$. For the denominator, we have

since $\sum_{a=1}^{k}\left\|{\pi_{\cdot,a}}\right\|_{1}=n$. Thus, we are able to obtain an upper bound on the first term as

where in the last inequality we use the orthogonality between $Z^{*}Z^{*T}$ and $11^{T}-Z^{*}Z^{*T}$. For its numerator, we have

Similar result holds for $|\hat{q}-q^{*}|$. Denote $\eta_{0}=\sqrt{\frac{k^{2}p^{*}}{n(p^{*}-q^{*})^{2}}}+\frac{3\left\|{\pi-Z^{*}}\right\|_{1}}{n/k}$, thus

By the assumption of $nI$ in Equation (18) and Proposition C.2, we have $n(p^{*}-q^{*})^{2}/(k^{2}p^{*})\asymp nI/k^{2}\rightarrow\infty$. Therefore, the first term in $\eta_{0}$ goes to 0. The second term in $\eta_{0}$ is at most $3c_{\text{init}}$ which implies $\eta_{0}\leq 4c_{\text{init}}$.

By the fact that the digamma function satisfies $\psi(x)\in(\log(x-1/2),\log x),\forall x\geq 1/2$, we have

Recall that we have shown $\sum_{i<j}\sum_{a=1}^{k}\pi_{i,a}\pi_{j,a}$ lies in the interval of $(n^{2}/(2k),n^{2}/2)$. By Equation (18), there exists a sequence $\eta^{\prime}=o(1)$ such that $\alpha_{p},\beta_{p}\leq\eta^{\prime}(p^{*}-q^{*})n^{2}/k$. Then we have

Recall that we assume $c_{0}p^{*}<q^{*}<p^{*}$. Thus $(\eta_{0}+\eta^{\prime})(p^{*}-q^{*})/p^{*}\leq 5c_{\text{init}}c_{0}$. When $c_{\text{init}}$ is sufficiently small, we have $(\eta_{0}+\eta^{\prime})(p^{*}-q^{*})/p^{*}\leq 1/2$. Then using the fact $-x\geq\log(1-x)\geq-2x,\forall x\in(0,1/2)$. We have

Analogously we can obtain the same upper bound on $\hat{t}-t^{*}$, and then

Identical analysis can be applied towards bounds on $|\hat{\lambda}-\lambda^{*}|$. Note that

similarly for $\alpha_{q},\beta_{q}$. Omitting the immediate steps, we end up with

The proof is complete after we unify and rephrase all the aforementioned results. ∎

Let $A\in^{n\times n}$ such that $A=A^{T}$ and $A_{i,i}=0,\forall i\in[n]$. Assume $\{A_{i,j}\}_{i<j}$ are independent random variable, and there exists $p\leq 1$ such that $9n^{-1}\leq\frac{2}{n(n-1)}\sum_{i<j}\text{Var}(A_{i,j})\leq p$, and then we have

with probability at least $1-e^{3}5^{-n}$.

Then by applying Grothendieck inequality we obtain

where $c$ is a positive constant smaller than 2. This concludes with

Assume $0<q<p<1$. Let $X\sim\text{Ber}(q)$ and $Y\sim\text{Ber}(p)$. Recall the definition $\lambda=\log\frac{1-q}{1-p}/\log\frac{p(1-q)}{q(1-p)}$, $t=\frac{1}{2}\log\frac{p(1-q)}{q(1-p)}$ and $I=-2\log[\sqrt{pq}+\sqrt{(1-p)(1-q)}]$. Then the following two equations hold

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

The following lemma on the operator norm of sparse networks is from . In the original statement of Lemma 12 in , “with probability $1-o(1)$” is stated. However, its proof in gives explicit form of the probability that the statement holds, which is at least $1-n^{-1}$.

[Lemma 12 of ] Suppose $M$ is random symmetric matrix with zero on the diagonal whose entries above the diagonal are independent with the following distribution

holds with probability at least $1-n^{-1}$.

by implementing Bernstein inequality. Applying Bernstein inequality again we have

Under the assumption that $0<q<p=o(1)$. For $I=-2\log\left[\sqrt{pq}+\sqrt{(1-p)(1-q)}\right]$ we have

Consequently, $(p-q)^{2}/(4p)\leq I\leq(p-q)^{2}/p$.

It is a partial result of Lemma B.1 in . ∎

Define $\lambda=\log\frac{1-q}{1-p}/\log\frac{p(1-q)}{q(1-p)}$. For any $p,q>0$ such that $p,q=o(1)$ and $p\asymp q$, there exists a constant $0<c<1/2$ such that