# Community detection for binary graphical models in high dimension

Julien Chevallier, Guilherme Ost

<figure><img src="image_1.png" /></figure>

## To cite this version:

Julien Chevallier, Guilherme Ost. Community detection for binary graphical models in high dimension. 2024. hal-04803948

HAL Id: hal-04803948

https://hal.science/hal-04803948v1

Preprint submitted on 26 Nov 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
COMMUNITY DETECTION FOR BINARY GRAPHICAL MODELS IN HIGH DIMENSION

BY JULIEN CHEVALLIER¹,ª ID, GUILHERME OST²,ᵇ ID

¹Univ. Grenoble Alpes, CNRS, Inria, Grenoble INP*, LJK, 38000 Grenoble, France,
ªjulien.chevallier1@univ-grenoble-alpes.fr

²Institute of Mathematics, Federal University of Rio de Janeiro, ᵇguilhermeost@im.ufrj.br

Let $N$ components be partitioned into two communities, denoted $\mathcal{P}_+$ and $\mathcal{P}_-$, possibly of different sizes. Assume that they are connected via a directed and weighted Erdös-Rényi random graph (DWER) with unknown parameter $p \in (0, 1)$. The weights assigned to the existing connections are of mean-field type, scaling as $N^{-1}$. At each time unit, we observe the state of each component: either it sends some signal to its successors (in the directed graph) or remain silent otherwise. In this paper, we show that it is possible to find the communities $\mathcal{P}_+$ and $\mathcal{P}_-$ based only on the activity of the $N$ components observed over $T$ time units. More specifically, we propose a simple algorithm for which the probability of exact recovery converges to 1 as long as $(N/T^{1/2}) \log(NT) \to 0$, as $T$ and $N$ diverge. Interestingly, this simple algorithm does not required any prior knowledge on the other model parameters (e.g. the edge probability $p$). The key step in our analysis is to derive an asymptotic approximation of the one unit time-lagged covariance matrix associated to the states of the $N$ components, as $N$ diverges. This asymptotic approximation relies on the study of the behavior of the solutions of a matrix equation of Stein type satisfied by the simultaneous (0-lagged) covariance matrix associated to the states of the components. This study is challenging, specially because the simultaneous covariance matrix is random since it depends on underlying DWER random graph.

**1. Introduction.** In this paper, we consider a system of $N$ interacting $\{0, 1\}$-valued chains $X = \{X_{i,t}, t \in \mathbb{Z}, 1 \le i \le N\}$ that are coupled via a random matrix $\theta = (\theta_{ij})_{1 \le i,j \le N}$ of i.i.d entries distributed as $\text{Ber}(p)$ with $0 \le p \le 1$. Let $X_t = (X_{1,t}, \dots, X_{N,t})$ denote the configuration of the system $X$ at time $t$. The dynamics of the model is as follows. Conditionally on $\theta$, the system $X$ evolves as a stationary Markov chain on the state space $\{0, 1\}^N$ in which the conditional distribution of $X_t$ given that $X_{t-1} = x$ is that of $N$ independent $\text{Ber}(p_{\theta,i}(x))$ random variables, $i = 1, \dots, N$, where

$$ (1) \quad p_{\theta,i}(x) = \mu + (1 - \lambda) \left( \frac{1}{N} \sum_{j \in \mathcal{P}_+} \theta_{ij} x_j + \frac{1}{N} \sum_{j \in \mathcal{P}_-} \theta_{ij} (1 - x_j) \right), \quad x = (x_1, \dots, x_N). $$

The sets $\mathcal{P}_+$ and $\mathcal{P}_-$ form a partition of $[N] := \{1, \dots, N\}$ and the values $0 < \lambda < 1$ and $0 \le \mu \le \lambda$ are fixed. Neither the sets $\mathcal{P}_+$ and $\mathcal{P}_-$ nor the values $\mu$ and $\lambda$ depend on the random environment $\theta$.

In this model, the transition probabilities $p_{i,\theta}(x)$ defined in (1) are increasing with respect to $x_j$ for $j \in \mathcal{P}_+$ and decreasing for $j \in \mathcal{P}_-$, so that the components $j \in \mathcal{P}_+$ have an excitatory role, whereas those belonging to $\mathcal{P}_-$ have an inhibitory role. In that respect, the sets $\mathcal{P}_+$ and

<footnote>*Institute of Engineering Univ. Grenoble Alpes</footnote>

<footnote>MSC2020 subject classifications: Primary 62J05; secondary 60J10, 60K35, 62F12.</footnote>

<footnote>Keywords and phrases: Community detection, exact recovery, graphical models.</footnote>
$\mathcal{P}_-$ partitioning the set of components $[N]$ of the system $X$ introduce a community structure into the model: components in the same partition element have the same role. Hereafter, for $a \in \{-, +\}$, we denote $r_a^N = |\mathcal{P}_a|/N$ the fraction of components of the system that belong to the community $\mathcal{P}_a$, and assume that there are values $0 < r_+, r_- < 1$ satisfying $r_+ + r_- = 1$ such that $|r_+^N - r_+| \vee |r_-^N - r_-| \le KN^{-1}$ for some constant $K$. In particular, we do not assume that the communities are balanced, i.e., the choice $r_+ \neq r_-$ is allowed.

**Problem formulation.** The goal of this paper is to address the following community detection problem. By observing a sample $X_1, \dots, X_T$ of the system $X$, can we find the communities $\mathcal{P}_+$ and $\mathcal{P}_-$? We do not assume any prior knowledge on $\mu, \lambda, p, r_+, r_-$ and $\theta$.

**1.1. Related literature.** The model considered here was recently introduced in (Chevalier, Löcherbach and Ost, 2024). In their main result, the authors show that it is possible to estimate the connectivity parameter $p$ from the sample $X_1, \dots, X_T$ of the system $X$ with rate $N^{-1/2} + N^{1/2}/T + (\log(T)/T)^{1/2}$, as long as the asymptotic fraction of components in each community (the values $r_+$ and $r_-$) is known a priori. As a consequence, the proposed estimator is shown to be consistent provided that $N^{1/2}/T \to 0$. In this paper, we complement these results, by showing that one can also find the communities defined by the sets $\mathcal{P}_+$ and $\mathcal{P}_-$. More specifically, we propose a simple algorithm for which we can show (see Corollary 2.5) that the probability of *exact recovery* converges to 1 as long as $(N/T^{1/2}) \log(NT) \to 0$. Interestingly, this simple algorithm does not require the knowledge of the asymptotic fraction of components that belong to each community.

Community detection problems have been extensively investigated over the past decades. A large part of the known results were established for the stochastic block model (SBM) - also known as planted partition model - (Abbe, Bandeira and Hall, 2016; Mossel, Neeman and Sly, 2015, 2016) and its generalizations such as the degree-corrected SBM (Zhao, Levina and Zhu, 2012; Lei and Rinaldo, 2015) and mixed-membership SBM (Airoldi et al., 2008). We refer the reader to (Abbe, 2018) for a recent review of the classical results in this framework. The problem of detecting the communities in such stochastic block models is different from ours for at least two reasons. First, for such models, one has access to the adjacency matrix with a block structure in terms of the edge probabilities and uses a function of this matrix to find the communities. In our framework, the adjacency matrix counterpart is the random matrix $A^N$ defined in (3). This matrix is a normalized and signed version of the random matrix $\theta$ in which all columns in $\mathcal{P}_-$ have a negative sign. Since the matrix $A^N$ depends on the matrix $\theta$ that is not observed, one cannot use the matrix $A^N$ in our detection problem. Instead, we can only use the information present in the sample $X_1, \dots, X_T$ which only gives us a skewed view of $A^N$ (this is better justified in (4)). Second, given the community assignments, there is no interaction between the components of those models: they behave independently of each other. Here, on the contrary, the components of the model do interact and these interactions are essentials for finding the communities.

More recently, the problem of community detection have been also addressed in the context of graphical models (Berthet, Rigollet and Srivastava, 2019) and G-latent models (Bunea et al., 2020). In both these papers, the observations are assumed to be i.i.d. copies of some $N$-dimensional random vector $Z$. In (Berthet, Rigollet and Srivastava, 2019), the authors suppose that the distribution of $Z$ is that of an Ising model without external field in which the interaction strength is $\beta N^{-1}$ between sites of the same block and $\alpha N^{-1}$ otherwise, where the real values $\alpha$ and $\beta$ satisfy $\beta > \alpha$. In particular, the block structure in this model is imposed directly on the interactions, similar to our model. In the framework of G-latent models considered in (Bunea et al., 2020), the community structure is imposed directly on the covariance matrix of the random vector $Z$, unlike the model considered in (Berthet, Rigollet and Srivastava, 2019) and ours. In the present paper, we go beyond the assumption that the observations
are i.i.d by addressing the community detection problem for dependent variables. In this setting, part of the literature has focused on dynamic SBM (Pensky and Zhang, 2019; Lin and Lei, 2024; Matias and Miele, 2016) - a framework fundamentally different from ours - where one observes sequences of adjacency matrices that slowly vary their community structures over time and the goal is to estimate these time varying communities. Closer to our setting, we are only aware of the works of (Peixoto, 2019) and (Hoffmann et al., 2020). Both papers adopt a Bayesian approach for detecting communities from time series data. Their focus is on discussing the properties of the proposed methods using synthetic and empirical data, rather than developing theoretical guarantees leading to exact recovery. Establishing such results in a similar context is the goal of our paper.

1.2. Our contribution. The key step in our analysis is to study the asymptotic behavior of the one unit time-lagged covariance matrix $\Sigma^{(1)}$ whose $(i, j)$-entries are defined as

$$ (2) \qquad \Sigma_{ij}^{(1)} = \mathrm{Cov}_{\theta}(X_{i,1}, X_{j,0}), \quad i, j \in [N], $$

where $\mathrm{Cov}_{\theta}(X_{i,1}, X_{j,0})$ denotes the covariance between the random variables $X_{i,1}$ and $X_{j,0}$, conditionally on $\theta$. Introducing the random matrix $A^N$ whose $(i, j)$-entries are given by

$$ (3) \qquad A_{ij}^N = \begin{cases} N^{-1}\theta_{ij}, & \text{if } j \in \mathcal{P}_+, \\ -N^{-1}\theta_{ij}, & \text{if } j \in \mathcal{P}_-, \end{cases} \quad i, j \in [N], $$

