Learning Large-Scale MTP$_2$ Gaussian Graphical Models via Bridge-Block Decomposition
Authors: Xiwen Wang
Jiaxi Ying
Daniel P. Palomar
Venue: NeurIPS 2023 poster (2023)
Keywords: MTP2 Gaussian Graphical Model
High-dimensional precision matrix estimation
Bridge-block decomposition.

Abstract: This paper studies the problem of learning the large-scale Gaussian graphical models that are multivariate totally positive of order two ($\text{MTP}_2$). By introducing the concept of bridge, which commonly exists in large-scale sparse graphs, we show that the entire problem can be equivalently optimized through (1) several smaller-scaled sub-problems induced by a \emph{bridge-block decomposition} on the thresholded sample covariance graph and (2) a set of explicit solutions on entries corresponding to  \emph{bridges}. From practical aspect, this simple and provable discipline can be applied to break down a large problem into small tractable ones, leading to enormous reduction on the computational complexity and substantial improvements for all existing algorithms.  The synthetic and real-world experiments demonstrate that our proposed method presents a significant speed-up compared to the state-of-the-art benchmarks.

Learning Large-Scale MTP, Gaussian Graphical Models via Bridge-Block Decomposition
Xiwen Wang', Jiaxi Ying 1.2 Daniel P. Palomar'
The Hong Kong University of Science and Technology
HKUST Shenzhen-Hong Kong Collaborative Innovation Research Institute (xwangew. jx.ging @connect.ust.hi palomar@ust.hk
Abstract
This paper studies the problem of learning the large-scale Gaussian graphical models that are multivariate totally positive of order two (MTP2). By introducing the concept bridge, which commonly exists large-scale sparse graphs, we show that the entire problem can be equivalently optimized through (1) several smaller- scaled sub-problems induced by a bridge-block decomposition on the thresholded sample covariance graph and (2) set of explicit solutions entries corresponding to bridges. From practical aspect, this simple and provable discipline can be applied to break down a large problem into small tractable ones, leading to enormous reduction on the computational complexity and substantial improvements for all existing algorithms. The synthetic and real-wor world experiments demonstrate that our proposed method presents significant speed-up up compared to the state-of-the -ftthe-art benchmarks.
Introduction
In recent decades, the surge in data availability has posed challenges for analyzing and understanding large-scale datasets. key challenge is examining pairwise relationships among numerous vari- ables, which can be represented using Gaussian graphical models (GGMs) [1 that depict variable connections as graphs. The precision matrix, which is the inverse of the covariance matrix, helps determine the non-zero patterns of GGMs 3 A traditional approach estimate the precision matrix is the graphical lasso 4 [}. which is formulated as a regularized log-determinant program. The solutions of graphical lasso possess an appealing property known as sparsity, which severs as common assumption particularly in large graphical models and has been shown to offer numerous benefits by significant amount of research. For example, the sparsity can reduce the model size and improve the interpretability of the model by highlighting the most important variables and their relationships, allowing better understanding the underlying causal structure 6-8).
In this paper, we solve the large-scale sparse precision matrix estimation problem for GGMs that fol- low multivariate totally positive of order two (MTP2) Gaussian distribution 10] equivalently, possess a precision matrix that is a symmetric M-matrix (II), where all off-diagonal elements are non-positive. The so-called MTP; property a special form of positive dependence and plays an essential role in applications where all interaction potentials are considered non-negative or the focus on emphasizing positive associations between variables [12155]. The MTP, property has been applied in various fields. Here are some examples. In finance, MTP, structures are often exploited as the asset returns are often considered positively correlated [16 In machine learning, MTP, GGMs are recognized as attractive Markov random field and have been used in applications such as taxonomic reasoning [19 and psychometrics
13 The MTP, precision matrix, also referred to
Corresponding author.
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

along with the combinatorial graph Laplacian [23 123 have found broad applications across a variety of fields due to their unique characteristics. These applications include, but are not limited to, graph Fourier transform 28 electrical circuits analysis [29]. image and video coding 30
as the generalized graph Laplacian 20 22
and structured graph financial time-series analysis 31
learning [33 (4].
The objective of this paper build a theoretical foundation and devise effective approaches for learning large-scale, sparse MTP₂ GGMs. The contributions of this paper are threefold.
This is the first work in the literature that introduces the notion of bridges in the context of learning MTP₂ GGMs to the best of our knowledge. Building upon this notion, we prove that the presence of explicit solutions for some entries depends whether an edge functions as bridge. Meanwhile, we demonstrate that the optimal solution exhibits decomposed structure through vertex-partition known as bridge-blod decomposition.
Based theoretical findings, we propose efficient framework that decomposes a large problem into several small tractable ones, each of which can be solved using any existing algorithm. With some negligible extra cost for bridge-block decomposition, the proposed method results in dimension reduction that significantly cuts down the computational complexity.
The synthetic and real-world experiments demonstrate that our proposed methods prompt consid- erable speed-up for all existing methods and enable the solving of large-scale MTP₂ GGMs that were previously considered infeasible.
1.1 Notation and Organization
Vectors and matrices are written as lower and upper case bold letters, respectively. An undirected graph is denoted G (V, E), where V is the set of nodes with size [V] p and & V x V is the set of edges. Note that we shall not distinguish between (i, E G and (i.j) E & Some graph terminology frequently used throughout the paper introduced as follows.
Support graph Given symmetric matrix A SP, the support graph of A, denoted as supp (A). is defined as undirected graph with the vertex set V = (l,.... P)} and the edge set & E V x V such that (i,j) E & if and only Aij ff 0 for every two different vertices i.j V.
Neighbors The neighbors of a node refer to the subset of vertices that connected to this node, i.e. (i) = { (i, E E}.
Path and Cycle A path from node i E V to node iT E V is a sequence (i.e., ordered set) of edges denoted di. each one incident the next, i.e., {(ir, de4 )) S &, A cycle is path from a node to itself when it does not duplicately include the same edge.
Partition A partition
(V....., Vx is a collection of non-empty, disjoint vertex sets Vi,..., Vx called cluster such that - V V.
Component A component of G refers to subset of vertices that connected to each other by edges, but are not connected to any other vertices in the graph.
The paper is organized as follows. In Section 2 we introduce the problem formulation and related works. In Section we present our main result along with potential impact. Section numerical experiments to show that the proposed methods can facilitate the performance of existing algorithms and solve large-scale MTP, GGMs problems In Section 5 we draw the conclusions.
4 presents the
Problem Formulation and Related Works
2.1 Problem Formulation
Let y be a random vector that follows a Gaussian distribution with zero mean and a precision matrix O, a.k.a. the inverse of covariance matrix Σ. A Gaussian graphical model represents the conditional dependency between random variables via a graph G, in which nodes correspond to y and edges between these variables represent conditional dependencies, which can be equivalently characterized by the non-zero patterns of the precision matrix + 2
Codes available in http::///gttuuuooo////ee1111///2 bbd

---

This paper considers estimating the precision matrix @ given n independent and identically distributed observations ........ Yn that follow an 2TP, Gaussian distribution. Equivalently. O is assumed to be symmetric M-matrix, i.e., iij S O for any + j. The problem is formulated as
minimize log det (e) + ( . + Eifj As 10yj, OEM²
where S is the sample covariance matrix, A.j 2 0 are the regularization coefficients, the objective to minimize the negative oigeliill of the data subject to a weighted fnorr penalty on the precision matrix, and M refers to a set of --maricces with dimension i.c.,
M = O E S" O
O and :;j S 0, }
(2)
We refer to O' as the optimal solution of (T) and to the support graph of @³, i.e., supp (0"), as the optimal graph. Particularly, we define the thresholded sample covariance matrix T as
Sig Agj if i and Sij > Aig.
To
otherwise,
and the support graph of T as the thresholded sample covariance graph (short for thresholded graph). Clearly, entries corresponding to Siy S Aæ will be thresholded to zero. Since the thresholded graph holds immense importance throughout the paper, unless explicitly mentioned otherwise, we denote G' = supp (T).
The involvement of MTP₂ constraints, which require that O, I 0 for all i + j. confers several advantages [[I]. For example, all M-matrices are inverse-positive, i.e.. O 2 and the non- smoothness from & norms can be eliminated via Σiy; A.y = (A, b), where A = (Ac)) and diag (A) = 0. Meanwhile, the MTP₂ constraints maintain the sparsity in the estimated precision matrix as an implicit regularizer [9] Throughout the paper, we require the following assumption.
Assumption 2.1. For any * j. we have Sij < SS
Under that assumption, the minimizer of Problem ( exists and is unique according to 19 Theorem 1].
Assumption 2.1 holds with probability if the number of observations n 2 2 13 19
2.2 Related Works
In literature, devising algorithms for learning MTP, GGMs has garnered considerable attentions. For instance, block coordinate descent (BCD) methods 19 26 351 update a single row/column of the precision matrix in a cyclic manner by solving nog-nega least squares problems. Projection- based methods, such as projected gradient descent (36] and projected quasi-Newton algorithms 377, iteratively take steps along the descent direction and then project the solutions back to the feasible region. Despite these efforts, existing research has not fully explored the potential of exploiting MTP, properties to reduce computational complexity. With complexities of Opp²) for BCD methods and (pp") for projection projection-based techniques, addressing large-scale problems continues to be challenging, particularly when manipulating full-dimensional matrices without problem reduction.
Instead of devising efficient algorithms, recently, leveraging the properties of the thresholded sample covariance graph has emerged as a popular approach for learning GGMs. Specifically, the existence of closed-fome solution has been established for graphical lasso when the thresholded sample covariance graph is acyelic (i.e., contains no cycles) [38] 39]. However, the applicability of their main results is limited by certain conditions that are difficult to verify. Interestingly. those conditions are unnecessary for the existence of closed-froo from solutions in MTP₂ GGMs 20). Despite the theoretical establishment on closed-form solutions, an exact acyelic structure is rarely observed in practice. Therefore, more research dives into the relationship between supp (T) and supp (0*). Research found that supp (0*) is subset of supp (T) [3), the number of connected components in supp (T) is identical to that in supp ( [20]. and there exist necessary conditions for the presence of edges in supp (0*) via analyzing the path products on thresholded matrix [13].
This paper advances prior research in two ways. Firstly, unlike previous studies that only provided an explicit solution for O.y in the case of acyclic thresholded graphs, we reveal that this explicit solution for G.y consistently applies to all (i. pairs acting as bridges, regardless of whether the thresholded graph is acyclic or non-acyclic. Secondly, we highlight that the optimal graph can be represented in an equivalent decomposed form through a vertex partition, termed as the bridge-block decomposition of the thresholded graph by leveraging MTP₂ properties.
3

---

Figure 1: An illustration of the bridge-block decomposition. Edges (5.6) and (9, 10) are identified as bridges since removing either of them increases the number con- nected components. After removing the bridge edges, the resulting decomposition phld ((1,2,3,4,6), 87,8, 9, 10) (11, 12, 13, 14, 15, 16)).
Proposed Methods
The goal of this paper is about developing an algorithm but shed light the remarkable properties concealed in the thresholded graph. This section is organized as follows. In Section[..] introduce bridge and bridge-block decomposition. Section 3.2 presents main result. Then in Section 3.3
we elaborate on the contributions and connections to existing research.
3.1 Bridge and Bridge-Block Decompusition
Bridge is one of the important concepts in graph theory. Technically, edge is called bridge if and only if its deletion increases the number of graph components. Therefore, an edge is a bridge only when it is not contained in any cycles [40]. Notably, when graph consists solely of bridges, it is referred to as acyelic. this paper, we denote B as the set of all bridges, i.e,
B { (i,j)) (i.j)in is bridge in G}-
Remark 3.1. Bridges are frequently observed, particularly in large-scale sparse graphs [1143] This attributed the fact that in sparse graphs, only the most significant relationships between variables are retained, and the removal of edges can create additional connected components, giving rise to the presence of bridges.
Remark 3.2. In practice, the set can be efficiently obtained via various bridge-finding algorithms 40 44). These algorithms employ a depth-first search approach, resulting in computational complexity of O ([V] El). In the sparse graphs that are of interest to us, the number edges |E] typically scales similarly to the number of nodes VV As result, the computational cost of identifying bridges in large-scale sparse graphs remains low.
Using Figure [ as an illustration, we perform the bridge-block decomposition 45 as follows. By removing all the bridges, we compute the clusters V& corresponding to the components of G -V, B). This process results a pertex-partition,, known as the bridge-block decomposition:
phbl
(V1, V2...., Vk}.
(5)
where K refers to the number of clusters, also the number of components in the graph (V, B). \\B).
Without loss of generality, we define the operator as (i) l if node belongs the k-th cluster. In the literature, the bridge-block decomposition is also known as the 2-edge-connected
decomposition [46] 47 fundamental concept in graph theory with numerous applications such as community search 48 social network mining [49], and transmission networks [50].
In practice, the of obtaining the bridge-block decomposition pbhdd which involves calculating the thresholded graph, bridges, and clusters, negligible. With the aforementioned preliminary knowledge, we formally introduce main result as follows.
3.2 Main Results
Considering square matrix SP, we define Av, SP as the principal sub-matrix of keeping the rows and columns indexed by Vk. in which pk Vill is the number of nodes in k-th cluster and

---

+ (1223) θ + diag (2,, Sa).
(4.5,6) + + diag 53.
3
(i,j) E B: O.j = T (SirS - T2y)
(i.j) $ and (b) + (j): O1 = 0.
Figure 2:
An example of obtaining optimal solution via Theorem 3.3 The thresholded graph of
nodes can be partitioned into 2 clusters V1 {1, 2, 3} and V/2 = - 5, The optimal solution
O 8 can be exactly computed via set of explicit solution {03,1, 11,3). solutions of
smaller-sized sub- problems i.e., θ E 3 and @ E $² and zeros in the rest positions.
have p Σx Ph. For each E VR, mark (i) as its corresponding index in Vg- Hence, define 3 as the optimal solution of k-th sub-problem i.e.,
= arg min OLEMPE
log det (0x) + (0x, Sv. - Av.)-
The main result is then given as follows.
Theorem 3.3. Under Assumption 2.1 given the bridge-block decomposition the the thresholded graph
supp (T) as pold, and the optimal solution of each sub-problem 446 as ② the optimal solution of Problem ie.. O+, can be obtained
+ (()))=((() + Gi
i Vk: i and i,j E Vk.
+ ())T(())(())
1
(7)
(S - -T2))
i((i,j) E B, otherwise.
which 5. Su E ((im)EB SaSmm T2
Th
and C im (i,m) &
Figure 2 depicts an example of how to apply Theorem 3.3 to obtain the optimal solution more efficiently. Theoretically, the key to show the optimality of via an explicit expression of the inverse of θ i... R = θ using following theorem.
Theorem 3.4. Given S.T T E SP, the bridge-block decomposition phhd -V1..... Vk3 of supp (T), and set of matrices ( E 1k=1 l
the inverse of @ in the form of is derived
if E Vk. (f. E otherwise,
Raj
Ti
(8)
Su,S Di (R)
where R2 [0p]= and for each i, different clusters, gij given as
9ú R
2T = Rug,ets/ Su where un i, u22+1 A j,
(9)
in which T - 1 is the number of bridges in dij, By. 0 dij G. and denotes a sequence of following the orders of diş while
'incident' bridges, i.e., diş n B = {(u11 U2) ........ (12T-11
u2r) (
preserving the elements that are bridges.
We note that all terms in (9) can be computed via ( and Theorem 3.4 itself generally holds for arbitrary S, T. and
{ 1K Am Detailed proofs and discussions are deferred to Appendix A
Theorem 3.4 is theoretically critical
provides explicit form of the inverse of O which was previously difficult to compute. is vital for demonstrating the optimality of θ by connecting the KKT system of the original large-scale problem the KKT systems of smaller sub-prohe sub-problems through the bridge-bloct block decomposition. Consequently, together with MTP, properties, we show that in the form of satisfies the KKT condition of Problem ① in Appendix
The role of MTP, constraints: Although Theorem 3.4 can be applied without requiring MTP, properties, is important to highlight that these properties serve sufficient conditions to demonstrate
5

---

the optimality 0 - R By incorporating the MTP₂ constraints, the non-smoothness in graphical lasso is eliminated, resulting in simplified optimality conditions shown below:
Vi: Ru: + S =
RRi S =
(10a)
vOiy # 0: Raj + Sij + Xij sign( = 0
- R.j + Sij Xij 0
(10b)
,0,, 0: --Ra + Siyl S As)
Rig + S.g - Au S (10c)
Simultaneously, becomes non-negative (2 0). As a result, The most challenging part of the KKT conditions- - -R;; Sij Aiy S 0 holds under MTP₂ properties. In Appendix C we will elaborate the details and explore additional conditions under which we can extend the applicability of Theorem 3.3 to the traditional graphical lasso.
In practical applications, it is often more efficient to employ the bridge-block decomposition technique instead of directly manipulating full-sized matrices. This approach involves breaking down the problem into smaller isolated sub-problems. By addressing these sub-problems individually, we can compute the optimal solution by utilizing the solutions obtained from each sub-problem, as outlined in Theorem 3.3 This approach offers several advantages:
Propused solving framework:
Significant reduction in computational cost. Foremost, the cost of bridge-block decomposition is cheap. Suppose that we use BCD solver of complexity (p*)), then the total cost is reduced to Ek= O(1) ()) < Op²"), where Σx [i
This can prompt an enormous difference.
Considerable reduction in memory cost. The memory cost is typically troublesome for large- scale problems as each full-dimensional matrix contains p² elements. Theorem 3.3 can avoid generating number of full-dimensional intermediate variables during computation.
Potential speed-up via parallel computing. The sub nub-problems can be optimized independently, which allows parallel computing for significant speed-up.
3.3 Connection Existing Research
From Theorem 3.3 we obtain a very interesting result that the (i. j-th entries of θ admits an explicit solution, i.c., O,5 = -Ty/ (S - 22 if edge (i,j) is a bridge in the thresholded graph supp (T). To the best of our knowledge, this result, shown in Corollary 3.5 is the first sufficient condition for edge belonging optimal graph supp (O*).
Corollary 3.5. edge (i. is bridge in the thresholded graph supp (T), then (i.j) is also bridge in the optimal graph supp (0*).
Corollary 3.5 can be obtained as follows. By the definition of a bridge, from nodes to j the path dij {i.j)} is the only path in the thresholded graph, and removing would result in increase in the number of components. Given that supp θ*) S sT)p( and θ4 * 0, it follows that dy (i,j)} remains the unique path from nodes i to in the optimal graph. This indicates that the edge (i,j) continues to function as a bridge.
In the literature, only some necessary conditions for edges in the thresholded graph supp (T) to be retained in the optimal graph supp (0+) have been established. For instance, [13] and 2 indicate that supp (0)) C supp (T), which implies that (i, & supp (+*) ) if (i,j) & supp (T))
Another explicit solution is also applicable when k-th cluster only contains node, i.e., V/2 (i). then θ = 1/5,, As consequence, Theorem 3.3 generalizes to the following corollary for acyclic graph which admits bridge-block decomposition into clusters.
Corollary 3.6. Suppose that the thresholded graph supp (T) is acyclic, then the optimal graph supp (+*) is also acyclic and O* admits a closed-form solution as
(1+EmENrHH EmEN(()
if= = j.
S..
+
Ty
if (1.)) EB otherwise.
(1))
Su,S 13
0
The proof directly follows Theorem Compared to our results in Theorem
3.3
Corollary 3.6 coincides with 20, Theorem 3] and [38
3.3
explicit solutions in literature heavily depend on the graph

---

structure. Hence, our theory appears as the first to reveal that, fundamentally, it is the edge property. ie., whether the edge is bridge, that decides the existence of explicit solutions.
In conclusion, our proposed methods generalize the existing results with more profound understand- ings, offering significant potential for solving large-scale precision matrix estimation problems.
4 Numerical Simulations
We conduct experiments synthetic and real-world data to evaluate how the proposed method accelerates the convergence of existing algorithms compared to algorithms that do not exploit it. All experiments were conducted on 2.60GHZ Xeon Gold 6348 machines and Linux OS. All methods are implemented in MATLAB and the state-of-the-art methods we consider includes
BCD: Block Coordinate Descent 19
of complexity Opp*).
PGD: Projected Gradient Descent 36 of complexity O(p²).
PQN-LBFGS: A projected quasi-New Newton method with limited memory BFGS [51] The complexity is (( + )) where m is the number of iterations stored for approximating the Hessian.
FPN: Fast projected Newton-like method [37 of complexity Opr").
Note that all methods converge the optimal solution. The target here is not to compare these methods but to evaluate how our proposed framework promotes the efficiency of these methods.
4.1 Synthetic Data Experiments
We use the processes described in [9] to generate the data. We begin with an underlying graph that has an adjacency matrix A E S", and define @ = 5I A, where 5 = 1.05 Amax (A) and A. his ensures that O is a positive definite matrix with off-diagonal elements being negative, making @ randomly generated M-matrix. Next, we normalize by substituting with DOD, where is suitably chosen diagonal matrix, ensuring that diag( O @-1)) 1. We then sample = 10p data points from a Gaussian distribution (0, θ O-1) and calculate the sample covariance matrix as S.
Amax (A) represents the largest eigenvalue
Following (25). construct the regularization matrix A. Briefly, given initial estimate θ set A, = =x *j + e) when * j and A. = 0 when = Here, > determines the sparsity level and is a small positive constant, such as 10-3 Generally, a large penalty is expected on θ; when @000 is small. We efficiently compute (0) using the following equation:
0 = Tiy' (SnS - - Vi *
(12)
and then appropriately adjust the value of X that supp (T) has only one connected component, while supp (@*) can roughly recover the underlying structure.
We consider two scale-frec random graphs underlying with their typical structures in Figure 4
Barabasi-Albert (BA) graph 52 53 of order one. BA models generate random scale-free networks via preferential attachment, suitable for modeling various networks like the Internet, world wide web, protein interactions, citations, and social/online networks.
2. Stochastic Block Model (SBM) of networks [54), a.k.a. the community graph 555 Stochastic block models serve as fundamental tools in network science, creating random networks based on community structures, where nodes within the same group more likely to form connections.
In the experiments, we consider = 5000 and p 20000 and calculate the relative errors as follows:
RE ( |/ (0) - / (0*)
(13)
where we compare the relative errors at each iteration against the computational time. Here, θ represents the optimal solution, and signifies the objective function of ( For methods employing bridge-block decomposition (bbd), we first calculate 6
(*) for each cluster at the s-th iteration, followed by calculating the intermediate solution (s) using 77. The computational time at the s-th
iteration is then the sum of the costs to obtain pbbd calculate all @ and compute (*) via
7

---

10
10
Relabbva Ervora
Realativa Emar
12
10
time (seconds)
time (seconds) BA graph, 20000
graph, 5000
107
Rellaivve Errors
Rullative Errors
10'
CPU time (seconds) SBM graph, 5000
CPU time (seconds) SBM graph, p 20000
BCD
PGD
FPH
PQN-LBFGS
BCD bbd)
PGD bbd) -0- FPN bod) --- PQN-LBFG5 bbd)
Figure 3: Relative errors of the objective values versus the computational time. Colors are used for distinguishing the state-of-the-art methods. Solid circle symbols stand for methods accelerated by bridge-block decomposition (bbd) and square symbols stand for methods without acceleration.
The results for - 5000 and p 20000 are displayed in Figure with all results averaged over five realizations. We highlight that the extra computational cost compared methods without acceleration is negligible shown in Table From the figures, we observe that:
Significant Speed-up: Our proposed framework substantially accelerates state-of-ther the-ttrrrtt methods by one to four orders of magnitude.
When = 20000, existing methods cannot optimize the problem within 10 seconds. In contrast, our framework greatly speeds up convergence, making possible to solve large-scale problems.
Solving Otherwise Infeasible Problems:
Higher Speed-up for Sparser Graph: Proposed method achieves greater acceleration on the BA graph than the SBM graph, as the former exhibits stronger sparsity pattern.
Higher Speed-up for Methods of Higher-complexity: Figure demonstrates that our method has more significant impact improving the efficiency of the BCD method, which has higher computational complexity and thus benefits more from dimension reduction.
To further shed light the factors that determine the magnitudes of improvement, in the next sub-section, we conduct additional experiments that dive into how bridge-block decomposition boosts the performance of existing methods subject different structures of the thresholded graph.
Table 1: Average computational time (seconds) of extra cost.
Conducting Bridge-block decomposition
Computing O from {
p 5k
p 20k
p 5k
p 20k
BA
0.2
2.2
0.0
0.0
SBM
0.3
2.3
0.1
0.2

---

Ratto
Proberrs Dimamalora (p)
(1) BA graph
(2) SBM graph
Figure 5: Ratios of improvements for BCD method with respect to computational time.
Figure 4: Structures of BA graph and SBM graph.
4.2 Further Experiments on Ratios of Improvement
investigate scenario where the underlying graph a community graph with a block-tridiagonal adjacency matrix, and each cluster contains an equal number of nodes (pi [Vll, vk). Without loss of generality, we assume that the internal edges within clusters form cycle. To generate S and A, we follow the procedures in Section 4.1 We then modify A (11" A) (1 a) A, ensuring [0, 1] and A is the adjacency matrix of the considered graph. that supp (T) = supp (A), where E
We experiment with various values of K (the number of clusters) and pl (the number of nodes in each cluster). For each configuration, we calculate the ratio of improvement, defined as the quotient of the time required achieve RE (+) < 10-6 without using bridge-block decomposition to the time when applying the decomposition. We conduct multiple trials of the BCD method and display the results in Figure 5 The results Figure Suggesst that the number of sub sub-graphs K is primary factor influencing the ratio values. The improvement is deemed tremendous, especially for large-scale settings. The involvement of bridge-block decomposition enables solving an MTP; GGMs with hundreds of millions of variables as long we can conduce bridge-blow block decomposition. Consequently, many previously unattainable real-life applications can now be optimized.
4.3 Real-word Data Experiment
We consider learning the MTP₂ GGM for the Crop Image dataset available from the UCR Time Series Archive [56). This dataset comprises p 24000 pixels, where each pixel corresponds time series record capturing spectral information in 46 instances. These instances represent geometrically and radiometrically corrected satellite images. By examining this dataset, we can observe the temporal changes in the observed areas through the recorded measurements. To facilitate our analysis, the data has been pre-processed by [57 using an indexing technique. This preprocessing step enables us to work with compressed informative indexes instead of handling large image files.
Our goal is to perform graph-based clustering the time series data that contains 46 observations, denoted as {y1,.... } using MTP₂ GGMs, where y, E 1220000 The data includes 24 classes corresponding to different indexed land covers, such corn, barley, or lake. It is important note that these labels are not involved in the clustering process but help evaluate the quality of the final graph-basse based clusters by computing the graph modularity [58]. Via analyzing the data, Appendix D we discuss that why learning MTP2 GGMs is statistically meaningful in the context of clustering and show that the MTP₂ assumptions approximately hold for this dataset.
The initial estimate (0) is computed using the method proposed in [59]. The regularization matrix is determined using the approach Section 4. with e = 0.01 and = 0.2. We then compute the sample correlation matrix as S E $24000 This results in precision matrix estimation problem involving 5.76 x 10 variables to optimize, which is quite challenging for state-of-the-an methods without bridge-block decomposition.
Figure shows the results of empirical convergence, and Figure visualizes a local structure of supp ) (modularity = 0.6849). In line with previous findings, our suggested approach greatly speeds up convergence of all current algorithms by minimum three orders of magnitude. This demonstrates immense potential for managing large-scale sparse graphical models.

---

Barley
Hardwood
18
Corn
Pasture
Ralatbee
14'44
Temporal Meadow
CPU time (seconds)
ECD
POR-LEPGS
Repeseed
(a) Visualizations.
(b) Experimental results Crop data (p 2-400).
Figure 6: Visualizations and experimental results of sparse MTP2 GGMs for Crop data set. (a) A local structure of optimal graph with 2008 nodes. Nodes with matching labels are assigned the same color and connected by matching edge color, while different groups of nodes are connected by gray edges. (b) Results of convergence to learn sparse MTP2 GGMs for Crop data set.
While our paper's ' primary focus not to uncover insights into the estimated MTP, GGMs for a deeper understanding of the underlying data, interesting phenomena still be observed in the structures presented in Figure 6a] Notably, we find that the majority of edges exist within the same crop type, while connections between nodes associated with different crops are relatively sparse. This observation is valuable for clustering processes and aligns with our expectations, as stronger positive dependencies are often observed within the same class, while dependencies between different classes tend to be weaker.
Moreover, graphical representation reveals more nuanced patterns. For instance, in Figure 6a of our manuscript, we observe dense network of edges between two crop types, 'temporary meadow' and 'pasture', indicating significantly stronger conditional dependency compared to other categories. This observation further supports our understanding. Hence, the inferred conditional dependency structure encapsulates the inherent interrelationships among different crops. As result, our proposed framework based on bridge-block decomposition effectively facilitates learning in large-scale MTP2 graphical models.
5 Conclusions and Discussions
Real-world sparse Gaussian graphical models often comprise subsets of variables that densely connected to one another, while variables different clusters maintain weak connections. However, standard estimation algorithms do account for this property. In this paper, by introducing the concept of bridge, we leverage these characteristics to reveal an interesting finding: the optimal solution for Problem (T) equivalently admits decomposed form via the bridge-block decomposition of the thresholded graph. We provide theoretical insights into the separability of optimal solutions and the existence of explicit expressions based on edge properties, specifically whether an edge is bridge. This method surpasses conventional approaches that depend on unique graph structures, such determining a graph acyclic. From the practical aspect, our proposed method offers a handy architect for accelerating any existing algorithms by handling number of small and tractable sub- coob--roblems.. Although our method is mainly developed for sparse large-scale graphs and therefore seems inapplicable for dense graph, as elaborated in the appendix C it can also provide quick way to obtain a solution, which can serve as either a starting point (warm-start) for numerical algorithms or alternative for acquiring approximate solutions. Overall, our simple and provable approach demonstrates exceptional performance and paves novel way for more extensive designs of large-scale MTP, Gaussian graphical models.
sarge-scale learning problem via
10

---

Acknowledgements
This work was supported by the Hong Kong Research Grants Council GRF 16207820, 163 16310620, and 16306821, the Hong Kong Innovation and Technology Fund (TTF) MHP/009/20, and the Project of Hetao Shenzhen-Hong Kong Science and Technology Innovation Cooperation Zone under Grant (-DBB-KCZYB-022000 We would also like to thank the anonymous reviewers for their valuable feedback on the manuscript.
References
[1] H. Rue and L. Held, "Gaussian Markov random fields: theory and applications, Monographs on Statistics and Applied Probability, vol. 104, 2005.
[2] S. L. Lauritzen, Graphical models. Clarendon Press, 1996, vol. 17.
[3] J. Whittaker, Graphical models in applied multivariate statistics. Wiley Publishing, 2009.
[4] O. Banerjee, L. El Ghaoui, and A. d' A'premoo,,, "Model selection through sparse maximum like- lihood estimation for multivariate Gaussian or binary data," The Journal of Machine Learning Research, vol. 9, pp. 485-516, 2008.
[5] J. Friedman, T. Hastie, and R. Tibshirani, "Sparse inverse covariance estimation with the graphical lasso," " Biostatistics, vol. 9, no. 3, pp. 432-441, 2008
[6] M. Yuan and Y. Lin, "Model selection and estimation in the Gaussian graphical model," Biometrika, vol. 94, no. 1, pp. 19-35, 2007.
[7] H. Liu, K. Roeder, and L. Wasserman, "Stability approach to regularization selection for high dimensional graphical models," Advances in Neural Information Processing Systems, vol. 23, 2010.
[8] Danaher, P. Wang, and D. M. Witten, "The joint graphical lasso for inverse covariance estimation across multiple classes," Journal the Royal Statistical Society: Series B (Statistical Methodology), vol. 76, no. 2, pp. 373-397, 2014.
[9] E. Bølviken, "Probability inequalities for the multivariate normal with non-negative partial correlations," Scandinavian Journal of Statistics, pp. 49-58, 1982.
[10] D. M. Malioutov, J. K. Johnson, and A. S. Willsky, "Walk-sums and belief propagation in Gaussian graphical models," The Journal of Machine Learning Research, vol. 7, pp. 2031-2064, 2006.
[11] R. J. Plemmons, "M-matrix characterizations. I-nonsingular M-matrices," Linear Algebra and its Applications, vol. 18, no. 2, pp. 175-188, 1977.
[12] Fallat, S. Lauritzen, Sadeghi, C. Uhler, Wermuth, and P. Zwiernik, "Total positivity in Markov structures," The Annals of Statistics, pp. 1152-1184, 2017.
[13] S. Lauritzen, C. Uhler, and P. Zwiernik, "Maximum likelihood estimation in Gaussian models under total positivity." The Annals of Statistics, vol. 47, no. 4, pp. 1835-1863, 2019.
[14] Röttger, S. Engelke, and Zwiernik, "Total positivity multivariate extremes," The Annals of Statistics, vol. 51, no. 3, pp. 962-1004, 2023.
[15] D. Rossell and P. Zwiernik, "Dependence in elliptical partial correlation graphs," Electronic Journal of Statistics, vol. 15, no. 2. pp. 4236-4263, 2021.
[16] A. Müller and M. Scarsini, "Archimedean copulae and positive dependence," Journal of Multi- variate Analysis, vol. 93, no. 2, pp.
2005.
[17] R. Agrawal, U. Roy, and C. Uhler, "Covariance matrix estimation under total positivity for
portfolio selection." Journal of Financial Econometrics, vol. 20, no. 2, pp. 367-389, 2022.
[18] R. Zhou, J. Ying, and D. P. Palomar, "Covariance matrix estimation under low-rank factor model with nonnegative correlations," IEEE Transactions on Signal Processing, vol. 70, pp. 4020-4030, 2022.
[19] M. Slawski and M. Hein, "Estimation of positive definite M-matrices and structure learning for attractive Gaussian markov random fields," Linear Algebra and its Applications, vol. 473, pp. 145-179, 2015.
11

---

[20] E. Pavez, H. E. Egilmez, and A. Ortega, "Learning graphs with monotone topology properties and multiple connected components." IEEE Transactions on Signal Processing, vol. 66, no. 9, pp. 2399-24492 2413, 2018.
[21] J. Ying. J. V. d. M. Cardoso, and D. P. Palomar, "A fast algorithm for graph learning under attractive Gaussian Markov random fields," in 2021 55th Asilomar Conference on Signals, Systems, and Computers, 2021, pp. 1520-1524.
[22] Z. Deng and A. M.-C. So, "A fast proximal point algorithm for generalized graph Laplacian learning," in ICASSP 2020-2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2020, pp. 5425-5429.
[23] J. Ying, J. V. d. M. Cardoso, and D. P. Palomar, "Does the & noom learn a sparse graph under Laplacian constrained graphical models?" arXiv preprint arXiv: 2092,5,,00. 2020.
[24]
"Nonconvex sparse graph learning under Laplacian constrained graphical model," in
Advances in Neural Information Processing Systems, vol. 33, 2020, pp. 7101-7113.
[25]
"Minimax estimation of Laplacian constrained precision matrices," in Proceedings of The 24th International Conference on Artificial Intelligence and Statistics, vol. 130, 2021, pp. 3736-3744.
[26] H. E. Egilmez, E. Pavez, and A. Ortega, "Graph learning from data under Laplacian and structural constraints," IEEE Journal of Selected Topics in Signal Processing, vol. 11, no. 6, pp. 825-841, 2017.
[27] J. Ying, X. Han, R. Zhou, X. Wang, and H. C. So, "Network topology inference with sparsity and Laplacian constraints," arXiv preprint arXiv:2309.0099: 2023.
[28] D. I. Shuman, S. K. Narang, P. Frossard, Ortega, and Vandergheynst, "The emerging field of signal processing on graphs: Extending high-dimensional data analysis to networks and other irregular domains," " IEEE signal processing magazine, vol. 30, no. 3, pp. 83-98, 2013.
[29] Dorfler and F. Bullo, "Kron reduction of graphs with applications to electrical networks," IEEE Transactions on Circuits and Systems I: Regular Papers, vol. 60, no. 1, pp. 150-163, 2012.
[30] H. E. Egilmez, Y.H. Chao, A. Ortega, B. Lee, and S. Yea, "GBST: Separable transforms based on line graphs for predictive video coding," in 2016 IEEE International Conference on Image Processing (ICIP). IEEE, 2016, pp. 2375-2379.
[31] J. V. D. M. Cardoso, J. Ying, and D. P. Palomar, "Graphical models in -tiled markets," in Advances in Neural Information Processing Systems, vol. 34, 2021, pp. 19 989-220011.
[32]
"Learning bipartite graphs: Heavy tails and multiple components," Advances in Neural Information Processing Systems, vol. 35, 2022, pp. 14 044-14057.
[33] S. Kumar, J. Ying, J. V. d. M. Cardoso, and D. P. Palomar, "A unified framework for structured graph learning via spectral constraints," Journal of Machine Learning Research, vol. 21, no. 22, pp. 1 1-60, 2020.
[34]
"Structured graph learning via Laplacian spectral constraints," in Advances in Neural Information Processing Systems, vol. 32, 2019, pp. 11 647-11 658.
[35] E. Pavez and A. Ortega, "Generalized Laplacian precision matrix estimation for graph signal processing," in 2016 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2016, pp- 6350-6354.
[36] J. Ying, J. V. D. M. Cardoso, and D. P. Palomar, "Adaptive estimation of graphical models under total positivity," " in International Conference on Machine Learning, 2023, pp. 40
[37] J.-F. Cai, J. V. d. M. Cardoso, D. P. Palomar, and J. Ying, "Fast projected Newton like method for precision matrix estimation under total positivity," arXiv preprint arXiv:2112.. 12.01939, 2021.
[38] S. Sojoudi, "Fquivalence of graphical lasso and thresholding for sparse graphs," The Journal of Machine Learning Research, vol. 17, no. 1. pp. 3943-3963, 2016.
[39] S. Fattahi and S. Sojoudi, "Graphical lasso and thresholding: Equivalence and closed-form solutions," " Journal of Machine Learning Research, 2019.
[40] M. Schmidt, "A simple test on 2-vertex-and 2-edge-connectivitty" Information Processing Letters, vol. 113, no. 7, pp. 241-244, 2013.
12

---

[41] S. Hojsgaard and S. L. Lauritzen, "Graphical Gaussian models with edge and vertex symmetries." Journal of the Royal Statistical Sociery: Series B (Statistical Methodology), vol. 70, no. 5, pp. 1005-1027, 2008.
[42] M. L.. Stein, "Space-time covariance functions," Jourrnal of the American Statistical Association,
vol. 100, 469, pp. 310-321, 2005.
[43] C. Li and H. Li, "Network-constraineed regularization and variable selection for analysis of genomic data," " Bioinformatics, vol. 24, no. 9, pp. 1175-1182, 2008.
[44] R. Tarjan, "A note on finding the bridges of a graph." Information Processing Letters, vol. 2, no. 6, pp. 160-161, 1974.
[45] J. Westbrook and R. E. Tarjan, "Maintaining bridge-connected and biconnected components on-line," Algorithmica, vol. 7, no. 1, pp. 433-464, 1992.
[46] L. Chang, J. X. Yu, L. Qin, X. Lin, C. Liu, and W. Liang, "Efficiently computing k-edge connected components via graph decomposition." " in Proceedings of the 2013 ACM SIGMOD International Conference on Management of Data, 2013, pp. 205-216.
[47] R. Zhou, C. Liu, J. X. Yu, W. Liang, B. Chen, and J. Li, "Finding maximal k-edge-connected subgraphs from a large graph." " Proceedings of 15th International Conference on Extending Database Technology, 2012, pp. 480-491.
[48] Y. Fang, X. Huang, L. Qin, Y. Zhang, W. Zhang, R. Cheng, and X. Lin, "A survey of community search over big graphs," The VLDB Journal, vol. 29, no. 1. pp. 353-392, 2020.
[49] R. Agrawal, S. Rajagopalan, R. Srikant, and Y. Xu, "Mining newsgroups using networks arising from social behavior," in Proceedings of the 12th International Conference on World Wide Web, 2003, pp. 529-535.
[50] I.. Lan and A. Zocca, "Refining bridge-block decompositions through two-stage and recursive tree partitioning," " arXiv preprint arXiv:2110.06900 0000698, 2021.
[51] D. Kim, S. Sra, and S. Dhillon, "Tackling box-constrained optimization via a new projected quasi-Newton approach," SIAM Journal on Scientific Computing, vol. 32, no. 6, pp. 3563, 2010.
[52] A.-L. Barahási and R. Albert, "Emergence of scaling in random networks," Science, vol. 286, no. 5439, pp. 509-512, 1999.
[53] R. Albert and A.-L. Barabási, "Statistical mechanics of complex networks," " Reviews of Modern Physics, vol. 74, no. 1, p. 47, 2002.
[54] P. W. Holland, K. B. Laskey, and S. Leinhardt, "Stochastic blockmodels: First steps." Social Networks, vol. 5, no. 2, pp. 109-137, 1983.
[55] S. Fortunato, "Community detection in graphs," Physics Reports, vol. 486, no. 3-5, pp. 75-174, 2010.
[56] H. A. Dau, A. Bagnall, K. Kamgar, C.-C. M. Ych, Y. Zhu, S. G.aarghabi, C. A. Ratanamahatana, and E. Keogh, "The UCR time series archive," IEEE/CAA Journal of Automatica Sinica, vol. 6,
no. 6, pp. 1293-1305, 2019.
[57] C. W. Tan, G. I. Webb, and F. Petitjean, "Indexing and classifying gigabytes of time series under time warping," in SIAM International Conference on Data Mining, 2017, pp. 1-10.
[58] M. E. Newman, "Modularity and community structure in networks." Proceedings of the National Academy of Sciences, vol. 103, no. 23, pp. 8577-8582, 2006.
[59] F. Nie. X. Wang, M. Jordan, and H. Huang, "The constrained Laplacian rank algorithm for graph- based clustering," in Proceedings of the AAAI Conference on Artificial Intelligence. vol. 30, no. 1, 2016.
13

---

Appendix
In what follows, we present technical proofs omitted in the main content. Section A] we introduce important properties of matrix R and prove Theorem 3.4, In Section B we prove the main result, ie.. Theorem 3.3 In Section C we introduce some extension of the proposed framework to broaden its applicability. In Section D we will comment on real-world experiments and explain why MTP₂
assumption approximately holds for Crop Image dataset.
A Proof of Theorem 3.4
This section contains detailed description of how we construct the matrix R and prove Theorem 3.4 based on the properties of R. The proof relies on the concept of a 'generalized path product', which involves a set of edges labeled as a 'bridge path'. The matrix R has several notable properties that are essential in demonstrating Theorem 3.3
A.1 Definitions and Important Properties
To begin, let us define the 'bridge path' Bij for nodes i and that belong to different clusters, i.e., (t) A vi) based on the bridge set B in the thresholded graph G. The bridge path Biy is defined as the set of bridges in any path di from node to node j. where each element in B,g a bridge in B and the elements in B.y follow the order of dis while preserving the elements in B. Symbolically, we can express this as:
B₁ dij n B {(u1, u2) (uz, u4) .....: (U2T (,,T=1,
(14)
where T denotes the number of clusters in any path from nodes i to j. It is important to note the following:
Each element in B.y is bridge in the thresholded graph.
The elements in By follow the orders of dij. This means that for any path from nodes to j. the path will visit (Ugg-1. U2:), then visit (u25-11 U225) whenever we have i < consequently, we have
(i) = (un). & (uz) ab (us) ...... 1b (u2r-2) = vp (uzr-1), & (w2r) +b (f). (15)
For any path dij. the bridge path Bi, dip n B is unique. Biy is defined as empty set if there is path from nodes to j.
is important mention that the thresholded graph acyelic, the bridge path can reduced to a traditional path. In this scenario, all edges in the path d,j become bridges, and the bridge path Biş is equivalent to the path dey. which can be expressed as Bij = dej-
We define the matrix R E SP based on the definition of bridge-path, as shown below:
if E Ve. if (i.j) E B, otherwise,
Rig
Ty VSnSjj 96 (R)
(16)
where R2 [0p] and for each i, in different clusters, Pyy is given as
Jū R I Ruc, / / Sun, Sucei where uo A i, u2r+1 A j.
(17)
where u denotes the bridge path Bij and 9ij (R) = 0 if B.j = 0. Remark A. Clearly, we have Raj Rys for any and
Example: For instance, using Figure 7 as example, we consider a path from node 1 to 10 as
d1,10 (1, (4.5) (5.6), (6, (8,)) (9. 10)
(18)
bridge
bridge
the bridge path then given as
B1.10 ((, 10))
(19)
14

---

10
15
8
Figure 7: An example of computing bridge path. For example, from node 1 to 10, all the paths will cover the bridges (5, 6) and (9, 10). Hence, we have Br {(5, 6), (9, 10)}.
and 10 is computed as
R1,10 = S1.1 91,10 (R)
S10.10
R1.5
R5.6
R6.9
Rg
R10,10
- SS,
100
VS,
v Se,
v
S1.9
10
(20)
Particularly, we note that
10), (5, 6)} do would always visit the bridge The bridge path B0, o can not be written (0),
(5, 6) before visiting (9, 10).
The bridge path B1,10 not be written {66, 5), (9, 10)) d leave first cluster at node then enter second cluster at node 6. Therefore, we have (11 = (5)).
It is easy to observe that B10,1 {(10, 9), 5)). Based definitions, we can obtain Ri0,1 R1,10-
Specially, the matrix S is given as a correlation matrix, i.e., diag(S) = 1, then R1,10 R1,5 R556 R0,9 R9 R10,10. which can be seen continued produce on path, known as path product.
The matrix R shares some important properties introduced follows.
Lemma A.2. For any i E VE. Rs together with θ would satisfy the following equations
-[Bal1(0) St, = - -Ree) +iij Aij =
i-
(21a)
fi j. (+t)a(() # 0,
(21b)
-1Rel + S Ais 5O
ii j. (041-00.0)) = 0,
(21c)
Proof. Let's denote (())))())) the dual variables associated with the constraints e (xl=(1)=)) 0.
With Rg = 0-1 the KKT conditions include
stationarity:
Rs Sv. - Av. + +Iv = 0,
(22a)
primal feasibility:
2 Okla()) S0, Vz(1) # T())
(22b)
dual feasibility:
r [())-(() 2 VT(()) * T(j)
(22c)
complementary slackness:
3 O)l(()))) Tx(),(0) = 0, VT + x(1)
(22d)
We can eliminate the dual variables as follows:
For (O1 < 0, the complementary slackness leads to Tx(0)*() = 0, which implies
-(Ra)+(0, + +St.j Ar, 0.
For (Ox))(0)a()) = 0. we have T (())))))) 2 0, which implies - -(Rula(0, + St, - Aas S 0.
Hence, Lemma A.2 is obtained.
15

---

In particular, we note that (1 = St, holds for any E V according to (2la) As a result, the following property can be easily verified via the definitions.
Corollary A.3. Suppose that two nodes and j are in different clusters, i.e., vp (1) + vb (j), and the
bridge path is not empty and denoted by Bij with form of (44) then we have
Ri,., Ray
Raj
Vt 0, ..... 2T +1,
(23)
S Suc.
Proof. Corollary A.3 can be easily verified via definitions. t = 0 or t = 2T + Corollary A.3 holds Ra St,i and Ri Sy.- Let be number that satisfies 0 I 2T. For any Ugy Us+1) E Bij. we have
Rag.ay
Rui
Ris Sunu. = VSis
v Sun
v Sunnio Suyur
S 1 Suniis
Raaua
Ras
Rugr-
Vunw
VS
Vannu Petts
S Suzt 11227- Purrour
Rag
(24)
Clearly, Corollary A.3 is obtained.
The aforementioned properties of matrix R are important for establishing the proof of Theorem 3.4 and Theorem 3.3
A.2 Proof of Theorem 3.4
Proof. Let F OR, where O and have the forms of and (). respectively. The target is to show Fü 1, Vi and Fij = 0, vi + using the following results we have derived.
(a) For any i Vx. have
O [+
+ 1. Sz O (())(()))))))
E(Em)EB S,,mm Tim
if =j= j.
+
(25)
i j.
and Ray = (1))) according to the definitions.
(b) For any nodes and in different clusters, suppose that Biy is not empty and edge (u, v) is bridge in B.j. then S.yS... = R4u Rug and S. = R,,, R.g hold according to Corollary A.3
(c) The matrix Fr, which is defined as F. = R2 Gr is equal to the identity matrix, i.c. F. I. In other words, we have
(Alea =1, 1, and (P)))()) O, Vi and i, Vk-
(26)
Based on (a), (b), and (c), we first show that Fs holds for any f1,...,ph, p). We suppose that node i in kette cluster, i.e., = tb (t). Then, by using Vm & (Vr (N) Oim = from definitions in (7). have
Fu = Σ
Rim Om =
Σ
Rum Om + R +
Σ
R O4m
m=1
mE K0
(K.))B
(0)
Σ
(0l + RK ( (x)+
mE
T2
Tim (-5--10 Tim
+ Si
im
Sa
Sa.
T2
Sam T22
(Lm)ES
(i.m)B
Σ
O + [Ral @kl
[Palea(0)m20 =.
(27)
16

---

Therefore, Ta holds for any E PP).
Let proceed with the proof that Fij Em=1 Rim Ojj 0 for any i * Without loss of generality, we suppose that node is in the k-th cluster, i.e., h: == (j). According to the definitions, have
Vm & (VEUN (i)) Omy = O.
(28)
Here, there would be three cases to analyze.
The first case happens where i Vk. which means that nodes and j are in the same cluster. In this case, we have
P
Fij
Σ m=1
Ram mij
Σ
Rim mij + Riy +
Σ
Rim mij
mEVaMA M5)
(m.j)BB
Re (mm))))mmm))))aammm.m
meV
T
(
Σ
Σ (mjj)B
Rim Oij
Smm S33
T
(mJ)eB
mj
(c) Felra()) + Ris
1
Σ
mj
Σ (m.J)EB
Ram Omi
Sjs
S Smm Sis
T
(m.J)
(b.c)
0 Rij
Σ
T2
m.j
RijRij Sjj
Tmj
Sjj (m.j)E
Smm Sjz
T
S Sjn Sjz
T ma
(m.j)
= 0.
(29)
The second case happens where & Vk. but there is no path between nodes and In such a case, for any node m in Vy J. (j), there is also no path from nodes to m. Consequently, based on the definition of R, we have
Vm E VA UN(G) R.m = 0.
(30)
As result, we obtain
FúJ
Σ 1==1
Ram Omj
Σ
Rim Omj = 0.
(31)
MEVILNN)
The third case happens where i é V4 and there exists a path from nodes i to j. Then we have (i) # vi (f) and the bridge path Bij denoted
Bi {w1, , (us, u4) (u2r-1, u2r)
(32)
where (u,) & (i) and & (uer) (f) k. In particular, emphasize that ugr V/ as (u2r-1, u2r) is a bridge to enter the cluster that contains Therefore, we have
P
hh)
R REATT
P Rugr Omj O,
Rim Omj
(33)
Sugr 212
m=
m=1
where the last equality is from 29 by using u2r E VR. In conclusion, via enumerating all the possible cases, we have F I, which means that Re =I I.
Remark A.4. Theorem 3.4 generally holds for any T. S and (Ô: E s.ykk -k=1- It does not utilize
MTP, properties from O. Hence, Theorem 3.4 itself may extended graphical models beyond MTP2 assumption.
17

---

B Proof of Main Results
To prove Theorem 3.3 require some preliminaries on MTP₂ properties.
Lemma B.I.
Theorem 11 real symmetric mairix X E S", which admits the form:
Vi j: Xij 0,
(34)
is a non-singular M-matrix if and only if one of the following statements holds: (1) Every eigenvalue of X is positive; (2) X inverse-positive That is, X -- exists and x-13 0.
direct outcome of Lemma B.1 is that the matrix R is a non-negative matrix.
Corollary B.2. The mulrix R given by 8 satisfies R 2 0.
Corollary B.2 holds according the fact that
vk E Kl. E
R 2 0.
(35)
Now, we begin proving Theorem 3.3
Proof of Theorem
3.3 First of all, we show that the matrix O given as
e ()) +
i = j E VR= if * and i.j E Vk.
o (
O.,
(36)
(S, - T3)
if (f, B, otherwise.
primal feasible and positive definite. Based on Assumption 2.1 which states that S. < v S
we know that T2j S (Si, A.,)² S St < Sz. Hence, for any i * j. we have O: S due to (S T?) s 0 and the primal feasibility of each @ Additionally, since we have R > according to Corollary|B.. we can conclude that θ is inverse-positee
according to Theorem 3.4 which implies that it is mositive-defiiite matrix. In other words, @ is both primal feasible and inverse-positive, and we can express this as θ > 0 based on Lemma B.1
In second step, we show that θ = -R-1 satisfies the optimality condition listed as follows, i.e.,
- [ -[0-1),, + Sij = 0
if
(37a)
if if + 0, if = - 0. For 37a we let h: = v (i). Then, according to (1 in Lemma A.2 we have
o } + Sij - -As
=0 -A==0
(37b)
|0-1] "|s Si -A 00
(37c)
- - 10-"1,; + Sij - R (x),.()) + Sij
(38)
For 37b as Oi * 0, the pair (i.j) either satisfies (i,j) E B or 4b (1) = 1b (j)
If (1,j) E B.
then by using Raj Tej. we have
e '1y + Sij Auj -
Ríj Tτ 0.
(39)
If nodes i and j are in the same cluster, then according (1bb in Lemma A.2 we have O.j
( (+)))())))) * 0 such that
-10- + S.j - Aoj = -IR + Sij - Ag 0.
(40)
For 37c,. the pair (i.j) either satisfies ub (i) * (j) or /b (1) = vb (j). For first case, we have
Si Aij I 0 as (i,j) & supp (T). Consequently, we obtain
- 11, + Sij Ai S Rij S 0,
(41)
as R 2 0. For second case where sb (i) vb (j), we have O.₃ = (0x) = 0. Then, applying
21c in Lemma A.2 we have
o
"o S.j A,j = - RRx + Sij Ag S
(42)
Hence, together with the fact that t is primal feasible and positive definite, we conclude that θ satisfies the KKT conditions convex Problem (T). Therefore, θ is the optimal solution of Problem (I and Theorem 3.3 is obtained.
18

---

C Extensions and Discussions
C.1 Extensions traditional graphical lasso
The traditional graphical lasso is formulated as
minimize - log det (0) (0)+ + Eizg Ag 10:yl 9-0
(43)
The KKT conditions of graphical lasso include:
- RR + Stj =
if if * j. 02 + 0, if = j. 05 = 0.
(44a)
Ris + Sis Auj sign (0t)
(44b)
Aig S - -R33 + Sij- S Aij
(44c)
In the realm of graphical lasso, the thresholded sample covariance matrix usually defined as
if | S Aij. To = Sig - As sign (S.z) if Sipl > Aij.
0
(45)
Our main results Theorem 3.3 can be generalized to graphical lasso with some modifications.
Theorem C.1. Suppose that P is the bridge- blockk decomposition of the thresholded graph, O is the optimal solution of k-th sub-problen, i.e.,
θ arg min log det (0k) + (Ox, Sv.) + +E[+ [Av 10xll 0
(46)
inj
The decomposed form, which is given as
+ ()) +
ifi -i=je Vk: if j and i.j E Vk.
+ (xlm(f,
(47)
(SnSS - T2j)
if(i,) E B, otherwise,
T S. Ci.meEB SuSSmm -1 and C: = if (i, m) & B, achieves optimality if C > and the following inequality is satisfied:
in which C =
Raj + Siyl S Aij.
Vi Vi j and i,j belong to different clusters.
(48)
important to note that the decomposed form 47) aligns with the structure of with the key distinction lying in the computation of the thresholded matrix T. The proofs of Theorem C.I can be conducted similar to Section B
The condition (48) can be guaranteed under some circumstances. For example, we can set large values for Ai whenever * and i, are different clusters. In this case, the thresholded graph is sparse and the condition (48) is satisfied. Nevertheless, in the majority of practical situations, verifying this condition poses a substantial challenge. Conversely, these conditions are inherently satisfied when MTP₂ constraints are applied.
C.2 Extensions to dense graphs
The proposed method operates under an implicit assumption of sparsity, which is a common charac- teristic of high-dimensional applications where the focus is on identifying only the most significant relationships among variables. However, is possible encounter scenarios where the thresholded graph is dense, making it difficult to identify any bridges for conducting tridge-block decomposition. In such cases, proposed framework can still be useful, and this section will primarily focus on exploring its potential applications under these circums We will discuss various strategies and techniques for adapting the proposed method to deal with dense graphs, enabling broader applicability of the proposed bridge-blod block decomposition approach.
19

---

One of typical approaches to make the proposed method applicable is to tune the values of regulariza- tion coefficients. As the thresholded matrix is computed as
T max (0, S A),
(49)
we can always increase the value of Ai to make the thresholded graph supp(T') sufficiently sparse. When tuning the hyper-parameters is not recommended, a more general strategy use the solution of proposed framework as an initial point of numerical algorithms.
Warm-start algorithm: One of the main possible strengths of the proposed decomposed-fore solution is that it can be treated as an initial point (warm-start) for the numerical algorithms specialized for solving Problem To apply the decomposition approach, we first need derive suitable vertex partition of the graph based on the given data matrix S and regularization matrix A. There are several possible vertex partitioning techniques that can be used, including:
k-edge-connected decomposition: This method divides a connected graph into its k-edge-connected sub- sub-graahs, which can be useful for identifying bridges bottlenecks in the graph structure.
k- vertex--oomeeted decomposition: This method decomposes a graph into its sub-graph -rappss which can be useful for identifying clusters or communities of vertices that are tightly connected to each other.
k-core decomposition: This method identifies the maximal subgraphs of the graph in which each vertex has degree of at least k. which can be useful for identifying the most densely connected regions of the graph.
Community detection: This method identifies groups of vertices that are densely connected within the group and sparsely connected vertices outside the group, which can be useful for identifying clusters or communities of vertices with similar properties or behavior.
Once we obtain a vertex partition P (Vi...., VR). ). where K refers the number of clusters, can compute the initial point using the following form:
+ ())) + +
if ffi=jEl Vh: if i and i,j E Vk, if (i.5)= E Ec. otherwise,
+ (x)x())x(())
O
(50)
(S - T3)
in which & refers to set of external edges among clusters, i.e..
E
E supp (T), + # (1),
(51)
G is computed
G Ssi (i,m)E&.
im
(52)
S
T
and O is computed as the optimal solution of k-th sub-problem:
θ k arg
min EMMM
log det (0x) (Ox, Sv. Av
(53)
Remark C.2. The initial point is the optimal solution if the vertex-partition P is bridge-block decomposition the optimal graph.
Remark C.3. When the vertex partition P is given P ((1) ( (p)). i.e. Vi = [k) for ... the initial point (50) reduces to a closed-form expression:
1
The 1 (1+ ((.))))) Sa8mm T Tis (S 10--1-7 T7))
if Fi=j. j. if (i.j) E supp (T), otherwise.
S..
O.j
(54)
In practice utilizing the warm-start strategy can also enhance the performance of certain numerical algorithms.
20