we establish in Theorem 2.1 that there are explicitly known constants $c_1 > 0$ and $c_2$ (both depending on $\lambda, p$, and $r_+$) such that for all entries $(i, j)$,

$$ (4) \qquad \Sigma_{ij}^{(1)} = c_1 A_{ij}^N + c_2 N^{-1} + \epsilon_{ij}^N, $$

where the error terms $\epsilon_{ij}^N$ are shown to have the first moment uniformly bounded by $K\sqrt{\log(N)/N^3}$ with $K$ being some positive constant depending only on $\lambda, p$ and $r_+$. The term $c_2N^{-1}$ can be interpreted as a bias term present in all the entries of the matrix $\Sigma^{(1)}$. It is worth mentioning that the sign of this bias term $c_2N^{-1}$ is the same as the sign of the $r_+ - r_-$. In particular, the bias is 0 if and only if the communities are balanced, i.e. $r_+ = r_- = 1/2$.

The matrix $A^N$ defined in (3) is a normalized and signed version of $\theta$ in which the entries associated with columns in $\mathcal{P}_-$ have a negative sign. Hence, the columns of $A^N$ carry information about the communities $\mathcal{P}_+$ and $\mathcal{P}_-$. The approximation (4) suggests that the same is true for the columns of the random matrix $\Sigma^{(1)}$ for large values of $N$. This motivates us to consider the vector $\sigma^N$ defined by, for all $j \in [N]$,

$$ \sigma_j^N = \sum_{i=1}^{N} \Sigma_{ij}^{(1)}. $$

These quantities are asymptotically of order one and allow us to recover the communities $\mathcal{P}_+$ and $\mathcal{P}_-$. Indeed, starting from (4) and using the fact that $N^{-1}\sum_{i=1}^N \theta_{ij} \to p$ as $N \to \infty$ (recall that the $\theta_{ij}$'s are i.i.d. Ber($p$)), we get

$$ (5) \qquad \sigma_j^N = \begin{cases} c_1 p + c_2 + \sum_{i=1}^N \epsilon_{ij}^N & \text{if } j \in \mathcal{P}_+, \\ -c_1 p + c_2 + \sum_{i=1}^N \epsilon_{ij}^N & \text{if } j \in \mathcal{P}_-. \end{cases} $$

Hence, if the aggregated errors $\sum_{i=1}^N \epsilon_{ij}^N$ are small enough (in fact, we can show that their first moment are uniformly bounded by $K\sqrt{\log(N)/N}$ using the a priori control of the error terms), a value of $\sigma_j^N$ close to $c_1p + c_2$ (respectively $-c_1p + c_2$) indicates that $j \in \mathcal{P}_+$ (resp. $j \in \mathcal{P}_-$). Of course, one cannot use $\sigma^N$ in practice, but there is a natural and convergent
estimator of $\sigma^N$, denoted $\hat{\sigma}^N$ and defined in (11). The proof of its convergence uses (Ost and Reynaud-Bouret, 2020, Theorem 3). Using $\hat{\sigma}^N$, exact recovery is feasible with high probability under the asymptotic regime of Corollary 2.5.

The asymptotic approximation (4) is achieved by first showing that the one unit time-lagged covariance matrix $\Sigma^{(1)}$ can be computed from the simultaneous (0-lagged) covariance matrix $\Sigma^{(0)}$ whose $(i, j)$-entries are defined as

$$ (6) \qquad \Sigma_{ij}^{(0)} = \mathrm{Cov}_{\theta}(X_{i,0}, X_{j,0}), \quad i, j \in [N]. $$

Next, we show that simultaneous covariance matrix $\Sigma^{(0)}$ satisfies some Stein matrix equation (see Proposition 6.1). To conclude, we then study the asymptotic behavior of the solutions of the matrix equation satisfied by $\Sigma^{(0)}$ as $N \to \infty$. This analysis is challenging, specially because the matrix $\Sigma^{(0)}$ depends on random matrix $\theta$ and, hence, it is itself a random matrix. The analysis involves two main ingredients: identifying an "almost" eigenvector of the Kroenecker product $A^N \otimes A^N$, and a control over empirical means based on the matrix $\theta$ which uses the maximal inequality (Boucheron, Lugosi and Massart, 2013, Theorem 2.5).

Part of our motivation comes from Neuroscience. Nowadays, the simultaneous activity of large neuronal networks can be routinely recorded and understanding which structural properties of these networks can be inferred from the recordings is an important question. In several cases, neuronal recordings have been modeled by discrete-time $\{0, 1\}$-valued stochastic chains that share some similarities with our model (see for example (Galves and Löcherbach, 2013; Duarte et al., 2019) and references therein). In these cases, the time interval in which the neuronal network has been recorded is binned in small intervals and in each of these intervals, and for each recorded neuron, we assign the symbol 1 if the neuron has fired an action potential in that time interval and the symbol 0 otherwise. In this context, our results show that in an idealized network of neurons in which the dynamics is described by our model, one can discriminate exactly the excitatory neurons from the inhibitory ones using only the recorded activity.

The paper is organized as follows. In the Section 2 below, we define our model rigorously and introduce the minimal set of notation in order to state our main results. Section 3 introduces more notation and results already used and proved in Chevallier, Löcherbach and Ost (2024). Section 4 contains the main asymptotic approximations (especially Proposition 4.2) used to prove our main results. Section 5 is devoted to the proofs of our main results. Section 6 is devoted to the proof of Proposition 4.2. In Section 7, we illustrate the performance of our estimator through some simulations. Finally, Section 8 contains technical lemmas based on the maximal inequality.

## 2. Model definition, notation and main results.

*Model definition.* We consider a system of $N$ interacting chains $X = \{X_{i,t}, t \in \mathbb{Z}, 1 \le i \le N\}$ taking values in $\{0, 1\}$ denoting the presence or the absence of a signal at a given time. This system evolves in a random environment which is given by the realization of a directed Erdös-Rényi random graph via the selection of $N^2$ i.i.d. random variables $\theta_{ij}, 1 \le i, j \le N$, distributed as a Bernoulli of parameter $p$ with $0 \le p \le 1$. Conditionally on $\theta$, the evolution of the system $X$ is that of a stationary Markov chain on the hypercube $\{0, 1\}^N$ with transition probabilities (which depend on $\theta$) given as follows. Writing $X_t = (X_{1,t}, \dots, X_{N,t})$ and $x, y \in \{0, 1\}^N$, we have, for all $t \in \mathbb{Z}$,

$$ (7) \qquad I_{\theta}(X_t = y | X_{t-1} = x) = \prod_{i=1}^{N} (p_{\theta,i}(x))^{y_i} (1 - p_{\theta,i}(x))^{(1-y_i)}, $$
where $p_{\theta,i}(x) = \mathbb{P}_{\theta}(X_{i,t} = 1 | X_{t-1} = x)$ is defined in (1). The existence and uniqueness of a stationary version of the Markov chain having transition probabilities as defined in (7) and (1) is granted by Theorem 3.3 of (Chevallier, Löcherbach and Ost, 2024). Only the stationary regime is considered throughout this paper.

**Notation.** Throughout the paper, we use the letters $t, s$ to denote some time indices whereas the letters $i, j$ denote some spatial values, that is the index of the corresponding component. For a vector $v$ in $\mathbb{R}^N$, the notation $\bar{v}$ denotes the spatial mean $\bar{v} = N^{-1} \sum_{i=1}^N v_i$. In agreement with the left-hand side of (7), we write $\mathbb{P}_{\theta}$ to denote the probability measure under which the environment $\theta$ is kept fixed and the process $X$ is distributed as the unique stationary version of the Markov chain having transition probabilities as defined in (7) and (1). The expectation taken with respect to $\mathbb{P}_{\theta}$ is denoted by $\mathbb{E}_{\theta}$. The variance and covariance computed from $\mathbb{E}_{\theta}$ are denoted $\text{Var}_{\theta}$ and $\text{Cov}_{\theta}$ respectively. Moreover, we write $\mathbb{P}$ to denote the probability measure under which the random environment $\theta = (\theta_{ij})_{i,j \in [N]}$ is distributed as a collection of i.i.d. random variables with distribution $\text{Ber}(p)$ and the conditional distribution of the process $X$ given $\theta$ is that of the process $X$ under $\mathbb{P}_{\theta}$, i.e., the following identity holds $\mathbb{P}(X \in \cdot | \theta) = \mathbb{P}_{\theta}(X \in \cdot)$. Finally, we denote $\mathbb{E}$ the expectation taken with respect to the probability measure $\mathbb{P}$, and $\text{Var}$ and $\text{Cov}$ the variance and covariance, respectively, computed from $\mathbb{E}$. The letter $K$ denotes a constant that may change from line to line.

For any $N$-by-$N$ matrix $M$, we shall denote $\|M\|_p$, with $1 \le p \le \infty$, the matrix norm induced by the vector norm $\|\cdot\|_p$, and $\|M\|_{\max} = \max_{i,j \in [N]} |M(i,j)|$. Recall that, although the matrix norm $\|\cdot\|_{\max}$ is not sub-multiplicative, it satisfies

$$ (8) \qquad \|AB\|_{\max} \le \min \{ \|A\|_{\infty} \|B\|_{\max}, \|A\|_{\max} \|B\|_{1} \}, $$

for any $N$-by-$N$ matrices $A$ and $B$. We write $I_N$ and $J_N$ to denote, respectively, the identity matrix and the matrix full of ones, both $N$-by-$N$.

Hereafter, for any subset $V$ of $[N]$, we write $1_V$ to indicate the $N$-dimensional vector having value 1 in each coordinate belonging to $V$ and value 0 in the remaining coordinates. To alleviate the notation, we will simply write $1_N$ instead of $1_{[N]}$, to indicate the vector full of ones.

**Main results.** Let us denote

$$ (9) \qquad m = m(\mu, \lambda, p) = \frac{\mu + (1 - \lambda)pr_{-}}{1 - p(1 - \lambda)(r_{+} - r_{-})}, $$

which was proved in Theorem 2.1 of (Chevallier, Löcherbach and Ost, 2024) to be the limit, in probability, of the spatio-temporal mean $(TN)^{-1} \sum_{i=1}^N \sum_{t=1}^T X_{i,t}$ as $\min\{T, N\} \to \infty$.

Recall the random matrices $\Sigma^{(1)}$ and $A^N$ defined in (2) and (3) respectively. The matrix $\Sigma^{(1)}$ is the one unit time-lagged covariance matrix and $A^N$ is a normalized and signed version of random matrix $\theta$. Remark that $\Sigma^{(1)}$ also depends on $N$ but we chose to not highlight this dependence in the notation. Moreover, note that this matrix is not symmetric. Our first main result establishes a relation between the matrices $A^N$ and $\Sigma^{(1)}$.

**THEOREM 2.1.** There exists a positive constant $K > 0$ depending only on $\lambda, p$ and $r_+$ such that the following inequality holds: for all $N \ge 2$,

$$ \mathbb{E} \left[ \left\| \Sigma^{(1)} - (1 - \lambda)m(1 - m) \left( A^N + \frac{(1 - \lambda)^2 p^3 (r_+ - r_-)}{1 - p^2 (r_+ - r_-)^2} \frac{J_N}{N} \right) \right\|_{\max} \right] \le K \sqrt{\frac{\log(N)}{N^3}}, $$

where the matrices $\Sigma^{(1)}$ and $A^N$ are defined in (2) and (3) respectively, $J_N$ is the matrix full of ones and the real value $m$ is defined in (9).
The proof of Theorem 2.1 is given in Section 5. For each $j \in [N]$, let us denote

$$ (10) \qquad \sigma_j^N = \sum_{i=1}^{N} \Sigma_{ij}^{(1)}. $$

Write $\sigma^N = (\sigma_1^N, \dots, \sigma_N^N)$. The vector $\sigma^N = (1_N)^\top \Sigma^{(1)}$ is obtained by summing over the columns of the matrix $\Sigma^{(1)}$. Our second main result reads as follows.

**THEOREM 2.2.** There exists a constant $K > 0$ depending only on $\lambda, p$ and $r_+$ such that we have for all $N \ge 2$,

$$ \mathbb{E} \left[ \left\| \sigma^N - (1-\lambda)pm(1-m) \left( \frac{(1-\lambda)^2 p^2 (r_+ - r_-)}{1 - p^2 (r_+ - r_-)^2} + 1_{P_+} - 1_{P_-} \right) \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}, $$
where the vector $\sigma^N$ is defined in (10) and the real value $m$ is defined in (9).

**REMARK 2.3.** The constants $K$ appearing in Theorems 2.1 and 2.2 diverge when $\lambda \to 0$ (the ergodicity of dynamics breaks) or $\lambda \to 1$ (the dependence between the components of the model becomes very weak). In such regimes, the upper bounds provided by Theorems 2.1 and 2.2 deteriorate. For fixed values of $\lambda$ (i.e. as functions of $p$ and $r_+$ only), these constants $K$ remain bounded.

The proof of Theorem 2.2 is also given in Section 5. Denoting $F_{\varepsilon}^N = \{\|\sigma^N - (1-\lambda)pm(1-m)\left(\frac{(1-\lambda)^2p^2(r_+-r_-)}{1-p^2(r_+-r_-)^2} + 1_{P_+} - 1_{P_-}\right)\|_{\infty} \le \varepsilon\}$ for $\varepsilon > 0$, it follows from Markov inequality and Theorem 2.2 that $\mathbb{P}(F_{\varepsilon}^N) \ge 1 - K\varepsilon^{-1}\sqrt{\log(N)/N}$. Moreover, writing $i \sim j$ to indicate that the components $i$ and $j$ belong to the same community and $i \not\sim j$ otherwise, one can easily check that on $F_{\varepsilon}^N$, the two communities $P_+$ and $P_-$ are well separated via $\sigma^N$ in the sense that

$$ i \sim j \Rightarrow |\sigma_i^N - \sigma_j^N| \le 2\varepsilon \quad \text{and} \quad i \not\sim j \Rightarrow |\sigma_i^N - \sigma_j^N| \ge 2(1-\lambda)pm(1-m) - 2\varepsilon, $$

(see Figure 1 for a graphical representation). This implies that taking $\varepsilon$ small enough, one could use the vector $\sigma^N$ to exactly recover the sets $P_+$ and $P_-$, with high probability.

In practice, however, the vector $\sigma^N$ cannot be used to that end because the matrix $\Sigma^{(1)}$ is not known. To circumvent this issue, one can try to find these sets by using the empirical estimate of $\sigma^N$.

Let $X_1, \dots, X_T$ be a sample of a stationary Markov chain with transition probabilities given by (7) and (1), associated with some realization $\theta$ of the random environment which is not observed. Given this sample, one can then estimate the vector $\sigma^N$ by the vector $\hat{\sigma}^N = (\hat{\sigma}_1^N, \dots, \hat{\sigma}_N^N)$ defined as,

$$ (11) \qquad \hat{\sigma}_j^N = \sum_{i=1}^{N} \left[ \frac{1}{T-1} \sum_{t=2}^{T} X_{i,t} X_{j,t-1} - \frac{Z_{i,T} Z_{j,T}}{T} \right], \quad j \in [N], $$

where $Z_{i,T} = \sum_{t=1}^{T} X_{i,t}$ denotes the number of signals emitted by the $i$-th chain over the discrete time interval $\{1, \dots, T\}$.

**THEOREM 2.4.** There exists a constant $K > 0$ depending only on $\lambda$ such that for all $N \ge 2$ and $T \ge K \log(N)$,

$$ \mathbb{E} [\|\hat{\sigma}^N - \sigma^N\|_{\infty}] \le K \frac{N \log(TN)}{\sqrt{T}}. $$
<figure><img src="image_3.png" /><figcaption>FIG 1. Representation of the concentration of the values of $\sigma_j^N$ and $\hat{\sigma}_j^N$ around $\sigma^+$ for all $j \in \mathcal{P}_+$ and $\sigma^-$ for all $j \in \mathcal{P}_-$ under the event $G_{\varepsilon}^N \cap F_{\varepsilon}^N$. The exact expressions of $\sigma^+$ and $\sigma^-$ appear in Theorem 2.2. The distance between each ticks on the axis is equal to $\varepsilon$. The upper row correspond to the values of $\hat{\sigma}^N$ and the lower row correspond to the values of $\sigma^N$. The dashed lines connect $\sigma_j^N$ and $\hat{\sigma}_j^N$ for each component $j \in [N]$.</figcaption></figure>

As the proof of Theorems 2.1 and 2.2, the proof of Theorem 2.4 is provided in Section 5.

For each $\varepsilon > 0$, let us denote $G_{\varepsilon}^N = \{\|\hat{\sigma}^N - \sigma^N\|_{\infty} \le \varepsilon\}$. Combining Markov inequality with Theorems 2.2 and 2.4, it then follows that

$$IP(G_{\varepsilon}^N \cap F_{\varepsilon}^N) \ge 1 - K\varepsilon^{-1} \left(\sqrt{\log(N)/N} + N\log(TN)/\sqrt{T}\right).$$

Moreover, on $G_{\varepsilon}^N \cap F_{\varepsilon}^N$, we have that

$$i \sim j \Rightarrow |\hat{\sigma}_i^N - \hat{\sigma}_j^N| \le 4\varepsilon \quad \text{and} \quad i \not\sim j \Rightarrow |\hat{\sigma}_i^N - \hat{\sigma}_j^N| \ge 2(1-\lambda)pm(1-m) - 4\varepsilon,$$

(see Figure 1 for a graphical representation). With $\varepsilon$ small enough, exact recovery is feasible with the standard k-means algorithm for instance. Let us denote $\hat{P}_{-}^N$ and $\hat{P}_{+}^N$ the two subsets which are the output of k-means algorithm applied to $\hat{\sigma}^N$. Under the assumption that $\varepsilon < (1 - \lambda)m(1 - m)p/2$, it is easy to check that the algorithm (with $k = 2$) initialized with the means $m_1 = \min_j \hat{\sigma}_j^N$ and $m_2 = \max_j \hat{\sigma}_j^N$ returns the partition made of $S_1 = \mathcal{P}_-$ and $S_2 = \mathcal{P}_+$ at its first step. Furthermore, if $\varepsilon < (1 - \lambda)m(1 - m)p/4$, then this partition is stable and so $\hat{P}_{-}^N = \mathcal{P}_-$ and $\hat{P}_{+}^N = \mathcal{P}_+$.

Hence, in the asymptotic framework where $T = T_N \to \infty$ and $\varepsilon = \varepsilon_N \to 0$ as $N \to \infty$, we have the following corollary.

**COROLLARY 2.5.** Let $\hat{\sigma}^N$ be the vector defined in (11). If as $N \to \infty$,

$$\varepsilon_N^{-1} \left( \sqrt{\log(N)N^{-1}} + N \log(NT_N) / \sqrt{T_N} \right) \to 0,$$

then the assumption of Theorem 2.4 is eventually valid and in turn the probability of exact recovery goes to 1:

$$\lim_{N \to \infty} IP \left( \left\{ \hat{P}_{+}^N = \mathcal{P}_{+} \right\} \cap \left\{ \hat{P}_{-}^N = \mathcal{P}_{-} \right\} \right) = 1,$$

where $\{\hat{P}_{-}^N, \hat{P}_{+}^N\}$ denotes the partition computed via the k-means algorithm applied to $\hat{\sigma}^N$.

**3. Some matrix notation and previous results.** The aim of this section is to introduce some notation and results from (Chevallier, Löcherbach and Ost, 2024). Let $L^{N, \bullet^-}$ be the vector of row-wise sums over the columns corresponding to $\mathcal{P}_-$ of the matrix $A^N$, i.e. $L^{N, \bullet^-} = A^N 1_{\mathcal{P}_-}$. Then, for all $x \in \{0, 1\}^N$, the vector $p_{\theta}(x) = (p_{\theta,1}(x), \dots, p_{\theta,N}(x))$ can be written as,

$$(12) \qquad p_{\theta}(x) = \mu 1_N + (1 - \lambda) (A^N x - L^{N, \bullet^-}).$$
Now, let us denote $m^N$ the (conditional) mean vector given by

$$ (13) \qquad m^N = (m_1^N, \dots, m_N^N) = \mathbb{E}_{\theta}(X_0) = \mathbb{E}_{\theta}[p_{\theta}(X_0)], $$

where we used the stationarity of the process, and denote its spatial mean $\overline{m^N} = N^{-1} \sum_{i=1}^{N} m_i^N$.

With this notation, let us remark that $\Sigma^{(1)}$ rewrites as

$$ (14) \qquad \Sigma^{(1)} = \mathbb{E}_{\theta}[(X_1 - m^N)(X_0 - m^N)^{\top}]. $$

Taking expectation in both sides of Equation (12), one has

$$ (15) \qquad m^N = \mu 1_N + (1 - \lambda) (A^N m^N - L^{N, \bullet -}), $$

which in particular gives that the centered version of $p_{\theta}(X_0)$ is

$$ (16) \qquad p_{\theta}(X_0) - m^N = (1 - \lambda) A^N (X_0 - m^N). $$

Moreover, using the matrix $Q^N = (I_N - (1 - \lambda)A^N)^{-1}$ and the notation for its row-wise sums $\ell^N = Q^N 1_N$ and $\ell^{N, \bullet -} = Q^N 1_{P_-}$, one has the following expression of $m^N$, see (Chevallier, Löcherbach and Ost, 2024, Equation (13)),

$$ (17) \qquad m^N = \mu \ell^N + 1_{P_-} - \ell^{N, \bullet -}. $$

Finally, let us mention that all coordinates of the vector of temporal means $m^N$ are expected to converge to the asymptotic spatio-temporal mean $m$ defined in (9). This statement is formalized in Equation (45).

**4. Asymptotic approximations.** Let us denote $v^N$ the variance vector given by

$$ (18) \qquad v^N = (m_1^N(1 - m_1^N), \dots, m_N^N(1 - m_N^N)). $$

According to the end of previous section, $v^N$ is expected to converge to $m(1 - m)1_N$. This is formalized in Equation (22) below.

Let us denote $\Sigma^{(0)} = \mathbb{E}_{\theta}[(X_0 - m^N)(X_0 - m^N)^{\top}]$ the simultaneous covariance matrix whose $(i, j)$-entry is

$$ (19) \qquad \Sigma_{ij}^{(0)} = \mathrm{Cov}_{\theta}(X_{i,0}, X_{j,0}), \quad i, j \in [N]. $$

Remark that $\Sigma^{(0)}$ depends on $N$ but we chose to not highlight this dependence in the notation even if we are interested in its asymptotic when $N \to \infty$. By the way, the main objective of this section is to provide an asymptotic approximation of the matrix $\Sigma^{(0)}$ as $N \to \infty$ (see Proposition 4.2).

The reason to study the simultaneous covariance matrix is that it serves as a proxy to get the asymptotics of the one unit time-lagged covariance matrix $\Sigma^{(1)}$. Indeed, we have the following relation:

**PROPOSITION 4.1.** $\Sigma^{(1)} = (1 - \lambda) A^N \Sigma^{(0)}$.

**PROOF.** By conditioning with respect to $\mathcal{F}_0$ in Equation (14) and then using (16), one has

$$ \Sigma^{(1)} = \mathbb{E}_{\theta}[(p_{\theta}(X_0) - m^N)(X_0 - m^N)^{\top}] = \mathbb{E}_{\theta}[(1 - \lambda)A^N(X_0 - m^N)(X_0 - m^N)^{\top}], $$

which ends the proof. $\square$

For a $N$-dimensional vector $u = (u_1, \dots, u_N)$, we write $\mathrm{diag}(u)$ to denote the $N$-by-$N$ diagonal matrix having value $u_i \in \mathbb{R}$ associated with its $i$-th diagonal entry, i.e. $\mathrm{diag}(u)_{ii} = u_i$. Let us denote $a^N = 1_N(1_N)^{\top} - I_N$ the $N$-by-$N$ matrix full of ones except on the diagonal.
**PROPOSITION 4.2.** There exists a constant $K > 0$ depending only on $\lambda, p$ and $r_+$ such that for all $N \ge 1$,

$$ \mathbb{E} \left[ \left\| \Sigma^{(0)} - \left( \operatorname{diag}(v^N) + \frac{(1-\lambda)^2 p^2 m (1-m)}{1-p^2(r_+-r_-)^2} \frac{a^N}{N} \right) \right\|_{\max} \right] \le K \sqrt{\frac{\log(N)}{N^3}}. $$

The proof of Proposition 4.2 is interesting per se and is postponed to Section 6. It relies on the fact that $\Sigma^{(0)}$ satisfies some kind of “Stein matrix equation” (see Proposition 6.1). To use this structure, we introduce the *vectorization operator*, the *Kronecker product* and their useful properties in Section 6.

Finally, the last main ingredients to describe the asymptotics of $\Sigma^{(1)}$ are the asymptotics of the vector of row-wise and column-wise sums of $A^N$, denoted by $L^N = A^N 1_N$ and $C^N = (1_N)^\top A^N$ respectively, and the asymptotics of $v^N$.

**PROPOSITION 4.3.** There exists a constant $K > 0$ depending only on $\lambda$ such that for all $N \ge 1$,

$$ (20) \qquad \mathbb{E} \left[ \| L^N - p(r_+ - r_-) 1_N \|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}, $$

$$ (21) \qquad \mathbb{E} \left[ \| C^N - p(1_{P_+} - 1_{P_-}) \|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}, $$

$$ (22) \qquad \mathbb{E} \left[ \| v^N - m(1-m) 1_N \|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}. $$

Proposition 4.3 is a particular case of Lemma 8.2 whose proof is given in Section 8.

**5. Proof of main results.** In this section, we prove our main results, namely, Theorems 2.1, 2.2 and 2.4. We start with the proof of Theorem 2.1.

**PROOF OF THEOREM 2.1.** First of all, remind that $\Sigma^{(1)} = (1-\lambda)A\Sigma^{(0)}$ (see Proposition 4.1). By Proposition 4.2 and using (8) with the fact that $\|A^N\|_{\infty} \le 1$, we have,

$$ (23) \qquad \mathbb{E} \left[ \left\| \Sigma^{(1)} - (1-\lambda) A^N \left( \operatorname{diag}(v^N) + \frac{(1-\lambda)^2 p^2 m (1-m)}{1-p^2(r_+-r_-)^2} \frac{a^N}{N} \right) \right\|_{\max} \right] \le K \sqrt{\frac{\log(N)}{N^3}}. $$

On the one hand, Equation (22) implies that

$$ \mathbb{E} \left[ \| \operatorname{diag}(v^N) - m(1-m) I_N \|_1 \right] \le K \sqrt{\frac{\log(N)}{N}}, $$

so that using (8) with the fact that $\|A^N\|_{\max} \le 1/N$, we get

$$ (24) \qquad \mathbb{E} \left[ \| A^N \operatorname{diag}(v^N) - m(1-m) A^N \|_{\max} \right] \le K \sqrt{\frac{\log(N)}{N^3}}. $$

On the other hand, for all $i, j \in [N]$, $(A^N a^N)_{ij} = \sum_{k=1:k \ne j}^n A_{ik}^N = L_i^N - N^{-1}\theta_{ij}$. Hence, Equation (20) implies that

$$ (25) \qquad \mathbb{E} \left[ \| A^N a^N - p(r_+ - r_-) J_N \|_{\max} \right] \le K \sqrt{\frac{\log(N)}{N}}. $$

Finally, it suffices to combine Equations (23)-(25) to conclude the proof. $\square$
We now prove Theorem 2.2.

PROOF OF THEOREM 2.2. By Theorem 2.1 and the definition of the vector $\sigma^N = (1_N)^\top \Sigma^{(1)}$, one can deduce that

$$ \mathbb{E} \left[ \left\| \sigma^N - (1-\lambda)m(1-m) \left( C^N + \frac{(1-\lambda)^2 p^3 (r_+ - r_-)}{1-p^2(r_+ - r_-)^2} 1_N \right) \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}, $$

and the desired result follows from (21). □

The proof of Theorem 2.4 is given below.

PROOF OF THEOREM 2.4. First of all, observe that

$$ \begin{aligned} \sigma_j^N &= \sum_{i=1}^N \Sigma_{ij}^{(1)} = \sum_{i=1}^N \mathrm{Cov}_{\theta}(X_{i,1}, X_{j,0}) \\ &= N \mathrm{Cov}_{\theta}(\overline{X}_1, X_{j,0}) = N \left( \mathbb{E}_{\theta}(\overline{X}_1 X_{j,0}) - \overline{m}^N m_j^N \right), \end{aligned} $$

and

$$ \hat{\sigma}_j^N = \frac{N}{T-1} \sum_{t=2}^{T} \overline{X}_t X_{j,t-1} - N \frac{\overline{Z}_T}{T} \frac{Z_{j,T}}{T}. $$

Hence,

$$ \|\sigma^N - \hat{\sigma}_j^N\|_{\infty} = N \max_{j \in [N]} \left| \frac{1}{T-1} \sum_{t=2}^{T} \overline{X}_t X_{j,t-1} - \mathbb{E}_{\theta}(\overline{X}_1 X_{j,0}) + \overline{m}^N m_j^N - \frac{\overline{Z}_T}{T} \frac{Z_{j,T}}{T} \right|. $$

By the triangle inequality and the facts that $\overline{m}^N \le 1$ and $T^{-1}Z_{j,T} \le 1$, we have that

$$ \max_{j \in [N]} \left| \overline{m}^N m_j^N - \frac{\overline{Z}_T}{T} \frac{Z_{j,T}}{T} \right| \le \max_{j \in [N]} \left| m_j^N - \frac{Z_{j,T}}{T} \right| + \left| \overline{m}^N - \frac{\overline{Z}_T}{T} \right|, $$

which together with the previous inequality gives

$$ (26) \quad \begin{aligned} \|\sigma^N - \hat{\sigma}_j^N\|_{\infty} \le N \max_{j \in [N]} & \left| \frac{1}{T-1} \sum_{t=2}^{T} \overline{X}_t X_{j,t-1} - \mathbb{E}_{\theta}(\overline{X}_1 X_{j,0}) \right| \\ & + N \max_{j \in [N]} \left| m_j^N - \frac{Z_{j,T}}{T} \right| + N \left| \overline{m}^N - \frac{\overline{Z}_T}{T} \right|. \end{aligned} $$

Now, we use Theorem 3 of (Ost and Reynaud-Bouret, 2020) to control each terms in the right hand side above.

First, for each $j \in [N]$, let $f_j : \{0,1\}^N \to \{0,1\}$ be the projection on the j-th coordinate: $f_j(x) = x_j$ and denote

$$ V(f_j) = \frac{1}{T} \sum_{t=1}^{T} (f_j(X_t) - \mathbb{E}_{\theta}(f_j(X_t))) = \frac{Z_{j,T}}{T} - m_j^N. $$

By Theorem 3 of (Ost and Reynaud-Bouret, 2020) with $m = M = 1$, $\mathcal{F} = \{f_j : j \in [N]\}$, $F = [N]$ and $\log(1/(1-\lambda))/2$ as their $\theta$, we have that there exists constants $K_1, K_2 > 0$ depending only on $\lambda$ such that

$$ \mathbb{P} \left( \cup_{j \in [N]} \left\{ V(f_j) > \sqrt{K_1 \frac{(1 + \log(TN))\xi}{T}} \right\} \right) \le K_2 T^{-1} + 2Ne^{-\xi}, $$
for all $\xi > 0$ and $T$ such that

$$T \ge 1 + \frac{2}{\log((1 - \lambda)^{-1})}(2\log(T) + \log(N)).$$

Since $|V(f_j)| \le 1$ for all $j \in [N]$, the probability inequality above implies that for any $\xi > 0$,

$$(27) \quad \mathrm{IE} \left[ \max_{j \in [N]} \left| \frac{Z_{j,T}}{T} - m_j^N \right| \right] \le \sqrt{K_1 \frac{(1 + \log(TN))\xi}{T}} + K_2 T^{-1} + 2Ne^{-\xi}.$$

Since, $|\overline{m^N} - \overline{\frac{Z_T}{T}}| \le \max_{j \in [N]} |Z_{j,T}/T - m_j^N|$, Inequality (27) implies that

$$(28) \quad \mathrm{IE} \left[ \left| \overline{m^N} - \overline{\frac{Z_T}{T}} \right| \right] \le \sqrt{K_1 \frac{(1 + \log(TN))\xi}{T}} + K_2 T^{-1} + 2Ne^{-\xi}.$$

Second, let us define $g_j : \{0,1\}^N \times \{0,1\}^N \to \{0,1\}$ for each $j \in [N]$ as $g_j(x, y) = x_j \bar{y}$ and denote

$$\begin{aligned} W(g_j) &= \frac{1}{T-1} \sum_{t=2}^{T} (g_j(X_{t-1}, X_t) - \mathrm{IE}_{\theta}(g_j(X_{-1}, X_t))) \\ &= \frac{1}{T-1} \sum_{t=2}^{T} (\overline{X}_t X_{j,t-1} - \mathrm{IE}_{\theta}(\overline{X}_1 X_{j,0})). \end{aligned}$$

By Theorem 3 of (Ost and Reynaud-Bouret, 2020) with $M = 1$, $m = 2$, $\mathcal{F} = \{g_j : j \in [N]\}$, $F = [N]$ and $\log(1/(1 - \lambda))/2$ as their $\theta$, there exist constants $K_3, K_4 > 0$ depending only on $\lambda$ such that

$$\mathrm{IE} \left( \bigcup_{j \in [N]} \left\{ W(g_j) > \sqrt{K_3 \frac{(1 + \log(TN))\xi}{T}} \right\} \right) \le K_4 T^{-1} + 2Ne^{-\xi},$$

for all $\xi > 0$ and $T$ such that

$$T \ge 2 + \frac{2}{\log((1 - \lambda)^{-1})}(2\log(T) + \log(N)).$$

Using $|W(g_j)| \le 1$ and arguing similarly as the previous case, one can show that for any $\xi > 0$,

$$(29) \quad \mathrm{IE} \left[ \max_{j \in [N]} \left| \frac{1}{T-1} \sum_{t=2}^{T} (\bar{X}_t X_{j,t-1} - \mathrm{IE}_{\theta}(\bar{X}_1 X_{j,0})) \right| \right] \le \sqrt{K_3 \frac{(1 + \log(TN))\xi}{T}} + K_4 T^{-1} + 2Ne^{-\xi}.$$

Putting together Equations (26)-(29) with $\xi = \log(TN)$, we show that the conclusion of the theorem is true for all $N \ge 2$ and $T \ge K \log(TN)$ for some $K$ depending only on $\lambda$. Hence, to conclude the proof, it remains to show that the conclusion of the theorem is also true for all $N \ge 2$ and $T \ge K' \log(N)$ for some $K'$ depending only on $\lambda$. To see that, take $T_0 = T_0(K)$ such that $T/2 \ge K \log(T)$ for all $T \ge T_0$, define $K' = (T_0 2K)/\log(2)$, and observe that if $T \ge K' \log(N) > T_0$ then $T = T/2 + T/2 \ge K \log(T) + (K'/2) \log(N) > K \log(T) + (T_0 K/\log(2)) \log(N) > K \log(T) + K \log(N)$. $\square$
**6. Asymptotic approximation of the simultaneous covariance matrix.** This section is devoted to the proof of Proposition 4.2. As mentioned above, the proof rely on some “Stein matrix equation” structure stated in Proposition 6.1 below.

For a $N$-by-$N$ matrix $M$, we denote $d_0(M)$ the same matrix where the diagonal terms are replaced by zeros, i.e. $d_0(M)_{ij} = M_{ij}1_{i \neq j}$.

**PROPOSITION 6.1.** The simultaneous covariance matrix $\Sigma^{(0)}$ is equal to $v^N$ on the diagonal and satisfies $\Sigma^{(0)} = (1-\lambda)^2 A^N \Sigma^{(0)} (A^N)^\top$ out of the diagonal, i.e.

$$ \Sigma^{(0)} = (1 - \lambda)^2 d_0 (A^N \Sigma^{(0)} (A^N)^\top) + \mathrm{diag}(v^N). $$

**PROOF.** Let $i \in [N]$. Under the stationary distribution, $X_{i,0}$ is a Bernoulli variable with parameter $m_i^N$. Hence, $\Sigma_{ii}^{(0)} = \mathrm{Var}_\theta(X_{i,0}) = m_i^N(1 - m_i^N) = v_i^N$.

Let $i, j \in [N]$ with $i \neq j$. By conditional independence between $X_{i,0}$ and $X_{j,0}$ given $X_{-1}$ and stationarity of the process, we have that

$$ \Sigma_{ij}^{(0)} = \mathrm{Cov}_\theta(X_{i,0}, X_{j,0}) = \mathrm{Cov}_\theta(p_{\theta,i}(X_{-1}), p_{\theta,j}(X_{-1})) = \mathrm{Cov}_\theta(p_{\theta,i}(X_0), p_{\theta,j}(X_0)). $$

Yet, the centered version of $p_{\theta,i}(X_0)$ is $p_{\theta,i}(X_0) - m_i^N$ so that using Equation (16) it is easy to get

$$ \Sigma_{ij}^{(0)} = (1 - \lambda)^2 \mathbb{E} [A^N (X_0 - m^N)(X_0 - m^N)^\top (A^N)^\top]_{ij} = (1 - \lambda)^2 [A^N \Sigma^{(0)} (A^N)^\top]_{ij}. $$

<figure><img src="image_4.png" /></figure>

**6.1. Vectorization and Kronecker product.** Two mathematical objects are known to be very useful to deal with “Stein matrix equations”: vectorization and Kronecker product. Here are their definition and useful properties.

The vectorization $\mathrm{vec}(M)$ of a $N$-by-$N$ matrix $M$ is the $N^2$-dimensional vector constructed column by column. More precisely, $M_{ij}$ associated with its $k$-th coordinate, where $k = i + (j-1)N$ with $i, j \in [N]$. With slight abuse of notation, will write $\mathrm{vec}(M)_{ij}$ to refer to the coordinate $k = i + (j-1)N$ of the vector $\mathrm{vec}(M)$. Recall that the vec operator is an invertible linear operator on the set of all matrices with real value entries. Given two $N$-by-$N$ matrices $A$ and $B$, we denote $A \otimes B$ the Kronecker product between $A$ and $B$. It has the following block representation:

$$ A \otimes B = \begin{bmatrix} A_{11}B & \cdots & A_{1n}B \\ \vdots & \ddots & \vdots \\ A_{m1}B & \cdots & A_{mn}B \end{bmatrix}. $$

More precisely, this is a $N^2$-by-$N^2$ matrix having value $A_{ik}B_{j\ell}$ associated with entry $(u, v)$, where $u = N(i-1) + j$ and $v = N(k-1) + \ell$ with $i, j, k, \ell \in [N]$. With slight abuse of notation, will write $(A \otimes B)_{ij,k\ell}$ to refer to the entry $(u, v) = (N(i-1) + j, N(k-1) + \ell)$ of the matrix $A \otimes B$.

Their most useful property regarding matrix equations is that, for all $N$-by-$N$ matrices $A, B, M$,

$$ (30) \qquad \mathrm{vec}(AMB) = (B^\top \otimes A) \mathrm{vec}(M). $$

Furthermore, remark that $\|M\|_{\max} = \|\mathrm{vec}(M)\|_{\infty}$.
PROOF OF PROPOSITION 4.2. The main idea is to apply the vectorization operator to the matrix equation of Proposition 6.1. In order to emphasize the difference between vectors of size N (for instance $m^N$, $v^N$) or size $N^2$ (for instance $\text{vec}(\Sigma^{(0)})$), we use bold notation for the latter. The same convention is used for $N^2$-by-$N^2$ matrices. Let us then introduce the notation:

$$ (31) \qquad \Sigma^{(0)} := \text{vec}(\Sigma^{(0)}), \quad v^N := \text{vec}(\text{diag}(v^N)), \quad a^N := \text{vec}(a^N) $$

and $\mathbf{I}$ denotes the $N^2$-by-$N^2$ identity matrix. Obviously, the abuse of notation described above for the coordinates of the vectorization and Kronecker product is inherited to the bold notation introduced here.

According to Equation (30), the starting idea to tackle the matrix equation satisfied by $\Sigma^{(0)}$ is to consider its vectorization and the Kronecker product $(1 - \lambda)^2(A^N \otimes A^N)$. Here is a slight modification in order to tackle the $d_0$ operator. Let $\mathbf{L} = (\mathbf{L}_{ij,k\ell})_{i,j,k,\ell \in [N]}$ be the $N^2$-by-$N^2$ matrix defined as

$$ \mathbf{L}_{ij,k\ell} = \begin{cases} (1 - \lambda)^2 (A^N \otimes A^N)_{ij,k\ell}, & \text{if } i \neq j, k, \ell \in [N] \\ 0, & \text{otherwise.} \end{cases} $$

In particular, remark that

$$ (32) \qquad \forall \mathbf{x} \in \mathbb{R}^{N^2}, i \in [N], \quad [\mathbf{L}\mathbf{x}]_{ii} = 0. $$

Moreover, $\|\mathbf{L}\|_{\infty} \le (1 - \lambda)^2 < 1$ so that $\mathbf{I} - \mathbf{L}$ is invertible.

**Step 1.** Here we prove that

$$ \Sigma^{(0)} = (\mathbf{I} - \mathbf{L})^{-1} \mathbf{v}^N = \sum_{k=0}^{\infty} \mathbf{L}^k \mathbf{v}^N. $$

*Proof of Step 1.* By Proposition 6.1 and linearity of vec, we have

$$ (33) \qquad \Sigma^{(0)} = (1 - \lambda)^2 \text{vec} \left( d_0 \left( A^N \Sigma^{(0)} (A^N)^\top \right) \right) + \mathbf{v}^N. $$

On the one hand, observe that for any $i, j \in [N]$ with $i \neq j$,

$$ \begin{aligned} (1 - \lambda)^2 \text{vec} \left( d_0 \left( A^N \Sigma^{(0)} (A^N)^\top \right) \right)_{ij} &= (1 - \lambda)^2 [A^N \Sigma^{(0)} (A^N)^\top]_{ij} \\ &= (1 - \lambda)^2 \sum_{k=1}^{N} \sum_{\ell=1}^{N} A_{ik}^N \Sigma_{k\ell}^{(0)} A_{j\ell}^N \\ &= \sum_{k=1}^{N} \sum_{\ell=1}^{N} \mathbf{L}_{ij,k\ell} \Sigma^{(0)}_{k\ell} = [\mathbf{L}\Sigma^{(0)}]_{ij}. \end{aligned} $$

On the other hand, for $i = j$, the definition of $d_0$ implies

$$ (1 - \lambda)^2 \text{vec} \left( d_0 \left( A^N \Sigma^{(0)} (A^N)^\top \right) \right)_{ii} = 0 = [\mathbf{L}\Sigma^{(0)}]_{ii}, $$

where the second equality comes from (32). All in all, we have proved that

$$ \Sigma^{(0)} = \mathbf{L}\Sigma^{(0)} + \mathbf{v}^N, $$

which, in turn, concludes Step 1.

The rest of the proof is devoted to the study of the iterates of the form $\mathbf{L}^k \mathbf{v}^N$. The next two steps provide useful estimates to tackle this study.
**Step 2.** In this step, we prove that there exists a constant $K > 0$ depending only on $\lambda$ such that for all $N \ge 1$,

$$ \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} m_k^N - p^2 m \right| \right] \le K \sqrt{\log(N)/N}. $$