---

-50
Average Test Log likelihood
-75
-100
MTP2
-150
0.25
0.50
0.75
1.00
Figure 8: Average lig-liielihood in out-of-sample data for different precision matrix estimation models as function of the sparsity promoting hyperparameter A.
Figure 9: Estimated Graphical Lasso of Crop data set at their highest average likelihood. Red edges represent positive conditional correlations, while blue edges represent negative ones.
D Additional Comments on Real-World Experiments
To show that MTP: assumption approximately holds for Crop dataset, we selected random subsets from the Crop dataset, For each subset, we computed the Graphical Lasso and the MTP, graphical model for different values of using the first 10 observations. The remaining 36 observations were used to calculate the out-of-sample log- ligliik,, which was then averaged across all datasets. This process allows us to evaluate how well these models generalize to unseen data.
As depicted in Figure of the attached PDF, the MTP₂ graphical model outperforms the Graphical Lasso, providing higher test log-likelihood. We present one instance of the estimated graphical lasso model in Figure 9 It reveals that most conditional correlations are positive (red edges, 90%), with few being negative (blue edges, 10%). This pattern implies strong positive dependence in the Crop data, aligning with the characteristics of MTP.
These results not unexpected, given that the Crop dataset comprises multiple clusters. Within the same cluster, we expect data points to exhibit greater similarity compared to those in different clusters. This situation signifies form of positive dependence, thereby justifying the plausibility of the MTP, assumption in this problem.
21