*Proof of Step 2.* First, since $\overline{m^N}$ is the mean of $m^N$, we can write

$$ \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} m_k^N = p^2 m + p^2 (\overline{m^N} - m) \\ + \frac{1}{N} \sum_{k=1}^{N} (\theta_{ik} \theta_{jk} - p^2) (m_k^N - \overline{m^N}) + \overline{m^N} \frac{1}{N} \sum_{k=1}^{N} (\theta_{ik} \theta_{jk} - p^2). $$

Next, use Cauchy-Schwarz inequality and the fact that $|\theta_{ik}\theta_{jk} - p^2| \le 1$ to deduce that

$$ \left| \frac{1}{N} \sum_{k=1}^{N} (\theta_{ik} \theta_{jk} - p^2) (m_k^N - \overline{m^N}) \right| \le N^{-1/2} \|m^N - \overline{m^N}\|_2, $$

so that, using the fact that $|\overline{m^N}| \le 1$,

$$ \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} m_k^N - p^2 m \right| \right] \le p^2 \mathbb{E} \left[ |\overline{m^N} - m| \right] \\ + N^{-1/2} \mathbb{E} \left[ \|m^N - \overline{m^N}\|_2 \right] + \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} - p^2 \right| \right]. $$

The last term on the right-hand side is controlled thanks to the maximal inequality (Lemma 8.1).

Now, using (Chevallier, Löcherbach and Ost, 2024, Proposition 4.7), one can show that

$$ (34) \qquad \mathbb{E} \left[ |\overline{m^N} - m| \right] + N^{-1/2} \mathbb{E} \left[ \|m^N - \overline{m^N}\|_2 \right] \le KN^{-1/2}, $$

for some constant depending only on $\lambda$, which ends the proof of this step.

**Step 3.** Here, we show that there exists a constant $K > 0$ depending only on $\lambda$ such that for all $N \ge 1$,

$$ \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} (m_k^N)^2 - p^2 m^2 \right| \right] \le K \sqrt{\log(N)/N}. $$

*Proof of Step 3.* We proceed similarly as in the proof of Step 2. First of all, check that we can write

$$ \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} (m_k^N)^2 = p^2 m^2 + \frac{p^2}{N} \sum_{k=1}^{N} ((m_k^N)^2 - m^2) \\ + \frac{1}{N} \sum_{k=1}^{N} (\theta_{ik} \theta_{jk} - p^2) ((m_k^N)^2 - (\overline{m^N})^2) + (\overline{m^N})^2 \frac{1}{N} \sum_{k=1}^{N} (\theta_{ik} \theta_{jk} - p^2). $$
First, observe that

$$ \frac{1}{N} \sum_{k=1}^{N} ((m_k^N)^2 - m^2) = \frac{1}{N} \sum_{k=1}^{N} (m_k^N + m)(m_k^N - \overline{m^N}) + (\overline{m^N} - m)(\overline{m^N} + m), $$