Reviews:
[{"review_id": "BKZvAE7EMM", "reviewer": "Reviewer_tnT4", "summary": "The paper focuses on the problem of learning large-scale Gaussian graphical models (GGMs) that are multivariate totally positive of order two (MTP2). The high-dimensional, sparse MTP2 GGMs are not easily manageable due to their size and complexity. The authors propose a novel approach, introducing the concept of a \"bridge\", to optimize the entire problem into several smaller, more manageable sub-problems and a set of closed-form solutions. The approach is based on the bridge-block decomposition of the thresholded sample covariance graph, which leads to reductions in computational complexity and improvements in existing algorithms. ", "questions": "Is there any way to characterize the effect of the proposed approach in terms of some global graphical properties, for example, the edge expansion? Asking this because bridge is related. ", "limitations": "See Above.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The proposed bridge-block decomposition framework on Gaussian graphical models seems novel. The problem is motivated nicely, and according to the authors, the framework could significantly reduce computational and memory cost. \n\nThe proposed method seems to subsume various network structures, including the BA graph and the SBM, which are common models used in network analysis.\n\nExperimental results are provided, and the computational results look promising. ", "weaknesses": "As the authors mentioned in the paper, the proposed method might not generalize to dense cases. Still I feel in many settings like BA and SBM, sparsity is a reasonable assumption.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## Answer to Questions \n\n> Is there any way to characterize the effect of the proposed approach in terms of some global graphical properties, for example, the edge expansion? Asking this because bridge is related.\n\n__Reply:__ We appreciate the reviewer's interesting and insightful question. Indeed, certain global graphical properties could potentially serve as indicators of a graph's connectivity strength, thereby providing a measure of our method's effectiveness. However, the computation of edge expansion presents a significant challenge, especially for high-dimensional graphs.\n\nAs an alternative, we opted to utilize another global graphical property known as __algebraic connectivity__. This metric, the second smallest eigenvalue of a graph's Laplacian matrix, offers a depiction of the graph's overall connectivity. It can be employed to characterize the impact of our proposed solution.\n\nTo demonstrate this, we executed tests on two random SBM graphs with distinct algebraic connectivities. We evaluate the how many times our proposed method can accelerate the convergence of existing algorithms. The corresponding results are displayed in the table below:\n\n|Algebraic connectivity |BCD|PGD|FPN|PQN-LBFGS|\n|:-:|:-:|:-:|:-:|:-:| \n| 4e-5 |$2064$ | $314$  |$213$  |$57$   |\n| 4e-4 |$26$  | $9$  |$13$  |    $5$  |  \n    \nIt is obvious that our proposed method is much more effective on graph with low algebraic connectivity.\n\nDetails regarding convergence (refer to Figure 3) and more in-depth analyses (refer to Figure 1) are available in the attached PDF. In detail, we generate numerous random graphs and for each, we examine the relationship between its algebraic connectivity and a theoretical number of how many time we can accelerate a BCD method. In Figure 1, each point plotted represents a specific graph, with the x-axis indicating the algebraic connectivity and the y-axis representing the acceleration factor.\n\n__Conclusion__: The results indicate a general trend: __as the algebraic connectivity decreases, the effectiveness of our proposed method enhances.__\n\nRemarks: The theoretical number of acceleration is computed as follows. Let's assume that a problem with dimension $p$ can be addressed by BCD algorithm in $c\\cdot p^4$ seconds, where $c$ is a constant. By implementing the bridge-block decomposition, this cost is reduced to $\\sum_kc\\cdot p_k^4$, with $p_k$ denoting the size of the $k$-th cluster. Therefore, we define the ratio $p^4\\big /\\sum_k p_k^4$ as the theoretical number by which our method can speed up the BCD algorithm.\n   \n## Reply to Comments\n\n> The proposed method seems to subsume various network structures, including the BA graph and the SBM, which are common models used in network analysis. As the authors mentioned in the paper, the proposed method might not generalize to dense cases. Still I feel in many settings like BA and SBM, sparsity is a reasonable assumption.\n\n__Reply:__  We appreciate the reviewer's comment and acknowledge that there are numerous applications beyond the scope of this paper where the graph is dense.\n\nIn response to this, we have added Section C in our Appendix, which discusses various alternative strategies for integrating our proposed method into the learning of dense MTP2 graphical models. For instance, though our decomposition formula only yields exact solutions for sparse graphs via bridge-block decomposition, it could be regarded as __an approximate solution__ for dense graphs. Additionally, the explicitly decomposed form could serve as an __effective warm start__ for other numerical algorithms."}, {"review_id": "xJisCZd2GZ", "reviewer": "Reviewer_vWNe", "summary": "This paper studies the problem of learning Gaussian Graphical Models (GGMs) satisfying a certain positive associativity condition among the variables, namely that the precision matrix has nonnegative off-diagonal elements. This condition is known as being \"multivariate totally positive of order two\", or MTP$_2$, and has applications in ML (where it corresponds to attractive Markov random fields), finance, and more.\n\nMTP$_2$ GGMs come with the benefit that the traditional optimization procedure used to estimate the precision matrix of a GGM, namely the graphical Lasso, takes on a particularly simple and smooth form. Prior work had shown various polynomial-time convergence guarantees for the graphical Lasso under the MTP$_2$ assumption, but these are still not very suitable for practical applications (scaling with the dimension $p$ as $O(p^3)$ or $O(p^4)$). Other prior work had shown a closed-form solution for the graphical Lasso under the assumption that the \"thresholded sample covariance graph\" (an object defined in terms of the sample covariance and the regularization parameters used in the graphical Lasso), or thresholded graph for short, is acyclic.\n\nThe main contribution of this paper is essentially to generalize the latter closed-form result in terms of the \"bridge-block decomposition\" of the thresholded graph. The bridge-block decomposition of a graph is essentially a partition of the graph into components connected only by \"bridge\" edges. Formally, a bridge edge is an edge such that deleting it increases the number of connected components in a graph; it is effectively \"the only edge\" bridging two different components (see Fig 1). The paper's main theorem (Thm 3.3) essentially says the following: to solve the graphical Lasso for an MTP$_2$ GGM, compute the bridge-block decomposition of the thresholded graph, run the graphical Lasso for each component separately, and stitch them together using a simple closed-form formula. Moreover, this theorem also readily recovers as a special case the prior closed-form result for acyclic thresholded graphs. This is because in an acyclic graph, every edge is a bridge, and the bridge-block decomposition is particularly simple.\n\nThus the main result amounts to a divide-and-conquer recipe for learning MTP$_2$ GGMs, and the authors show various numerical experiments suggesting the practical superiority of this method over all prior methods (which operate on the entire graph). A key benefit is that the subproblems may be solved using any graphical Lasso implementation whatsoever, and potentially in parallel.\n\n\n\n", "questions": "I think the biggest question I have is about the conceptual role of the bridge-block decomposition, as explained above. A few other questions:\n- One of the key parts of the eventual proof of Thm 3.3 is Lemma A.2, which is simply described as following from the KKT conditions. Presumably the authors mean the KKT conditions associated with Problem (5)? Since it plays a relatively important conceptual role, I feel this could definitely use more elaboration.\n- What exactly is the time required to compute the bridge-block decomposition? This is described as negligible in Sec 3.1, but is it e.g. $O(p^2)$? Similarly, what about the quantities used in the closed-form formula? Currently the only discussion of the overall asymptotic running time appears in the first bullet of the \"Proposed Framework\" list on page 5.\n- In the experimental section, it would be helpful to say more about it is natural to synthesize/define $\\Theta$ and $\\Lambda$ in the specific way described in the first two paragraphs of Section 4.1, especially for readers who are unfamiliar with the prior work mentioned in those paragraphs.\n\nA couple other nits regarding the presentation:\n- It is not immediately obvious how Corollary 3.5 follows from Thm 3.3, and this could use a line of explanation.\n- The notation $\\mathbb{S}^p$ is not formally defined in the paper, and MTP$_2$'s expansion is only mentioned in the abstract.", "limitations": "The paper could use a couple additional lines in the final section about the technical limitations of this work and what the major next steps could be. I am not aware of any significant potential negative societal impact of this theoretical work.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "Disclaimer: I am not very familiar with the literature in this area, and my review should be taken as that of a relative outsider.\n\nThe paper's main result is both an interesting structural result about MTP$_2$ GGMs as well as a genuinely practical algorithmic advance in learning such models. From a conceptual point of view, the idea of leveraging the bridge-block decomposition seems novel. The overall result seems like a useful and nice contribution to the literature on this problem, and to the extent that one considers MTP$_2$ GGMs significant, one should consider this result significant as well.\n\nThe paper is largely clear and easy to follow (modulo some occasionally confusing bits; see the Questions section). It does a good job of setting up the main problem as well as the necessary context. I did not manage to verify the proofs in detail, but they seemed fairly clean, relying on an analysis of the KKT conditions of the graphical Lasso as well as some clever algebraic manipulation.", "weaknesses": "I think the main things to really evaluate about this paper are its novelty and significance. As an outsider to this area, I find this hard to accurately gauge, but I think the paper scores well on these fronts.\n\nI do think the paper could benefit from a better conceptual overview of the main proof and the role of the bridge-block decomposition. The context and benefits are discussed adequately, but the key ideas in the proof do not come through very well, and the main proof seemed slightly magical to me. Why would one have expected the bridge-block decomposition to help? Was its role surprising?\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## Answer to Questions Part 1 \n> One of the key parts of the eventual proof of Thm 3.3 is Lemma A.2, which is simply described as following from the KKT conditions. Presumably the authors mean the KKT conditions associated with Problem (5)? \n\n__Reply:__ Thank you for your insightful comment. Here, we elaborate on their relationship.\n\nLet's denote $\\Gamma_{ij}$ as the dual variables associated with the constraints $\\Theta_{ij}\\leq 0$.  With $\\mathbf{R}=\\boldsymbol{\\Theta}^{-1}$, the KKT conditions include (1) $-\\mathbf{R}+\\mathbf{S}-\\boldsymbol{\\Lambda}+\\boldsymbol{\\Gamma}=\\mathbf{0}$, (2) $\\Theta_{ij}\\leq0,\\forall i\\neq j,$ (3) $\\Gamma_{ij}\\geq0,\\forall i\\neq j,$ (4) $\\Theta_{ij}\\cdot\\Gamma_{ij}=0,\\forall i\\neq j.$ \n\nWe can eliminate the dual variables as follows:\n  \n1. For $\\Theta_{ij}<0$, the complementary slackness leads to $\\Gamma_{ij}=0$, which implies $-R_{ij}+S_{ij}-\\Lambda_{ij}=0$;\n\n2. When $\\Theta_{ij}=0$, we have $\\Gamma_{ij}=R_{ij}-S_{ij}+\\Lambda_{ij}\\geq 0$, which indicates that $-R_{ij}+S_{ij}-\\Lambda_{ij}\\leq 0$.\n\nFollowing these steps, we arrive at the optimal conditions for original problem. As the sub-problems have the same form with the original problem, Lemma A.2 holds accordingly. We will incorporate these details in the revised version.\n\n## Answer to Questions Part 2\n\n> What exactly is the time required to compute the bridge-block decomposition? Is it e.g. O(p^2)? What about the quantities used in the closed-form formula?\n\n__Reply:__ We sincerely appreciate the reviewer's question. Please find our detailed answer in the Part 3 of our global response.\n\n## Answer to Questions Part 3\n\n> In the experimental section, it would be helpful to say more about it is natural to synthesize/define $\\boldsymbol{\\Theta}$ and $\\boldsymbol{\\Lambda}$.\n\n__Reply:__ We appreciate your valuable recommendation and will revise our paper accordingly. The methods we use to synthesize $\\boldsymbol{\\Theta}$ and $\\boldsymbol{\\Lambda}$ are guided by [1] and can be explained as follows.\n\nWe synthesize $\\boldsymbol{\\Theta}$ as $\\boldsymbol{\\Theta}=1.05\\cdot\\lambda_{\\max}(\\mathbf{A})\\cdot\\mathbf{I}-\\mathbf{A}$, where $\\mathbf{A}$ is the adjacency matrix of the underlying graph. This ensures that $\\boldsymbol{\\Theta}$ is a positive definite matrix with off-diagonal elements being negative, making $\\boldsymbol{\\Theta}$ a randomly generated M-matrix. Then, we normalize $\\boldsymbol{\\Theta}^{-1}$, thereby deriving a randomly generated correlation matrix.\n\nFollowing this, we set $\\Lambda_{ij}=\\chi\\big/(\\epsilon+\\Theta_{ij}^{(0)})$. Here, a high penalty is expected on $\\Theta_{ij}$ if its initial estimate $\\Theta_{ij}^{(0)}$ is small. By applying this way, we can recover the underlying structure by selecting an appropriate value of $\\chi$.\n\n[1] Martin Slawski and Matthias Hein. Estimation of positive definite M-matrices and structure learning for attractive Gaussian markov random fields. Linear Algebra and its Applications, 473:145–179, 2015.\n\n## Answer to Questions Part 4\n> It is not immediately obvious how Corollary 3.5 follows from Thm 3.3.\n\n__Reply:__ We appreciate your valuable suggestion. Here's a clarification on why a bridge $(i,j)$ in the thresholded graph will persist as a bridge in the optimal graph.\n\nAn edge $(i,j)$ is a bridge if and only if there exists exactly one unique path connecting nodes $i$ to $j$, denoted as $d_{ij}=\\{(i,j)\\}$. According to the definition of a bridge, the removal of it would lead to an increment of graph's components. The presence of any additional paths would contradict this definition.\n\nAs the optimal graph is a subset of the thresholded graph and $\\Theta_{ij}\\neq 0$, we find that $d_{ij}=\\{(i,j)\\}$ continues to be the unique path from nodes $i$ to $j$, which indicates that $(i,j)$ remains a bridge.\n\n## Answer to Questions Part 5\n\n> The notation $\\mathbb{S}^p$ is not formally defined in the paper, and MTP$_2$'s expansion is only mentioned in the abstract.\n\n__Reply:__ We appreciate your keen observation that calls for clarity. The symbol $\\mathbb{S}^p$ represents the set of symmetric matrices with a dimension of $p\\times p$. We will incorporate the full form of MTP2 into the body of the paper.\n\n## Reply to Comments\n\n> I do think the paper could benefit from a better conceptual overview of the main proof and the role of the bridge-block decomposition. Why would one have expected the bridge-block decomposition to help? Was its role surprising?\n\n__Reply:__ We sincerely value your comprehensive review of our paper and your engagement with the core concepts underlying our proof.\n\nOur research is influenced by the study of (Fattahi et al., 2019), which indicates that a closed-form solution exists for large-scale graphical lasso problems. Their intriguing findings hold potential for large-scale data sets, yet their theorem imposes some hard-to-validate conditions and requires the thresholded graph to be acyclic.\n\nIn the realm of the MTP2 graph learning problem that interests us, we aimed to conduct related research and unearth similar findings. Our investigations led us to discover that (1) the existence of a closed-form solution hinges not on the acyclicity of the thresholded graph, but on whether an edge is a bridge; (2) Given that bridges provide closed-form solutions, if our aim is to decompose the problem, we should form the sub-problem excluding the bridges, giving birth to our bridge-block decomposition strategy; (3) By referencing the proofs in (Fattahi et al., 2019), we eventually prove that all conditions that are typically difficult to verify are naturally satisfied in MTP2 graphical models.\n\nAs our approach does not impose additional conditions, and it does not require the thresholded graph to be acyclic, our method has a wider range of applicability.\n \n-- Salar Fattahi and Somayeh Sojoudi. Graphical lasso and thresholding: Equivalence and closed-form solutions. Journal of Machine Learning Research, 2019."}, {"review_id": "4UyNTm3qwr", "reviewer": "Reviewer_bj5B", "summary": "The paper studies the problem of estimating the precision matrix, which is the inverse of the correlation matrix, of a given Gaussian random vector $y$. The precision matrix $\\Theta$ is assumed to satisfy a technical condition called MTP2 which states that $\\Theta$ is symmetric and $\\Theta_{i,j} \\le 0$. This seems to be a well motivated assumption from various applications. The contribution of this paper is a technique for estimating $\\Theta$ as follows: given a predicted sparsity pattern on $\\Theta$ in the form of a graph $G$ , the natural optimization problem for estimating $\\Theta$ can be solved by first solving the optimization problem on smaller 'blocks' and combining them across 'bridges'. They are defined as follows. Bridges are single edge cuts in  $G$ and the resulting connected vertices are called blocks.", "questions": "What happens to the SBM experiments if the edge probabilities between the different blocks are increased? The quality of the method should deteriorate as the probabilities increase, since the block structure increasingly deteriorates. It would be interesting to see how much the proposed method can tolerate.", "limitations": "No ethical concerns.", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "The paper shows that given a block bridge decomposition, the optimization problem of estimating $\\Theta$ can be efficiently solved by first solving the problem on the individual blocks and then combining the solution across bridges. The paper gives an explicit formula for doing so. In the case where $G$ is sparse, this can represent significant computational savings over estimating the entire precision matrix at once. Furthermore, since the work provides a structural theorem, any optimization algorithm can be used in conjunction with their observation.", "weaknesses": "I am not familiar with the literature but it seems like a big assumption to know the threshold graph explicitly. What happens if this graph is unknown? It seems to be more natural that the graph is unknown and one must estimate it.\n\nAn intermediate setting which also seems interesting is in the case where we know a noisy approximation to the block bridge structure. How do the proposed methods perform under such noisy information? Are the derived formulas robust?\n\nWhat is the motivation for even assuming the block bridge structure? I can see social networks being one motivation but it would be more convincing if there were experiments on (real) social networks.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## Reply to Comments Part 1\n\n> I am not familiar with the literature but it seems like a big assumption to know the threshold graph explicitly. What happens if this graph is unknown? It seems to be more natural that the graph is unknown and one must estimate it.\n \n__Reply:__ We appreciate your feedback and would like to provide further clarification. Our problem is defined such that we already have access to the sample covariance matrix, denoted by $\\mathbf{S}$, and the regularization matrix, represented as $\\boldsymbol{\\Lambda}$. \n\n__Given these matrices__, we can __precisely__ compute the thresholded matrix, $\\mathbf{T}$, using the formula $\\mathbf{T}=\\max(\\mathbf{0},\\mathbf{S}-\\boldsymbol{\\Lambda})$. Hence, the knowledge of the thresholded graph isn't an assumption we're imposing, but rather a natural outcome resulting from the specific optimization problem at hand.\n\n## Reply to Comments Part 2\n \n> An intermediate setting which also seems interesting is in the case where we know a noisy approximation to the block bridge structure. How do the proposed methods perform under such noisy information? Are the derived formulas robust? \n\n__Reply:__ Thank you for posing this question. In practical applications, samples inevitably carry noise, making it challenging to obtain an accurate estimate of the sample covariance matrix.\n\nHowever, regardless of the sample size or the degree of sample noise, we can still apply our proposed method to the noisy sample covariance matrix. This is because the step of computing the bridge-block decomposition from the sample covariance matrix is deterministic, and __we can always obtain an optimal and exact solution__ via our proposed framework. As a result, from an optimization perspective, our method is robust.\n\n## Reply to Comments Part 3\n\n> What is the motivation for even assuming the block bridge structure? I can see social networks being one motivation but it would be more convincing if there were experiments on (real) social networks.  \n\n__Reply:__ Thank you for your insightful question. Our motivation in assuming the block bridge structure stems from our focus on sparse graph learning. In sparse graphs, bridge is a commmon feature, and here's why.\n\nOne motivating for learning sparse graphs is to enhance interpretability by preserving only the most significant relationships among variables. Each node connects with its most important neighbors, reducing the probabilities for cycle formations. Furthermore, the connectivity in sparse graphs is rather weak due to the limited number of edges, making them easy to separate. As a result, numerous edges in sparse graphs become bridges - these are edges whose removal would create additional components.\n\n\nIn the context of social networks, the stochastic block model is a common tool for representing such networks. Notably, social networks often exhibit strong intra-group connections and weaker inter-group ties. This structure may potentially harbor a considerable number of bridges.\n\nWhile we didn't have the opportunity to apply our methods to social network learning due to time limit, we believe there is potential merit that our methods could contribute to the more efficient learning of these networks.  \n\n## Answer to Questions\n\n> What happens to the SBM experiments if the edge probabilities between the different blocks are increased? The quality of the method should deteriorate as the probabilities increase, since the block structure increasingly deteriorates. It would be interesting to see how much the proposed method can tolerate.  \n\n__Reply:__ We appreciate your highlighting of this practical issue. We concur that the effectiveness of the method might decline as the graph grows denser.\n\nIn response to your queries, we conducted additional experiments to examine the extent to which we can speed up the convergence of the BCD method for various values of $p'$. Here, $p'$ represents the probability of forming an edge $(i,j)$, where $i$ and $j$ are any two distinct nodes in neighboring communities.\n\nWe utilized a SBM graph with $1500$ nodes for this study. The results are illustrated in Figure 2 of the attached PDF and are also encapsulated in the table below:\n \n|Edge formation probability|0.05|0.10|0.15|0.20|0.25|0.30|0.35|0.40|0.45|  \n|:-:|:-:|:-:|:-:|:-:|:-:|:-:|:-:|:-:|:-:| \n| Ratio of Improvement |$179.7$|$150.1$|$80.1$|$60.1$|  $31.1$|$16.6$|$3.6$|$3.5$|$3.2$| \n\nHere, Ratio of Improvement refers to how many times we can accelerate the convergence of BCD methods. As expected, a rise in $p'$ increases the chance of multiple edges linking the blocks, which can potentially hamper the efficiency of our approach.\n\nAs we address in Section C of our Appendix, while our method is primarily intended for sparsely connected graphs, it retains its usefulness even for dense graphs. Our method can serve as an approximate solution or provide a warm start for other numerical algorithms, thereby boosting their computational efficiency."}, {"review_id": "otfWB8ISVv", "reviewer": "Reviewer_ghhX", "summary": "This paper studies a graphical lasso problem where the precision matrix is restricted to be symmetric M-matrix and the associated GMRF graph has a special structure. Specifically, the authors consider the situation when the graph allows bridge-block decomposition so that vertices can be partitioned into k parts by cutting k-1 \"bridges\".  Under this situation, the authors show that the original problem can be decomposed into k subproblems, illustrated in theorems in section 3.2. Since the existing algorithm's complexity quickly increases as the dimensionality increases, decomposition into subproblems offers large computational benefits. Authors compare the performance of four different algorithms, using synthetic and real data, with and without leveraging the bridge-block decomposition and show that bridge-block decomposition provides huge computational benefits if the given graph has many edges that are “bridges”.", "questions": "- Does the main result, decomposition into smaller graphical lasso subproblems under block-bridge decomposition, holds without MTP2 constraint? \n- If not, please clarify the the role of the MTP2 constraint in the main results; does MTP2 constraint is necessary for the assumption 2.1. to be hold?\n- What happens if graph has many connected components (without bridges) in MTP2 constrained graphical lasso problem? Does similar decomposition holds?\n- Section 4.1: Sample size $n$ is 10 times the dimension $p$? This setting completely violates the main motivation of graphical lasso where $\\ell_1$ penalty has been proposed since MLE does not exist when $n<p$. \n- Section 4.3: Shouldn't be data $y_i$ is 24000-dimensional with $i=1,\\dots,46$ observations?\n", "limitations": "The main results are meaningful but simulation results and real data settings are less convincing. No potential negative societal impact.", "rating": 6, "confidence": 4, "soundness": 2, "presentation": 2, "contribution": 2, "strengths": "Authors provide a closed form of the elements of precision matrix corresponding to bridges, and make a connection with existing literature when underlying GMRF graph is an acyclic graph. This method is useful when one wants to solve MTP2 constrained graphical lasso problem with a penalty parameter that allows block-bridge decomposition. This work is resembles the existing literature on graphical lasso (without MTP2 constraint), such as Witten et al (2011); Mazumder and Hastie (2012), where graphical lasso problem can be decomposed to the smaller subproblems when graph has many connected components, and adds a contribution in the context of MTP2 constrained graphical lasso when graph has many bridges.\n\n\nWitten, D. M., Friedman, J. H., & Simon, N. (2011). New insights and faster computations for the graphical lasso. Journal of Computational and Graphical Statistics, 20(4), 892-900.\n\nMazumder, R., & Hastie, T. (2012). Exact covariance thresholding into connected components for large-scale graphical lasso. The Journal of Machine Learning Research, 13(1), 781-794.", "weaknesses": "Simulation settings and real data examples have some discrepancy from usual (graphical) lasso settings and its motivation. The main goal of graphical lasso is to discover the underlying conditional independency structure from the multivariate Gaussian data with various choice of penalty parameter (aka lasso path diagram). However authors first fix the graph in their settings such as preferential attachment graph and stochastic block model and choose penalty parameter according to the graph. It is understandable in the simulation study setting to show the computational benefits with block-bridge decomposition, but fixing the graph in the real data is not convincing, and also it is not clearly stated why MTP2 graphical lasso is appropriate to this crop image data problem (is it reasonable to assume that crop image data comes from high-dimensional multivariate Gaussian distribution? what is the interpretation of the resulting graph? why MTP2 constraint is necessary / plausible in this problem?)\n\nIn practice, graphical lasso is ran under various settings of penalty parameters and chosen appropriately such as cross validation. The existing graphical lasso decomposition (Witten et al (2011); Mazumder and Hastie (2012)) are useful since graphical lasso with high penalty parameter often leads to a graph with many connected components. In similar spirit, I suggest authors to illustrate how often the graph that allows block bridge decomposition appears by the choice of penalty parameter.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## Answer to Questions Part 1\n\n>Does the main result holds without MTP2 constraint? \n\n__Reply:__ We thank the reviewer for raising this point. The MTP2 constraints are essential prerequisites for our main results.\n\n## Answer to Questions Part 2 \n\n> If not, please clarify the the role of the MTP2 constraint in the main results. \n\n__Reply:__ We appreciate the reviewer for raising this insightful question. Please refer to part 1 of global response for our detailed answer.\n  \n> Does MTP2 constraint is necessary for the assumption 2.1. to be hold?  \n\n__Reply:__ The MTP2 constraint is not a prerequisite for assumption 2.1. Assumption 2.1 primarily ensures the existence and uniqueness of the optimal solution.\n\n## Answer to Questions Part 3 \n \n> What happens if graph has many connected components (without bridges) in MTP2 constrained graphical lasso problem? Does similar decomposition holds?  \n  \n__Reply:__ We appreciate your thought-provoking question. Yes, our theorem can also find the exact solution like the existing graphical lasso decomposition. Here's a simple comparison: \n\n|Graph Type|Description|\n|:-|:-|  \n|__Single-component__ graph with bridges|Our method works, the __existing method doesn't__.| \n|__Multi-component__ graph with bridges|Our method __works better__ than existing method.|  \n|__Multi-component__ graph without bridges|Our method has __equivalent__ effectiveness to the existing method.| \n \nIn theory, bridge-block decomposition refers to the components after removing all bridges. So, even in a graph with many components and no bridges, we can still find clusters. This makes bridge-block decomposition versatile, as it can deal with __both connected and disconnected__ sparse thresholded graphs.\n\n## Answer to Questions Part 4\n\n> Section 4.1: Sample size n is 10 times the dimension p? This setting completely violates the main motivation of graphical lasso. \n\n__Reply:__ Thank you for your thoughtful feedback. We understand the graphical lasso is often used when $p>n$. However, __the aim of our experiments__ is not merely to uphold the statistical rationale of the MTP2 graphical lasso, but to __demonstrate how our proposed theorem can significantly accelerate the learning of high-dimensional sparse MTP2 graphs from an optimization perspective__. Hence, our experimental design is crafted to ensure successful graph structure recovery across all experiments. \n \nIn response to your suggestions, we've added new experiments (see Figure 4 in the attached PDF). Here, we set the sample size to $n=0.1p$, where $p$ is the dimension. These experiments show that our method still effectively speeds up convergence due to the sparsity of the underlying structure.\n\n## Answer to Questions Part 5\n\n> Section 4.3: Shouldn't be data yi is 24000-dimensional with i=1,...,46 observations?\n\n__Reply:__ In the context of graphical models, every node corresponds to a data, while the edges represent the conditional dependencies between the data. As a result, our network comprises 24000 nodes. It is important to clarify that in this scenario, yi does not represent a feature, but rather a signal consisting of 46 observations. This yields a configuration where p=24000 and n=46.\n\n## Reply to Comments Part 1 \n\n> Fixing the graph is not convincing, and also why MTP2 graphical lasso is appropriate to this crop image data problem (what is the interpretation of the resulting graph? why MTP2 constraint is necessary / plausible in this problem?)\n\n__Reply:__ We appreciate your practical question. To alleviate your concerns, we have included more experiments in Part 4 of the global response to show that the MTP2 graphical Lasso is suitable for graph-based clustering with the CROP dataset.\n\nThe results show that the CROP data exhibit positive dependency suit for MTP2 assumption. To the estimated graph, edges represent positive conditional dependence between variables. In a clustering context, this often suggests that these interconnected nodes in the estimated graph tend to fall within the same cluster.\n\nGiven our observation that the data adheres to MTP2 properties, it would be advantageous to include MTP2 constraints __as prior information to improve the learning of graphical models__. The MTP2 structure also possesses excellent mathematical properties that confer __significant computational advantages__, enabling applying bridge-block decomposition to solve high-dimensional sparse graphs, which can not be solved via existing methods.\n\n## Reply to Comments Part 2\n\n>  The existing graphical lasso decomposition are useful since graphical lasso with high penalty parameter often leads to a graph with many connected components. I suggest authors to illustrate how often the graph that allows block bridge decomposition appears by the choice of penalty parameter.\n\n__Reply:__ Thank you for your useful comments on practical aspects. To respond, we've done more tests using the same steps as in our synthetic data experiments. We created a random SBM model with 2000 nodes. By adjusting $\\lambda$, we checked how the decomposition methods perform under different thresholded graphs. The results are shown in the table below:\n\n|$\\lambda$|0|0.01|0.03|0.1|0.15|0.18|0.3|0.5|1|\n|:-|:-:|:-:|:-:|:-:|:-:|:-:|:-:|:-:|:-:|\n|Number of clusters by __existing decomposition__|1|1|1|1|1|__15__|__1681__|__1974__|__2000__| \n|Number of clusters by __bridge-block decomposition__| 1|1|__9__|__402__|__406__|__411__|__1793__|__1978__|__2000__|\n\nWe found that:\n\n1. Our method becomes effective when $\\lambda\\geq 0.03$, while the existing decomposition method is applicable when $\\lambda\\geq 0.18$.\n2. Theoretically, our approach can handle both scenarios: when the graph has bridges or when it has multiple components. Conversely, the existing method can only cope with the latter case.\n\nThe results clearly illustrate that __our method exhibits a broader range of applicability__ compared to existing decomposition methods."}, {"review_id": "RoGzVqtJZ8", "reviewer": "Reviewer_pk8s", "summary": "The authors show that, in Gaussian MTP2 distributions, bridges in the graph structure have a closed form solution.  They use this observation to suggest practical solutions that can be applied whenever such models are being fit.", "questions": "- \"calculating the thresholded graph, bridges, and clusters, is negligible\" Sure, but I feel that a precise statement and citation is needed here.  In particular, you still have to compute the sample covariance matrix for this.\n\n- I actually had to read it a couple times to understand the exact procedure that was being proposed.  I suggest adding some clarifying details to help the reader.\n\n-  What is the primary limitation that makes this result only apply to MTP2 distributions?  Is it possible that something similar could be true more generally?  I didn't really get a sense of why MTP2 was required in the main text.\n\n-  Do you need Gaussian distributions here?  Can something be done for more general MTP2 distributions, e.g., those that can still be represented as pairwise graphical models?\n\nMisc. typos (only a few listed):\n- \"which severs as a common assumption\"\n- \"make possibilities of solving...\" (revise)\n- \"As thresholded graph plays\"\n- \"garnered considerable attentions\"\n- \"Bridge is one of the important concepts\"\n- \"closed-form solutions in literature heavily\"\n- \"seems inapplicable for dense graph\"", "limitations": "Some discussion of limitations, though it would have been nice to see a more clearly identified set of problems yet to be solved.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "A nice, more or less self-contained theoretical work that unifies and extends some existing lines of research into MTP2 Gaussian distributions.", "weaknesses": "Some room for improvement in terms of presentation and some minor typos.  It's harder for me to judge practical utility.  While it is true that MTP2 distributions have some practical utility (I've used them myself), it isn't clear that many real datasets meet that requirement (even approximately).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## Answer to Questions Part 1\n\n>\"calculating the thresholded graph, bridges, and clusters, is negligible\" Sure, but I feel that a precise statement and citation is needed here. In particular, you still have to compute the sample covariance matrix for this. \n\n__Reply:__ We are grateful for your valuable suggestion. We would like to elaborate on this matter and provide a more precise statement. For example, given sample covariance matrix, the computational cost associated with calculating the thresholded graph, bridges, and clusters is relatively low [1,2]. \n\nFor a detailed discussion on why the computational cost associated with determining the thresholded graph, bridges, and clusters is much lower compared to solving sub-problems, please refer to Part 3 of our global response.\n\n__References__:\n\n[1] Jens M Schmidt. A simple test on 2-vertex-and 2-edge-connectivity. Information Processing Letters, 113(7):241–244, 2013.\n\n[2] R Endre Tarjan. A note on finding the bridges of a graph. Information Processing Letters, 2(6):160–161, 1974.\n\n## Answer to Questions Part 2\n\n> I actually had to read it a couple times to understand the exact procedure that was being proposed. I suggest adding some clarifying details to help the reader.\n\n__Reply:__ We appreciate your suggestion that the procedure could benefit from additional clarification. Please find our detailed answer in the part 3 of the global response. \n\n## Answer to Questions Part 3\n\n>What is the primary limitation that makes this result only apply to MTP2 distributions? Is it possible that something similar could be true more generally? I didn't really get a sense of why MTP2 was required in the main text.  \n \n__Reply:__ We appreciate the reviewer for raising this interesting question. In part 1 of our global response, we highlight the importance of MTP2 in ensuring the validity of our proposed Theorem 3.3. This also explains why our results only apply to MTP2.\n\nIn response to whether our methodology can hold more broadly, in part 2 of our global response, we detail the circumstances under which our theoretical framework could be applied to the graphical lasso as well. Please check our detailed answer. \n\n## Answer to Questions Part 4\n> Do you need Gaussian distributions here? Can something be done for more general MTP2 distributions, e.g., those that can still be represented as pairwise graphical models? \n\n__Reply:__ Thank you for your insightful comment. While our paper primarily focuses on learning $\\mathrm{MTP}_2$ Gaussian graphical models, our method is not explicitly tied to Gaussian distributions. This is due to the deterministic nature of our optimization problem, as presented below:\n\n$ \\min_{\\boldsymbol{\\Theta}\\succ \\mathbf 0} - \\log \\det ( \\boldsymbol{\\Theta}) + \\langle \\boldsymbol{\\Theta}, \\mathbf{S} \\rangle + \\sum_{ij} \\Lambda_{ij} | \\Theta_{ij}|.$\n\nExtending our method to non-Gaussian cases may be straightforward. For instance, in the context of elliptical distributions, we could substitute the sample correlation matrix with the Kendall's tau correlation matrix, leading to positive partial correlation graphs [1]. It's worth highlighting that the MTP2 property is somewhat restrictive, and as demonstrated in [2], an elliptical distribution that is MTP2 across all dimensions is essentially Gaussian. Lastly, we would like to clarify that the sole assumption in our paper is that $S_{ij} < \\sqrt{S_{ii} S_{jj}}$ for any $i \\neq j$. \n\nWe will include a more detailed discussion on this aspect in our revised manuscript.\n\n__References__:\n\n[1] R. Agrawal, U. Roy, and C. Uhler. \"Covariance matrix estimation under total positivity for portfolio selection.\" Journal of Financial Econometrics, 20(2):367-389, 2022.\n\n[2] D. Rossell, and P. Zwiernik. \"Dependence in elliptical partial correlation graphs.\" Electronic Journal of Statistics, 15(2):4236-4263, 2021.\n\n## Reply to Comments\n\n>  While it is true that MTP2 distributions have some practical utility (I've used them myself), it isn't clear that many real datasets meet that requirement (even approximately)\n\n__Reply:__ Thank you for your thoughtful questions. From a practical perspective, MTP2 becomes a natural choice when the __variables are anticipated to display positive dependence__. Numerous real-world scenarios, particularly in sectors such as finance and social sciences, manifest this positive dependence.\n\nTo illuminate how data expressing positive correlation may approximate MTP2 properties, we have conducted supplementary experiments with the CROP Image Data set, which is used in our real-world experiment. The details are described in Part 4 of our global response. These experiments reveal that the CROP data set shows a form of positive correlation, thereby making it suitable for modeling with MTP2 graphical models. Please refer to our global response for our detailed replies.\n\n## Additional Comments\n__Reply__: Thank you for pointing out the typos. We will correct them in the revised version."}]

Author rebuttal:
## Part 1: Answers to Questions Regarding Roles of MTP2 

> We gather this question from Reviewers ghhX and pk8s, who sought the role of the MTP2 constraint in our main results and why these only apply to MTP2 distributions.

__Reply__: We thank the reviewers for bringing up this interesting point. Our research offers two main theoretical advancements. The first proposes an explicit form for the inverse of $\boldsymbol{\Theta}$, while the second confirms that $\boldsymbol{\Theta}=\mathbf{R}^{-1}$ meets optimality without extra conditions.

__The MTP2 properties are sufficient conditions for the second contribution to hold__. Hence, while the bridge-block decomposed form of $\mathbf{R}$ could be broadly applied to other graphical models, our main findings are only applicable to MTP2 graphical models.

Technically speaking, the MTP2 constraints eliminate the non-smoothness, thus __simplifying the KKT conditions__ as follows: 

||Graphical Lasso|MTP2|
|:-|:-|:-|
|$\forall i$| $-R_{ii}+S_{ii}=0,$|$-R_{ii}+S_{ii}=0,$|
|$\forall\Theta_{ij}\neq0$|$-R_{ij}+S_{ij}+\lambda_{ij}\text{sign}\left(\Theta_{ij}\right)=0,$| $-R_{ij}+S_{ij}-\lambda_{ij}=0$|
|$\forall\Theta_{ij}=0$|$\vert -R_{ij}+S_{ij}\vert\leq \lambda_{ij}$|$-R_{ij}+S_{ij}-\lambda_{ij}\leq0$|

Simultaneously, $\mathbf{R}$ becomes non-negative ($\geq\boldsymbol{0}$). As a result, The most challenging part $-R_{ij}+S_{ij}-\lambda_{ij}\leq 0$ of the KKT conditions holds under MTP2 properties. (The details refers to Section B of our Appendix.)

## Part 2: Answer to Questions Regarding Generalizing Our Results to Graphical Lasso

> Reviewers ghhX and pk8s provided feedback regarding the applicability of our findings to the graphical lasso.

__Reply__: We thank reviewers posing this thought-provoking question. In the context of the graphical lasso, our current, yet-to-be-published research indicates that our decomposed form achieves optimality if the following inequality is satisfied:

$| -R_{ij} +S_{ij}|\leq \lambda_{ij},\quad\forall (i,j)\notin \mathcal B\text{ and }i,j \text{ belong to different clusters,}$ 

where $\mathbf{R}$ is computed by Theorem 3.4. This condition could be easily established in certain cases, such as when $\lambda_{ij}\gg 1$. Nevertheless, in the majority of practical situations, verifying this condition poses a substantial challenge.

Conversely, these conditions are inherently satisfied when MTP2 constraints are applied. We hope that our theorem will instigate future studies to simplify these conditions for the graphical lasso. 

## Part 3: Answers to Questions Regarding the Exact Procedures of Proposed Framework

> We receive this question from Reviewers vWNe and pk8s, who seek to comprehend the exact procedures and processing time of the proposed method.

__Reply__: We'll use an SBM graph with p=5000 as a reference. The specifics are outlined in the following:

__Preprocessing__: In this phase, we have three steps:

|Step|Elements to Compute|Time (s)|
|:-:|:-|:-:|
|1| Derive the thresholded matrix $\mathbf{T}$ from $\mathbf{S}$ and $\boldsymbol{\Lambda}$ |0.01|
|2|Extract the set of all bridges $\mathcal B$ from the thresholded graph |0.21|
|3|Determine the bridge-block decomposition|0.12|

__Solving Subproblems__: For all clusters, we solve the associated sub-problems.

__Computing Optimal Solution via Theorem 3.3__: The optimal $\boldsymbol{\Theta}$ is then derived as follows:

|Conditions|Formulas|Cost|
|:-|:-|:-|
|$i,j\in \mathcal V_k$| Derive $\Theta_{ij}$ from $\widehat{\boldsymbol{\Theta}}_k$|Depends on specific algorithms (BCD: 967s; FPN: 101s).|
|$(i,j)\in\mathcal B$|$\Theta_{ij}= -T_{ij}/(S_{ii}S_{jj}-T_{ij}^2)$|0.12s|
|otherwise|$\Theta_{ij}=0$|0s|

The findings reveal that the additional costs, such as the time for pre-processing and the employment of Theorem 3.3, are considerably less compared to the time invested in solving sub-problems.

We find bridges using a bridge-finding algorithm from [1]. This method uses a depth-first search, which means the complexity is $\mathcal{O}(|\mathcal V| +|\mathcal E|)$. In the sparse graphs we're interested in, the number of edges $|\mathcal E|$ usually scales similarly to the number of nodes $|\mathcal V|$. Hence, the computational cost of bridges for high-dimensional sparse graphs is low.

[1] R Endre Tarjan. A note on finding the bridges of a graph. Information Processing Letters, 2(6):160–161, 1974.

## Part 4: Experiments on Testing MTP2 in CROP Dataset

> In our paper, we mainly focus on how proposed method acclerate the learning of MTP2 graphical models on CROP Data set. It is suggested by the reviewer ghhX that we should also justify why MTP2 constraint is plausible in this problem. To address this, we present additional experiments to tackle this concern.

__Reply__: We selected 20 random subsets from the CROP dataset. For each subset, we computed the Graphical Lasso and the MTP2 graphical model for different values of $\lambda$ using the first 10 observations. The remaining 36 observations were used to calculate the out-of-sample log-likelihood, which was then averaged across all datasets. This process allows us to evaluate how well these models generalize to unseen data.

As depicted in Figure 5 of the attached PDF, the MTP2 graphical model outperforms the Graphical Lasso, providing a higher test log-likelihood. We present one instance of the estimated graphical lasso model in Figure 6. It reveals that __most conditional correlations are positive__ (red edges, 90%), with a few being negative (blue edges, 10%). This pattern implies strong positive dependence in the CROP data, aligning with the characteristics of MTP2.

These results are not unexpected, given that the CROP dataset comprises multiple clusters. Within the same cluster, we expect data points to exhibit greater similarity compared to those in different clusters. This situation signifies a form of positive dependence, thereby justifying the plausibility of the MTP2 assumption in this problem.

Meta-review: This paper studies methods to more efficiently learn large Gaussian graphical models from data, specifically in the context of the MTP2 (attractiveness/positive correlation) assumption and by attempting to take advantage of inherent 'bridge' structure in the data. This is a fairly specific setting (for example, the MTP2 assumption is reasonably well-motivated but probably not exactly satisfied in most real life datasets), but after the discussion overall the reviewers appreciated the contribution including the experimental validation. Based on this, I recommend acceptance.
Decision: Accept (poster)