so that Cauchy-Schwarz inequality together with the fact that $\|m^N\|_\infty \vee m \le 1$ implies that

$$ (35) \qquad \left| \frac{1}{N} \sum_{k=1}^{N} ((m_k^N)^2 - m^2) \right| \le 2N^{-1/2} \|m^N - \overline{m^N} 1_N\|_2 + 2|\overline{m^N} - m|. $$

Second, using Cauchy-Schwarz inequality one can deduce (using the fact that $m_k^N + \overline{m^N} \le 2$),

$$ (36) \qquad \left| \frac{1}{N} \sum_{k=1}^{N} (\theta_{ik}\theta_{jk} - p^2)((m_k^N)^2 - (\overline{m^N})^2) \right| \le 2N^{-1/2} \|m^N - \overline{m^N} 1_N\|_2. $$

Third, combining the inequality $|\overline{m^N}| \le 1$ with the maximal inequality (Lemma 8.1), we also have that

$$ (37) \qquad \mathrm{IE} \left[ \max_{i,j \in [N]: i \ne j} (\overline{m^N})^2 \left| \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} - p^2 \right| \right] \le \sqrt{\frac{\log(\sqrt{2N})}{N}}. $$

Hence, combining Equations (35)-(37) with (34) the result follows.

**Step 4.** We now prove that there exists a constant $K > 0$ depending only on $\lambda > 0$ such that for all $N \ge 1$,

$$ \mathrm{IE} \left[ \left\| \mathbf{L} \mathbf{v}^N - \frac{(1-\lambda)^2 p^2 m (1-m)}{N} \mathbf{a}^N \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N^3}}. $$

*Proof of Step 4.* Let $i, j \in [N]$. First, observe that $[\mathbf{L} \mathbf{v}^N]_{ii} = 0 = \mathbf{a}_{ii}^N$ by (32) and the definition of $\mathbf{a}^N$. Then, if $i \ne j$,

$$ \begin{aligned} [\mathbf{L} \mathbf{v}^N]_{ij} &= \sum_{k=1}^{N} \mathbf{L}_{ij,kk} m_k^N (1 - m_k^N) \\ &= \frac{(1-\lambda)^2}{N^2} \left( \sum_{k=1}^{N} \theta_{ik} \theta_{jk} m_k^N (1 - m_k^N) \right), \end{aligned} $$

and then use Steps 2 and 3 to obtain that

$$ \mathrm{IE} \left[ \max_{i,j \in [N]: i \ne j} \left| [\mathbf{L} \mathbf{v}^N]_{ij} - \frac{(1-\lambda)^2 p^2 m (1-m)}{N} \right| \right] \le K \sqrt{\log(N)/N^3}. $$

**Step 5.** Here, we show that $\mathbf{a}^N$ is almost an eigenvector of $\mathbf{L}$ in the sense that: there exists a constant $K > 0$ depending only on $\lambda > 0$ such that for all $N \ge 1$,

$$ \mathrm{IE} \left[ \left\| (\mathbf{L} - p^2(r_+ - r_-)^2 \mathbf{I}) \mathbf{a}^N \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}. $$

*Proof of Step 5.* Let $i, j \in [N]$. First, observe that $[(\mathbf{L} - p^2(r_+ - r_-)^2 \mathbf{I})\mathbf{a}^N]_{ii} = 0$ by (32) and the definition of $\mathbf{a}^N$. Second, if $i \ne j$, let us denote
• $k \leftrightarrow \ell$ if and only if $(k, \ell) \in \mathcal{P}_+^2 \cup \mathcal{P}_-^2$ and $k \neq \ell$,

• $k \leftrightarrow \ell$ if and only if $(k, \ell) \notin \mathcal{P}_+^2 \cup \mathcal{P}_-^2$.

Let us remark that the two relations above are not complementary. Indeed, $k = \ell$ does not satisfy $k \leftrightarrow \ell$ nor $k \leftrightarrow \ell$. Then,

$$ [\mathbf{L}\mathbf{a}^N]_{ij} = \sum_{(k,\ell) \in [N]^2; \ell \neq k} \mathbf{L}_{ij,k\ell} = \frac{(1-\lambda)^2}{N^2} \left( \sum_{k \leftrightarrow \ell} \theta_{ik}\theta_{j\ell} - \sum_{k \leftrightarrow \ell} \theta_{ik}\theta_{j\ell} \right). $$

The cardinality of the summation sets $\{k \leftrightarrow \ell\}$ and $\{k \leftrightarrow \ell\}$ above are respectively $N^2[(r_+^N)^2 + (r_-^N)^2 - 1/N]$ and $2N^2r_+^N r_-^N$. Hence, $[\mathbf{L}\mathbf{a}^N]_{ij} = (1-\lambda)^2(b^N + D_1^N(i,j) - D_2^N(i,j))$ where

$$ \left\{ \begin{array}{l} D_1^N(i,j) = N^{-2} \sum_{k \leftrightarrow \ell} (\theta_{ik}\theta_{j\ell} - p^2), \\ D_2^N(i,j) = N^{-2} \sum_{k \leftrightarrow \ell} (\theta_{ik}\theta_{j\ell} - p^2), \\ b^N = [(r_+^N - r_-^N)^2 - \frac{1}{N}]p^2. \end{array} \right. $$

On the one hand, by assumption on $r_+, r_-$, we have

$$ |b^N - p^2(r_+ - r_-)^2| \leq KN^{-1}. $$

On the other hand, the maximal inequality (Lemma 8.1) yields

$$ \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} |D_1^N(i,j)| + |D_2^N(i,j)| \right] \leq 2 \sqrt{\frac{\log(\sqrt{2N})}{N}}. $$

All in all, we get that

$$ \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| [\mathbf{L}\mathbf{a}^N]_{ij} - p^2(r_+ - r_-)^2 \right| \right] \leq K \sqrt{\frac{\log(N)}{N}}, $$

which ends the proof of this step.

We are now in position to end the proof. According to steps 4 and 5, let us write

$$ \mathbf{L}\mathbf{v}^N = \frac{\zeta}{N}\mathbf{a}^N + \varepsilon^{N,1} \quad \text{and} \quad \mathbf{L}\mathbf{a}^N = \gamma\mathbf{a}^N + \varepsilon^{N,2}, $$

where $\zeta = (1-\lambda)^2p^2m(1-m)$, $\gamma = p^2(r_+ - r_-)^2$, $\mathbb{E}[\|\varepsilon^{N,1}\|_{\infty}] \leq K\sqrt{\log(N)/N^3}$ and $\mathbb{E}[\|\varepsilon^{N,2}\|_{\infty}] \leq K\sqrt{\log(N)/N}$. By standard linear algebra, we have $\mathbf{L}^0\mathbf{v}^N = \mathbf{v}^N$, and by induction for all $k \geq 2$,

$$ \mathbf{L}^k \mathbf{v}^N = \frac{\zeta}{N} \left( \gamma^{k-1} \mathbf{a}^N + \sum_{\ell=0}^{k-2} \gamma^{\ell} \mathbf{L}^{(k-2)-\ell} \varepsilon^{N,2} \right) + \mathbf{L}^{k-1} \varepsilon^{N,1}. $$

In particular, denoting $\delta = \max\{\gamma, 1-\lambda\} < 1$ and reminding that $\|\mathbf{L}\|_{\infty} \leq (1-\lambda)^2$, one has for all $k \geq 1$,

$$ (38) \qquad \left\| \mathbf{L}^k \mathbf{v}^N - \frac{\zeta}{N} \gamma^{k-1} \mathbf{a}^N \right\|_{\infty} \leq \zeta (k-1) \delta^k \frac{\| \varepsilon^{N,2} \|_{\infty}}{N} + (1-\lambda)^{k-1} \| \varepsilon^{N,1} \|_{\infty}. $$

To conclude the proof, it suffices to remind the result of step 1, remark that $\sum_{k=1}^{\infty} \gamma^{k-1} = 1/(1-p^2(r_+ - r_-)^2)$ and use the summability of the sequence $(k\delta^k)_{k \geq 0}$. $\square$
**7. Simulation study.** The numeric experiments were made using Julia programming language and are available in the package MeanFieldGraph.jl, as well as the material used to produce the figures below.

The simulation framework follows the one of Chevallier, Löcherbach and Ost (2024). In particular, we consider the parameter $\beta = \mu/\lambda$ instead of $\mu$ because the set of admissible values of $\beta$ is independent of $\lambda$ (which is not the case for $\mu$). If not specified otherwise, the following parameter values are used:

$$N = 50, \quad r_+ = .5, \quad \beta = .5, \quad \lambda = .5, \quad p = .5, \quad N_{\text{simu}} = 1000.$$

Furthermore, the fractions of excitatory and inhibitory components are chosen as $r_+^N = [r_+N]$ and $r_-^N = N - |P_+|$.

<figure><img src="image_5.png" /><figcaption>FIG 2. Estimated PER as a function of N and T (computed over $N_{\text{simu}}$ simulations). The green curve corresponds to the equation $T = N^2$ which is the critical condition of our asymptotic regime.</figcaption></figure>

Figure 2 illustrates the asymptotic regime of Corollary 2.5. It shows that there is a separation between couples $(T, N)$ for which exact recovery never occurs (black color) and those for which exact recovery always occurs (white color). Furthermore, the green curve demarcating this separation is compatible (modulo a log factor) with our condition $(N/T^{1/2}) \log(NT) \to 0$ implying exactly recovery.

Figure 3 gives an overview of the performance of the proposed classification method as one of the parameters $(N, r_+, \beta, \lambda$ or $p)$ vary. The performance is evaluated via two metrics: the probability of exact recovery (PER = $\mathbb{P}(\hat{P}_+^N = P_+ \cap \{\hat{P}_-^N = P_-\})$) and the mean misclassification proportion (MMP = $\mathbb{E}[\text{Card}(\hat{P}_+^N \Delta P_+)/N]$), both estimated via Monte Carlo method ($N_{\text{sim}}$ indicates the number of simulations).

Our objective of exact recovery is hard by essence especially when $N$ is large. Our method achieves its objective only for very large values of $T$. However, it performs well in terms of MMP. Moreover, the transition from almost never exact recovery to almost sure recovery occurs when the MMP is around 0.01 in most of the parameter sets.

In the whole paragraph below, the qualitative descriptions are written with respect to the following framework: the time horizon $T$ of the observation is kept fixed while one of the other parameters vary. The parameters $\beta$ and $r_+$ seem to have no influence on the performance: the influence of $r_+$ on the MMP is rather an artifact of its influence on the proportion of each communities. As $\lambda$ increases, the performance deteriorates as expected because the
<figure><img src="image_6.png" /><figcaption>FIG 3. Estimated PER and MMP as function of T. The increasing dashed lines correspond to the PER (the ribbon corresponds to the standard confidence interval). The decreasing solid lines correspond to the MMP (the width of the ribbon is equal to the standard deviation). The panels correspond to different choices of varying parameter (the non-varying parameters are chosen according to the default values given at the beginning of Section 7). The values of the varying parameter are given by the color legends.</figcaption></figure>

strength of interaction in the system decreases so that the sample contains less and less in-
formation on the underlying graph. Mathematically, it appears in the fact that the separation
between $\sigma^+$ and $\sigma^-$ reduces (see Figure 1) and the constants of Theorems 2.2 and 2.4 di-
verge. Despite the fact that these constants also diverge when $\lambda \to 0$ (see Remark 2.3), our
method performs well for small values of $\lambda$. As $p$ increases, the performance improves as ex-
pected because the strength of interaction in the system increases. Mathematically, it appears
in the fact that the separation between $\sigma^+$ and $\sigma^-$ increases (see Figure 1).

**8. Lemmas based on the maximal inequality.** The two lemmas stated here rely on the maximal inequality (Boucheron, Lugosi and Massart, 2013, Theorem 2.5.) which involves sub-Gaussian random variables. In this section, the notation $\text{SubG}(\nu)$ is used to denote the property “sub-Gaussian with variance factor $\nu$”.
LEMMA 8.1. Let $(\theta_{ij})_{1\le i,j\le N}$ be i.i.d. entries distributed as $\text{Ber}(p)$. Then,

$$ \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| \frac{1}{N} \sum_{k=1}^{N} \theta_{ik} \theta_{jk} - p^2 \right| \right] \le \sqrt{\frac{2 \log(2N)}{N}}, $$

and, for all $A \times B \subset [N]^2$,

$$ \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| \frac{1}{N^2} \sum_{(k,l) \in A \times B} (\theta_{ik} \theta_{jl} - p^2) \right| \right] \le 2 \sqrt{\frac{2 \log(2N)}{N}}, $$

PROOF. The proofs of the two statements are similar. We only prove the second one which is a bit more complex.

Let $i \neq j$. First of all, observe that we can write

$$ \frac{1}{N^2} \sum_{(k,l) \in A \times B} (\theta_{ik} \theta_{jl} - p^2) = \left( \frac{1}{N} \sum_{k \in A} \theta_{ik} \right) \sum_{l \in B} \frac{(\theta_{jl} - p)}{N} + \frac{p|B|}{N} \frac{1}{N} \sum_{k \in A} \frac{(\theta_{ik} - p)}{N}. $$

Next, using that $(|A| \lor |B|)/N \le 1, 0 \le p \le 1$ and $|\theta_{mn}| \le 1$, we obtain from the previous inequality the following upper bound:

$$ (39) \quad \mathbb{E} \left[ \max_{i,j \in [N]: i \neq j} \left| \frac{1}{N^2} \sum_{(k,l) \in A \times B} (\theta_{ik} \theta_{jl} - p^2) \right| \right] \le \mathbb{E} \left[ \max_{j \in [N]} \left| \sum_{l \in B} \frac{(\theta_{jl} - p)}{N} \right| \right] + \mathbb{E} \left[ \max_{i \in [N]} \left| \sum_{k \in A} \frac{(\theta_{ik} - p)}{N} \right| \right]. $$

Now, for all $\ell \in [N]$, the random variables $(\theta_{j\ell} - p)/N$ are bounded by $1/N$ so that they are SubG$(1/N^2)$. Since for any fixed $j$, the random variables $((\theta_{j\ell} - p)/N)_{\ell \in B}$ are independent, their sum is SubG$(|B|/N^2)$. Therefore, (Boucheron, Lugosi and Massart, 2013, Theorem 2.5.) implies

$$ \mathbb{E} \left[ \max_{j \in [N]: \ell \in B} \left| \frac{1}{N} \sum_{\ell \in B} (\theta_{j\ell} - p) \right| \right] \le \sqrt{2 \frac{|B|}{N^2} \log(2N)} \le \sqrt{\frac{2}{N} \log(2N)}. $$

Arguing similarly, one can also show that

$$ \mathbb{E} \left[ \max_{j \in [N]: \ell \in A} \left| \frac{1}{N} \sum_{\ell \in A} (\theta_{j\ell} - p) \right| \right] \le \sqrt{2 \frac{|A|}{N^2} \log(2N)} \le \sqrt{\frac{2}{N} \log(2N)}. $$

The result follows by combining these last two inequalities with (39). $\square$

The lemma below explicitly contains Proposition 4.3.

LEMMA 8.2. There exists a constant $K > 0$ depending only on $\lambda$ such that for all $N \ge 1$ and $b \in \{-, +\}$,

$$ (40) \qquad \mathbb{E} \left[ \left\| L^{N, \bullet b} - bpr_b 1_N \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}, $$

$$ (41) \qquad \mathbb{E} \left[ \| L^N - p(r_+ - r_-) 1_N \|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}, $$
$$
(42) \quad \mathbb{E} \left[ \left\| C^N - p(1_{\mathcal{P}_+} - 1_{\mathcal{P}_-}) \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}},
$$

$$
(43) \quad \mathbb{E} \left[ \left\| \ell^{N, \bullet b} - \left( 1_{\mathcal{P}_b} + \frac{b(1-\lambda)p r_b}{1-(1-\lambda)p(r_+ - r_-)} 1_N \right) \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}},
$$

$$
(44) \quad \mathbb{E} \left[ \left\| \ell^N - \frac{1}{1 - (1 - \lambda)p(r_+ - r_-)} 1_N \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}},
$$

$$
(45) \quad \mathbb{E} \left[ \left\| m^N - m 1_N \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}.
$$

$$
(46) \quad \mathbb{E} \left[ \left\| v^N - m(1-m)1_N \right\|_{\infty} \right] \le K \sqrt{\frac{\log(N)}{N}}.
$$

REMARK 8.3. The asymptotics above are related to (Chevallier, Löcherbach and Ost, 2024, Lemmas 4.2. and 4.4.). There are two differences: 1) the 2-norm there is replaced by the infinity one here, 2) the rate of order 1 there is replaced by a rate of order $(\log(N)/N)^{1/2}$ here (remind the natural upperbound $\|x\|_2 \le \sqrt{N}\|x\|_{\infty}$ for $x \in \mathbb{R}^N$).

PROOF. Let us start with the proof of (40). By definition, for all $i \in [N]$ and $b \in \{-, +\}$,

$$
L_i^{N, \bullet b} - bpr_b^N = \sum_{j \in \mathcal{P}_b} (A^N(i, j) - bp) = \frac{b}{N} \sum_{j \in \mathcal{P}_b} (\theta_{ij} - p).
$$

The independent random variables $\theta_{ij} - p$ are bounded by 1 so that they are SubG(1/4). In turn, $L_i^{N, \bullet b} - bpr_b^N$ is sub-Gaussian with variance factor $r_b^N/4N \le 1/4N$. We consider the maximum of $2N$ of such variables (the factor 2 comes from the positive and negative part) so that the maximal inequality (Boucheron, Lugosi and Massart, 2013, Theorem 2.5.) implies

$$
\mathbb{E} \left[ \left\| L^{N, \bullet b} - bpr_b^N 1_N \right\|_{\infty} \right] \le \sqrt{2 \frac{1}{4N} \log(2N)} = \sqrt{\frac{\log(2N)}{2N}}.
$$

hence the result by using the fact that $|r_b^N - r_b| \le KN^{-1}$. Then, Equation (41) comes from the fact that $L^N = L^{N, \bullet +} + L^{N, \bullet -}$ and the triangle inequality. Moreover, the proof of Equation (42) is similar and therefore omitted.

We now turn to (44). Let us denote $\varepsilon^N = L^N - p(r_+ - r_-)1_N$. Starting from $1_N = Q^N[I_N - (1-\lambda)A^N]1_N$, one gets that

$$
\ell^N = 1_N + (1 - \lambda)Q^N L^N = 1_N + (1 - \lambda)Q^N \left[ p(r_+ - r_-)1_N + \varepsilon^N \right],
$$

so that

$$
\ell^N = \frac{1}{1 - (1 - \lambda)p(r_+ - r_-)} (1_N + (1 - \lambda)Q^N \varepsilon^N).
$$

Hence,

$$
\left\| \ell^N - \frac{1}{1 - (1 - \lambda)p(r_+ - r_-)} 1_N \right\|_{\infty} \le (1 - \lambda) \|Q^N\|_{\infty} \|\varepsilon^N\|_{\infty}
$$

and (44) follows from (41). Along the same lines, one can prove that

$$
\ell^{N, \bullet b} = 1_{\mathcal{P}_b} + b(1-\lambda)p r_b \ell^N + (1-\lambda)Q^N \varepsilon^N,
$$

and in turn that (43) follows from (44).
Starting from (17) and using the fact that

$$ \frac{\mu}{1 - (1 - \lambda)p(r_+ - r_-)} 1_N + 1_{\mathcal{P}_-} - \left( 1_{\mathcal{P}_-} - \frac{(1 - \lambda)pr_-}{1 - (1 - \lambda)p(r_+ - r_-)} 1_N \right) = m1_N, $$

one easily deduces (45) from (43) and (44).

Finally, remark that

$$ ||(m^N)^2 - m^2 1_N||_\infty \le 2||m^N - m1_N||_\infty, $$

so that (46) follows from (45).

<figure><img src="image_7.png" /></figure>

REFERENCES

ABBE, E. (2018). Community Detection and Stochastic Block Models: Recent Developments. *Journal of Machine Learning Research* **18** 1–86.

ABBE, E., BANDEIRA, A. S. and HALL, G. (2016). Exact Recovery in the Stochastic Block Model. *IEEE Transactions on Information Theory* **62** 471–487. https://doi.org/10.1109/TIT.2015.2490670

AIROLDI, E. M., BLEI, D., FIENBERG, S. and XING, E. (2008). Mixed Membership Stochastic Blockmodels. In *Advances in Neural Information Processing Systems* (D. KOLLER, D. SCHUURMANS, Y. BENGIO and L. BOTTOU, eds.) **21**. Curran Associates, Inc.

BERTHET, Q., RIGOLLET, P. and SRIVASTAVA, P. (2019). Exact recovery in the Ising blockmodel. *The Annals of Statistics* **47** 1805 – 1834. https://doi.org/10.1214/17-AOS1620

BOUCHERON, S., LUGOSI, G. and MASSART, P. (2013). *Concentration Inequalities: A Nonasymptotic Theory of Independence*. Oxford University Press. https://doi.org/10.1093/acprof:oso/9780199535255.001.0001

BUNEA, F., GIRAUD, C., LUO, X., ROYER, M. and VERZELEN, N. (2020). Model assisted variable clustering: Minimax-optimal recovery and algorithms. *The Annals of Statistics* **48** 111 – 137. https://doi.org/10.1214/18-AOS1794

CHEVALLIER, J., LÖCHERBACH, E. and OST, G. (2024). Inferring the graph of dependencies density of binary graphical models in high dimension. *Arxiv*.

DUARTE, A., GALVES, A., LÖCHERBACH, E. and OST, G. (2019). Estimating the interaction graph of stochastic neural dynamics. *Bernoulli* **25** 771–792.

GALVES, A. and LÖCHERBACH, E. (2013). Infinite Systems of Interacting Chains with Memory of Variable Length—A Stochastic Model for Biological Neural Nets. *Journal of Statistical Physics* **151** 896–921.

HOFFMANN, T., PEEL, L., LAMBIOTTE, R. and JONES, N. S. (2020). Community detection in networks without observing edges. *Science Advances* **6** eaav1478. https://doi.org/10.1126/sciadv.aav1478

LEI, J. and RINALDO, A. (2015). CONSISTENCY OF SPECTRAL CLUSTERING IN STOCHASTIC BLOCK MODELS. *The Annals of Statistics* **43** 215–237.

LIN, K. Z. and LEI, J. (2024). Dynamic clustering for heterophilic stochastic block models with time-varying node memberships. *Arxiv*.

MATIAS, C. and MIELE, V. (2016). Statistical Clustering of Temporal Networks Through a Dynamic Stochastic Block Model. *Journal of the Royal Statistical Society Series B: Statistical Methodology* **79** 1119–1141. https://doi.org/10.1111/rssb.12200

MOSSEL, E., NEEMAN, J. and SLY, A. (2015). Reconstruction and estimation in the planted partition model. *Probability Theory and Related Fields* **162** 431–461. https://doi.org/10.1007/s00440-014-0576-6

MOSSEL, E., NEEMAN, J. and SLY, A. (2016). Belief propagation, robust reconstruction and optimal recovery of block models. *The Annals of Applied Probability* **26** 2211 – 2256. https://doi.org/10.1214/15-AAP1145

OST, G. and REYNAUD-BOURET, P. (2020). Sparse space-time models: Concentration Inequalities and Lasso. *Annales de l'Institut Henri Poincaré (B) Probabilité et Statistique* **56** 2377–2405.

PEIXOTO, T. P. (2019). Network Reconstruction and Community Detection from Dynamics. *Phys. Rev. Lett.* **123** 128301. https://doi.org/10.1103/PhysRevLett.123.128301

PENSKY, M. and ZHANG, T. (2019). Spectral clustering in the dynamic stochastic block model. *Electronic Journal of Statistics* **13** 678 – 709. https://doi.org/10.1214/19-EJS1533

ZHAO, Y., LEVINA, E. and ZHU, J. (2012). Consistency of community detection in networks under degree-corrected stochastic block models. *The Annals of Statistics* **40** 2266 – 2292. https://doi.org/10.1214/12-AOS1036