LogSpecT: Feasible Graph Learning Model from Stationary Signals with Recovery Guarantees
Authors: Shangyuan Liu
Linglingzhi Zhu
Anthony Man-Cho So
Venue: NeurIPS 2023 poster (2023)
Keywords: Graph Signal Processing
Spectral Template
Network Inference
Optimization
Linearized ADMM

Abstract: Graph learning from signals is a core task in graph signal processing (GSP). A significant subclass of graph signals called the stationary graph signals that broadens the concept of stationarity of data defined on regular domains to signals on graphs is gaining increasing popularity in the GSP community. The most commonly used model to learn graphs from these stationary signals is SpecT, which forms the foundation for nearly all the subsequent, more advanced models. Despite its strengths, the practical formulation of the model, known as rSpecT, has been identified to be susceptible to the choice of hyperparameters. More critically, it may suffer from infeasibility as an optimization problem. In this paper, we introduce the first condition that ensures the infeasibility of rSpecT and design a novel model called LogSpecT, along with its practical formulation rLogSpecT to overcome this issue. Contrary to rSpecT, our novel practical model rLogSpecT is always feasible. Furthermore, we provide recovery guarantees of rLogSpecT from modern optimization tools related to epi-convergence, which could be of independent interest and significant for various learning problems. To demonstrate the practical advantages of rLogSpecT, a highly efficient algorithm based on the linearized alternating direction method of multipliers (L-ADMM) that allows closed-form solutions for each subproblem is proposed with convergence guarantees. Extensive numerical results on both synthetic and real networks not only corroborate the stability of our proposed methods, but also highlight their comparable and even superior performance than existing models.

LogSpecT: Feasible Graph Learning Model from Stationary Signals with Recovery Guarantees
Shangyuan Liu The Chinese University of Hong Kong shangyuanl cuhk, edu. hk
Linglingzhi Zhu The Chinese University of Hong Kong 1lzzhu@se cuhk. edu. hk
Anthony Man-Cho So
The Chinese University of Hong Kong manchoso8se. cuhk, edu. hk
Abstract
Graph learning from signals is core task in graph signal processing (GSP). A significant subclass of graph signals called the stationary graph signals that broadens the concept stationarity of data defined on regular domains to signals graphs is gaining increasing popularity in the GSP community [10, 20, 25]. The most commonly used model to learn graphs from these stationary signals is SpecT' [31], which forms the foundation for nearly all the subsequent, more advanced models. Despite its strengths, the practical formulation of the model, known as rSpecT, has been identified to be susceptible to the choice of hyperparameters. More critically, may suffer from infeasibility as an optimization problem. this paper, we introduce the first condition that ensures the infeasibility of rSpecT and design a novel model called LogSpecT, along with its practical formulation rLogSpecT to overcome this issue. Contrary to rSpecT, our novel practical model rLogSpecT is always feasible. Furthermore, we provide recovery guarantees of rLogSpeeT from modern optimization tools related to epi-convergence, which could of independent interest and significant for various learning problems. To demonstrate the practical advantages of rLooSpecTT, highly efficient algorithm based on the linearized alternating direction method multipliers (L-ADMM) that allows closed-form solutions for each subproblem is proposed with convergence guarantees. Extensive numerical results on both synthetic and real networks not only corroborate the stability of proposed methods, but also highlight their comparable and even superior performance than existing models.
1 Introduction
Learning with graphs has proved its relevance in many practical areas, such as life science [36, 37], signal processing [13, 14, 38, 39], and financial engineering [1, 21], just name few. However, there are many cases that the graphs are not readily prepared and only the data closely related to the graphs be observed. Hence, core task in Graph Signal Processing (GSP) is learn the underlying graph topology based on the interplay between data and graphs [22].
Among these commonly used data properties, assumption that is referred to the graph signal stationarity [10] is gaining attention [30, 20, 25J. Graph signal stationarity the extension of the notion of stationarity from the regular space/time domain to the irregular domain (i.e. graphs). The latter traditional hypothesis in signal processing to capture a special type of statistical relationship between samples of temporal signal. It is cornerstone of many signal processing methods. Although the graph signal stationarity property is proposed more from a theoretic end, the
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

techniques based on this idealistic assumption pave the way for the development of the following more advanced and applicable models [35, 34]. Furthermore, several works have shown that some real datasets can be approximately viewed as stationary or partially explained by the stationarity assumption. For instance, [25] revealed that the well-known USPS dataset and the CMUPIE set of cropped faces exhibit near stationarity. [10] found that certain weather data could be explained by stationary graph signals. Additionally, [31] highlighted the relevance of the stationarity assumption in learning protein structure. [43] applied the stationary assumption to US Senate Roll Calls and achieved good performance.
The predominant methods to process stationary graph signals and learn topology under the stationarity assumption are the -eterraltemplaeeeeeee baseed models. The start of this line of works is [31], which proposed a vanilla model called SpecT [31] to learn graphs from stationary graph signals. Many extensions of this fundamental model have been made since then to accommodate the more realistic and complicated settings, e.g. joint inference of multiple graphs [32, 27], learning with streaming graph signals [33], learning with hidden nodes [4], learning from both smooth and stationary graph signals [5]. In practice, SpecT requires the unknown data covariance matrix. Hence, robust formulation called rSpecT is proposed [31], which replaces the unknown covariance matrix with its estimate and introduces hyperparameter to reflect the estimation inaccuracy. Unfortunately, the model is sensitive to this hyperparameter and improper tuning of it may jeopardize the model performance or lead to model infeasibility.
The current approach to selecting an appropriate value is heuristic. It finds the smallest one that allows for feasible solution. This method has two shortcomings. On one hand, It is computationally costly. The smallest value is the exact solution to a second-order conic programming. which is solved approximately with computational methods. Any approximate solution that is smaller than the exact one will make rSpeeT infeasible. On the other hand, such an approach lacks theoretical guarantees as how the smallest value behaves as the number of observations increases is unknown. This issue makes the existing recovery guarantees invalid. To the best of our knowledge, only a few works have focused on how to understand the performance of the robust formulation rSpecTT. [24] empirically showed that the feasible set of the robust formulation approximates the original one as more samples are collected. [23] provided the first recovery guarantee for rSpecT. Their work exploits the fact that the objective is the l fi-orrm and the constraints are linear [44]. However, the conditions needed in their work are not only restrictive but also hard check, as they require the full rankness of large-scale matrix related to random samples. Another model [34] tries to circumvent the model infeasibility issue by turning a constraint into penalty. However, this approach introduces another hyperparameter that is neither easy to tune nor amenable to theoretical analysis. Also, this hyperparameter should go to infinity heuristically when more and more samples are collected, which may exacerbate model instability.
In this paper, we propose an alternative formulation of SpecT that can bypass the aforementioned issues. We first provide condition that ensures the infeasibility of the fundamental model rSpecT. Then, we propose a novel formulation called LogSpecT and its practical formulation rLogSpecT based the spectral template and incorporating log barrier. The inclusion of a log barrier regularizer on the degrees is inspired by its extensive application in graph learning from smooth graph signals [15, 16, 19, 47] and the discussion in [7] emphasizing that restricting the degrees to preclude isolated nodes is crucial to learning graphs with special structures. It is worth noting that the approach in [7], which directly imposes constraints, may exacerbate the infeasibility issue in our problem. Although rl rLogSpecT still requires a parameter to reflect the estimation error of the unknown covariance matrix, it longer suffers from the infeasibility issue associated with the improper choice of the parameter. Furthermore, we provide the decaying rate of the parameter that can guarantee the recovery of rLooSPecT. This not only aids in the parameter tuning, but also establishes a connection between rLogSpecT and LogSpecT through non-asymptotic recovery analysis. These theoretical results are derived using contemporary optimization tools related to epi-convergence [29]. Unlike the current guarantees for rSpecT that are built on an l analysis model [44], our approach based on ppieeenveeeence is not confined to the specific types of optimization problems (e.g. the combination of the & f--ormm and log barrier in the objective function) and consequently admits broader applications that can be of independent interest.
In the algorithmic aspect, we design a linearized alternating direction method of multipliers (L- ADMM) to solve rLogSPeTT. The subproblems of L-ADMM admit closed-form solutions and can be implemented efficiently due to the linearization technique. Also, we provide the convergence result of the proposed method. Finally, we conduct extensive experiments on both synthetic data and

---

real networks. Notably, the infeasibility issue of rSpecT frequently occurs in experiments on real networks, thereby emphasizing its practical significance beyond theoretical concerns. Furthermore, we observe that even when the classical models are set to be feasible, our novel models (LogSpecT and SpeeT)) exhibit improved accuracy and stability, particularly in the case of synthetic BA graphs. This empirically illustrates that although our models are initially proposed to circumvent the infeasibility issue, they are superior alternatives to classical methods the general cases.
1.1
Notation
The notation we use in this paper is standard. We use [m]] to denote the set {1, 2,....... ,..,, for any positive integer m. Let the Euclidean space of all real matrices be equipped with the inner product (X,Y) tr(x Y) for any matrices X, Y and denote the induced Frobenius norm by -|-1| (or -1 when the argument
vector). For any X E Rm*"" we use (B p) to denote the closed Euclidean ball centering at X with radius Let || be the operator norm, || = Σ., [Xay) |x maxi.j [x,g] and let diag(X) be vector formed by the diagonal entries of X. For column vector æ, let Diag(æ) a) be the diagonal matrix whose diagonal elements are given by æ. Given closed and convex set C, we use (wc(w) == argmin / w|l2 to denote the projection of the point w onto the set We use 1 (resp. 0) to denote an all-one vector (resp. all-zero vector) whose dimension will be clear from the context. For set S and any real number a, let a& == {ar x E S},
S1 == {x1 x E S), and (s(-) (.) be the indicator function of S. For two non-empty and compact sets X' and D, the distance between them is defined as dist(X\)) supxex infyey |X -Y|l mexxxxiiffxey
2 Preliminaries
Let G - E) be a graph, where V = [m]] is the set of nodes and & C V V is the set of edges. Let S be the weight matrix associated with the graph where S,j represents the weight of the edge between nodes and j. In this paper, we consider undirected graphs without self-loops. The set of adjacency matrices for such graphs is S {S E mm S = S (iag(S) = N 0}. Suppose that the adjacency matrix S admits the sigen-decomposition S UAU where A is diagonal matrix and U is an orthogonal matrix. A graph filter is a linear operator : mm - Rmxm defined as h(S) Di S U (E h;,A)) UT, where p 0 is the order of the graph filter and {h+ are the filter coefficients. According to the convention, we have 5 I.
graph signal can be represented by a vector æ E R", where the i-th element r, is the signal value associated with node i. A signal is said to be stationary if it is generated from
x l(S))
where satisfies E[w] 0 and E[w"T T wi of ær, which is denoted by Coc, shares the same eigenvectors with S. Hence, we have the constraint
I. Simple calculations give that the covariance matrix
Cx =S
Based on this, the following fundamental model SpeeT is proposed to learn graphs from stationary signals without the knowledge of graph filters [31, 34]:
Sij = =1] (SpecT) min || s.t. Cx SC,. S Sesn Rmxm -,, S
where the constraint Σ/= Si = 1 is used to preclude the trivial optimal solution S = 0. When C.o is unknown and only n i.id samples of {æ,} from (1) are available, the robust formulation, which is based the estimate Cn of C.c and called rSpecT, is used:
min || s.t. ||C,S - SC S 6, S SESn S
TR T
>, S1, 115=1 (rSpecT)
For this robust formulation, the recovery guarantee is studied empirically [24] and theoretically in [23] under some conditions that are restrictive and hard to check.
It reduces to the classic point-to-set (resp. point-to-poiin distance when y (resp. X and y) is singleton.
3

---

3 Infeasibility of rSpecT and Novel Models
Even though rSpecT has gained much popularity, there are few works that discuss the choice of the hyperparameter 8 and how it affects the model feasibility. In this section, we present condition under which rSpecT is always infeasible and then propose an alternative formulation. To motivate our results, let us consider the following 2-node example.
Example 3.1. Consider a simple graph containing 2 nodes. Then, the set given by the second constraint of rSpecT is S (S |R=x -,== Stj = 1] [0, 1, which is a singleton. Suppose that the sample covariance matrix is C, [hn 1 h12 hz2) Then, the constraint
||C,, - -Cr S 6 is reduced to 2(h11 h22) < 5². Hence, when h1 # h22 and 5 < v2 h11 h22), rSpecT has no feasible solution.
Before delving into the general case, we introduce linear operator B p "pm"xm(m-D 1)/2 that maps the vector æ E 4 (pmm-1)m/2 to the vectorization of an adjacency matrix S E S of simple, undirected graph:
if i j, F Pij+ 7 (2m -)
(Bx) (Bx) mm)))(11)))
0,
if i= j
(3)
1 j- 22 (2m-1))
if i i
where i, [m]] We also define A,, = I -n C, 2 E Rm xm² that the first constraint of rSpecT can be rewritten as |C,S SC, < < A vec(S) S 8, where vec(.) is the vectorization operator. We now give condition that guarantees the infeasibility of rSpecT.
Theorem 3.2. Consider the linear system
A, = y 0S0,
E Di + O,
m--
(4)
where y E IR infeasible for all &n E [0, 8.).
m)
If (4) has no feasible solution, then there exists a on 0 such that rSpecT is
Remark 3.3. From Theorem 3.2 we can infer that for any fixed n, the linear system (4) has no solution when A, has full column rank (e.g., Example 3.1). This leads to the infeasibility rSpecT with 5n E [0, 6.). As remarked in [31], J, one should tackle the feasibility issue of rSpecT with caution.
The failure of rSpecT (SpecT) lies in the existence of the constraint (S1); to preclude the trivial solution S = O. For this reason, we resort to an alternative approach to bypassing the trivial solution. When the graphs are assumed to have isolated node, the log barrier is commonly applied [15]. In these cases, the zero solution is naturally precluded. This observation inspires us to propose the following novel formulation, which combines the log barrier with the spectral template (2) to learn graphs without isolated nodes from stationary signals:
which is used
min || -al log(S1)
(LogSpecTT)
s.t. SeS, SC. S.
where || || is convex relaxation of || llo promoting graph sparsity, 1 log(S1) is the penalty to guarantee the nonexistence of isolated nodes, and is the tuning parameter. As can be seen from the following proposition, the hyperparameter a in LogSpecT only affects the scale edge weights instead of the graph connectivity structure.
Proposition 3,4. Let Opt(C,a) a) be the optimal solution set of LogSpecT with input covariance matrix C and parameter a Then, for any y > 0, it follows that
Opt(C,a)= *pt
Opt(C,1).
Remark 3.5. The result of Proposition 3.4 spares us from tuning the hyperparameter when we are coping with binary graphs. In fact, certain normalization will eliminate the impact of different values of a and preserve the connectivity information. Hence, may simply set a: 1 in implementation.
Note that the true covariance matrix Coo LogSpecTT is usually unknown and an estimate C,, from i.id samples is available. To tackle the estimation inaccuracy, we introduce the following robust

---

formulation:
min ||| ad log(S1)
(rLogSpecT)
s.t. Se |SC,, - S
This formulation substitutes C. with C, and relaxes the equality constraint to an inequality con- straint with tolerance threshold ón. Contrary to rSpecT, we prove that rLooSSecT is always feasible.
Proposition 3.6. For 6,, 0 and C, with any fixed n, rLogSpecT always has nontrivial feasible solution.
Recovery Guarantee of rLogSpecT
In this section, we investigate the non-asymptotic behavior of rLogSpecT when and more i.i.d samples are collected. For the sake brevity, we denote f = Il1.1 - al log (-1). The theoretical recovery guarantee is as follows:
Theorem 4.1. If 6 2 2p - Coell with given = max am, am(log a - 1)), then there exist constants C1, Ci; 02 such that
i - Eni
(ii) dist(Sm**, (2e,) S Eni (iii) dist S EiEn + 02
where En ci (8n 2p Cn Coo []), Su.* the optimal solution set of rLogSpecT, f (resp. fi) denotes the optimal value of LogSpecT (resp. rLogSpecT), and
S = (S E S| SCxx C. ((S) S +
is the E-suboptimal solution set of LogSpecT.
Remark 4.2. (i) Compared with the conclusion (ii) in the Theorem 4.1, the conclusion (iii) links the optimal solution sets of rLogSpecT and LogSpecT instead of the sub-optimal solutions.
(ii) A byproduct of the proof (see Appendix E.2 for details) shows that the optimal node degree vector set (S) is singleton. However, there is guarantee of the uniqueness of S itself. We will discuss the impact of such non unn-uiqqueesss on model performance in Section 6.4.
Remark 4.3. The proof of Theorem 4.1 relies on two important optimization concepts: Epi- convergence and truncated Hausdorff distance. Epi-convergence is closely related to the asymptotic solution behavior of approximate minimization problems, and the truncated Hausdorff distance is used to characterize the epi-comvereence non-asymptotically: With the help of the Kenmochi condition, which allows us to explicitly calculate the truncated Hausdorff distance, we are able to study the non yynptssmptoicc behavior of the optimal value and optimal solutions of the models. We refer the reader to Appendix E and also [28, Chapter 7J, [29, Chapter 6.J] for more details.
Corollary 4.4. Under the assumptions in Theorem 4.1, follows that if 8n 0, then
lim dist(S".*, ,8)
(5)
Remark 4.5. From the strong law of large numbers, if we choose 8n = O(IICn -Coll), then almost surely one has lim || Cxall and consequently 8n O. Hence, Corollary 4.4 shows that (5) holds almost surely.
For the remaining part of this section, we study the choice of 8n under certain statistical assumptions. A large number of distributions (e.g., Gaussian distributions, exponential distributions and any bounded distributions) can be covered by the sub-Gaussian distribution (cf. [41, Section 2.5]), whose formal definition is as follows.
Definition 4.6 (Sub-Gaussian Distributions). The probability distribution of random vector up is called hab-Gaussian if there are C, u > 0 such that for every t > 0, PP[[[]]] > t) S Ce -et
Consider the case that w in the generative model (1) follows sub-Gaussian distribution. The following result is adapted from [40, Proposition 2.1].
Lemma 4.7. Suppose that w in the generative model (1) follows sub-Gaussian distribution. Then, x follows a sub-Gaussian distribution, and with probability larger than - O(2), one has ||C, Cooll OOvllo n//n.
5

---

Equipped with non-asympiotic results in Theorem 4.1, we choose on with a specific decaying rate. Corollary 4.8. If the input in (1) follows uab-Gaussian distribution and 8n O log n//). then the assumptions in Theorem 4.1 I hold with probability larger than 1
Corollary 4.8 together with Theorem 4.1 illustrates the non-asymptotic convergence of optimal function valuc/suboptimal solution set/optimal node degree vector of rLogSpect from i.id samples to the ideal model LogSpect with the convergence rate ((s))/0(ew)/O() where
En O(IC, C I Olylog n//n) for sub-Gaussian w. This matches the convergence rate O / log / of classic spectral template models, e.g., Proposition [31] and Theorem 2 in [23] for the SpecT model, which shows that LogSpecT is also competitive on recovery rate.
5 Linearized ADMM for Solving rLogSpecT
In this section, we design a linearized ADMM algorithm to solve rLogSpecT that admits closed-fome solutions for the subproblems. The ADMM-type algorithms have been successfully applied to tackle various graph learning tasks [46, 42].
Firstly, we reformulate rLooSSpecT such that it fits the ADMM scheme:
min (S. 11 a1 log q + (s(S) +B(O,6.. 8.2.4
(6)
s.t. C - SC,, = Z, S1.
The augmented Lagrangian function problem (6) is
2(S, A,)2) = (11T, S) - -x11 log + (s(S) + (300,..) + (A, C.S - SC,, - Z)
+ > (q 51) - PlIC,sS S - SC, Z + P 114 S1
2
Update primal variables and q Since Z and g are separable given S, we update the primal variables in two blocks: S and (Z, q). More specifically, the update of (Z, q) in the k-th iteration is
z (k+1) (+) == argmin (2.9) ((S Z, A (K) A (1))
which admits the following closed-form solution:
Proposition 5.1. The update (7) can be explicitly rewritten as
1, 8n
(k++1) g g2 4a/ p1
((++1) min (
2,
|1
where 2 = (()) g(k) + A
q 1- 1X
Update of primal variables S The ADMM update of does not admit closed-form solution. Hence, we solve the linearized version La of the augmented Lagrangian function L when all the variables except S are fixed. It is given by
(S( = (K) S)
PT |s - ((k) (S()
2
where (k) = 11 + C, - A (+)
x 1 plS C + C,S(M) - S(k)) Z (1+1) - C, 2 (k+1) - ( 11 and is the parameter for linearization. Therefore, the update of the k-tth iteration is
g(++1) = argmin ( = (s(((() - C(kk)/p))
(8)
where (X) is the projection of X onto the closed and convex set S and can be computed by
(Is(X)) = v max (o, X. + XzJ). if if
0,
if i=

---

Update of dual variables and parameters Finally, we update the dual variables A and X via gradient ascent steps of the augmented Lagrangian function:
(k+1)= == A (K) + p(C,, S(k+1) g(++)) C. - z((+1))
(9)
X [++1) = 2 p(q (k+1) 1).
(10)
Together with the updating rule for the augmented parameter p. we complete the whole algorithm. The details are presented in Appendix F.
Complexity of each iteration To illustrate the efficiency of L-ADMM, we calculate the complexity of each subproblem. In the update of Z and q, the most time consuming part is to calculate C2 When the number of nodes m is larger than the sample size n, we rewrite (()) as (1/xX)(XT S(A)), where X E R the data matrix. this case, the complexity is in the order of (mm²)) When m is smaller than n, we first calculate C, and then multiply it with ((A), This gives the complexity Onnm²). Similarly, the complexity needed to update S is O + mn²). The update of the dual variables can be decomposed as the summation of components calculated when updating the primal variables. In summary, the total complexity of each iteration is Onmm² + +mn?).
Convergence of L-ADMM For the convergence analysis of L-ADMM, we treat (Z,g) one variable, and the two constraints S, SC,, Z and S1 can be written into single one:
& C - C, & I; & J vec(S) (S) = [vec( (rec()):0 q].
Then can apply two-block proximal ADMM (i.e., L-ADMM) the problem by alternatingly updating S (see details in Appendix F.4 for this part) and (Z, q). Consequently, the convergence result in [45, Theorem 4.2] can be invoked directly to derive the following theorem.
Theorem 5.2. If > m + | & Im Im @ 1, then lima- f(S(A) = fñ.
Experiments
In this section, we evaluate our proposed algorithm and models via numerical experiments. The experiments are conducted on both synthetic and real networks. We apply the standard metrics [18]: F-measure, Precision, and Recall to assess the quality of learned graphs. The source code is available at: https://github: cot.//teeaass/lteeuttPooo00 LogSpecT.
6.1 Data Generation
Random Graphs. In the experiments, we consider two types of synthetic graphs, namely, the Erdōs-Rény: (ER) graph [9] and the Barabasi-Albert model graph (BA) [2]. The ER graphs generated by placing an edge between each pair of nodes independently with probability p = 0.2 and the weight on each edge is set to 1. The BA graphs are generated by having two connected nodes initially and then adding new nodes one at a time, where each new node is connected to exactly one previous node that is randomly chosen with a probability proportional its degree at the time.
Graph Filters and Signals. Three graph filters are used in the experiments. The first one is the low-pass graph filter (lowpass-EXP) h(5)) exp( 1 {)). The second one is the high- figh-passsa graph filter (highpass-EXX EXP) h(S) cxp(-S). The last one is quadratic graph filter (QUA) h(S) S + I. Note that the quadratic graph filter is neither low-pass nor high-pass. We will also use random graph filters developed from the above ones for experiments on real networks. They will be introduced when applied. The stationary graph signals are generated from the generative model (1). The input graph signal ~ N/ mm) is random vector following the normal distribution.
6.2 How to Infer Binary Graphs
Many practical tasks require binary graphs instead of weighted graphs. However, the results from LogSpect and rLogSpecT are generally weighted. In this section, we tackle the issue of how to convert a weighted graph to binary one. Firstly, we normalize each edge in the graph W by dividing it over the maximal weight, which yields collection of values ranging from 0 to 1. Our task is then to select a threshold & to round these values to and construct a binary graph W accordingly. Mathematically, this procedure is captured by the formula:
7

---

1, if W max wey) 2 C, (W"), = 0. W max Wey) < <.
To choose the threshold E, we either use a training- based strategy or searching based strategy. The training- based one searches the best on k graphs and applies the value to the newly learned graph. The searching-baser one simply searches for that yields the best performance for graph.
6.3 Convergence of L-ADMM
We present the convergence performance of L-ADMM in this section. Since there does not literally exist customized algorithm for rLogSpecT, we compare our L-ADMM with CVXPY [8]. The solver used in CVXPY is MOSEK. We conduct the experiments 10000 stationary graph signals generated from low-pass EXP graph filter. The parameter 5 is set as /lg n/n. For L-ADMM algorithm, the target accuracy is set as 10-6 and the initialization set as zero. The running time for L-ADMM takes around 20 seconds while the solver takes over 110 seconds. The primal residual and dual residual of L-ADMM in each iteration are plot in Figure The results corroborate the efficiency of our proposed algorithm. Furthermore, the locally linear convergence can be observed from the experiment. This implies that L-ADMM converges with fast speed to high accuracy. remains an open question to prove the property theoretically.
100-node BA graphs with
Residuals
rLogSpecT Performance
primal residual
0.8
E-meapure
0.6
10-1
N
0.2
lowpass-EXK
highpass
QUA
10a
10a
10a
10"
10s
16a
iteration
Sample Size
Figure 1: L-ADMM residual on 100-node graph
Figure 2: rL.ogSpecT with on = 0.2 /l///og//// n/n
El Stppoppass-EPP
Eit-qUA
Bi-highpass--:0
Nighighoazs
LogSpecT
3.0
6.8
Avgrage megaire
2.6
0.6
0.4-
0.4
0.2
0.0
20
30
40
50
60
0.0
20
30
50
60
Figure 3: Comparison of LogSpecT & SpecT on Synthetic Data
6.4 Experiments on Synthetic Networks
To evaluate the efficacy of LogSpecT and rLogSpecT, we compare them with SpecT and rSpecTT on synthetic data. Notice that the prior of no isolated nodes has been incorporated in practice by adding restrictions on the smallest degree. This approach introduces additional hyperparameter and may exacerbate the infeasibility issue robust formulation. Also, no recovery guarantees are developed for this trick. For these reasons, we still compare LogSpect & rLogSpecT with SpecT & rSpecT.

---

Performance of LogSp We first compare the performance of the ideal models: LogSpeeTT and SpecT on the two types of random graphs with three graph filters. We conduct experiments on graphs with nodes ranging from 20 to 60 and use the searching-basee strategy to set threshold s*. After 10 repetitions, we report the average results in Figure The left column presents the performance of LogSpecTT and the right column presents that of SpecTT Firstly, we observe that on ER graphs, both models achieve good performance. They nearly recover the graphs perfectly. BA graphs differs. In this case, LogSpecT can work efficiently. The learned graphs from LogSpecT enjoy much higher F-measure values than those from SpecT. Secondly, the comparison between different graph filters shows that different graph filters have few impacts on recovery performance for both LogSpecT and SpecT. Finally, we observe that the outputs of LogSpeeT on BA graphs tend to possess higher F-measure when the number of nodes increases. This suggests that SpeeT may behave better on larger networks. However, such phenomena cannot be observed from Spec'T results.
However, the performance
Performance of rLooSSecT. Since LogSpec'T persuasively outperforms SpecT,, we only present the performance of rLogSpecT when more samples are collected. As we have mentioned before, the optimal solution to LogSpecT is not necessarily unique. Thus, we do not expect to be able to show how PoegSpe optimal solution converges to LogSpecT's optimal solution. Moreover, the non-uniqueness may jeopardize the performance of rLogSpecT intuitively. However, the following experiment on BA graphs shows that the learned graphs from rLogSPecT with the Ainary-grahle transforming strategy tend to approach the ground truth when enough samples are collected.
In this experiment, the sample size n is chosen from 10 to 10 and 6, is set as 0.2 / /2ooo///// We rely on the training-basee strategy to obtain the best threshold e* from 10 randomly chosen training graphs. We then calculate the F-measure of the learned binary graph in the testing set. The result is reported in Figure 2. It shows that rLogSpecT works for all three graph filters on BA graphs and tends to possess higher F-measure values when more and more signals are collected. This indicates that similar to LogSpecT, rLogSpecTT is efficient for different types of graph filters BA graphs, including the high-pass ones. For more experiments and discussions, we refer readers to Appendix G.
6.5 Experiments on Real Networks
In this set of experiments, we compare the performance of LogSpecT (resp. rLooSSecTT) with SpecT (resp. rSpecTT) and other methods from the statistics and GSP communities on Protein database and Reddit database from [6]. The Protein database is a Bioinformatics dataset, where nodes are secondary structure elements (SSEs) and there is an edge between two nodes if they are neighbors in the amino-acid sequence or in 3D space. The Reddit database is a social network dataset. In this dataset, each graph corresponds to an online discussion thread, where nodes correspond to users, and there is an edge between two nodes if at least one of them responds to another's s comments. We choose graphs databases whose numbers of nodes are smaller than 50 and generate stationary graph signals from the generative model (1), with the graph filtering chosen as H (S) exp S. Our experiments include in total 871 testing graphs in the Protein database and 309 in Reddit database.
Infeasibility of rSpecT. For these real networks, we first check whether the infeasibility issue encountered by rSpecTT is significant. To this end, we adopt the random graph filters 4 + + I, where f,, = 1, 2, 3 are random variables following a Gaussian distribution with μ = 0 and o = 2. Then, we calculate the smallest 6, such that rSpecT is feasible. If the smallest ó, > O, rSpecT is likely to be encountered with the infcasibility issue. The results with different numbers of graph signals observed are shown in Table 1.
Table 1: Likelihood of infeasibility
Protein
Reddit
sample size
10
100
1000
10
100
1000
frequency
0,980
0.999
0.994
0.984
1
mean of δmin
38.514
26.012
14.813
1094.78
788
830.642
531.006
The experiment results present a decrease in the mean value of &min when more samples are used. However, a high frequency of the infeasibility issue occurring in these real datasets cannot be mitigated by the increasing sample size. These observations highlight the necessity of the careful treatment for choosing ó in rSpecT again.

---

Performance of different graph learning models. In this section, the stationary graph signals are generated by the low-pass-EXP ExP filter. We choose the famous statistical method called (thresholded) correlation [17] and the first GSP method that applies the log barrier to graph inference [15] as the baselines. The optimal threshold for the correlation method is selected from 0.1 to 0.6 and we search in (0.01, 0.1, 1, 10, 100, 1000) to obtain the best hyperparameter in Kalofolias' model. The parameter 8, in rLogSpecT is set / n n and in rSpec'T it is set as the smallest value that allows a feasible solution [31]. We also rely on the searching- based strategy to convert the learned weighted graphs from (r)SpecT & (r)).TSSeeTT.
The results of the ideal models with the true covariance matrix C.c applied are collected Figure We observe that on the real graphs, LogSpec'T achieves the best performance on average (median represented by the red line). Also, compared with SpecT, LogSpecT performs more stably. We remark that since the graph signals are not necessarily smooth, Kalofolias' model cannot provide guaranteed performance, especially on the Reddit networks.
Postenietaste
1.0
38
2.6
022
30
SpecT
Comvlation caiiiiio
Speet Commitiio canfl
137
327
100
10
100
16
1.2
100
102
109
Figure 4: Performance of LogSpecT and SpecTT Figure 5: Asymptotic behavior of rLogSpecT and rSpecT.
Figure 5 compares the performance of four different methods when different numbers of signals are observed." When the sample size increases, the models focusing on stationarity property can recover the graphs more accurately while correlation method and Kalofolias' method fail to exploit the signal information. This can also be inferred from the experiment results in Figure 4 since the models fail to achieve good performance from full information, let alone from partial information. The experiment also shows that when a fixed number samples are observed, the learned graphs from rLogSpecTT approximate the ground truth better than rSpeeTT. This further corroborates the superiority rLogSpecT on graph learning from stationary signals.
Conclusion
In this paper, we directly tackle the infeasibility issue encountered in graph learning from stationary graph signals [31] and propose the first condition that guarantees the model infeasibility. To overcome this, we propose an efficient alternative. The recovery guarantees of its robust formulation are analyzed with advanced optimization tools, which may find broader applications in learning tasks. Compared with current literature [23], these theoretical results require less stringent conditions. We also design an L-ADMM algorithm that allows for efficient implementation and theoretical convergence. Extensive experiments on both synthetic and real data are conducted in this paper. The results show that our proposed model can significantly outperform the existing ones.
We believe this work represents an important step beyond the fundamental model SpecT. Its general formulation allows for the transfer of the current extensions made on Spec'T. Testing our proposed models with these extensions is one future direction. Also, we notice that although the recovery guarantees for robust formulations are clear, the estimation performance analysis for the ideal-case models (i.e. SpecT and LogSpecT) is still incomplete. Investigating the exact recovery conditions is another future direction.
"The model in [34] is substitute of rSpecT to approximate SpecT and its hyperparameter is hard to determine. Hence, we omit the performance of that model.
10

---

Acknowledgments and Disclosure of Funding
This work was supported in part by the Hong Kong Research Grants Council (RGC) General Research Fund (GRF) Project CUHK 14203920.
References
[1] Daron Acemoglu, Asuman Ozdaglar, and Alireza Tahbaz-Salchi. Systemic risk and stability in financial networks. American Economic Review, 105(2))564-608, 2015,
[2] Álzzrt-L.ászíó Barabási and Réka Albert. Emergence scaling in random networks. science, 286(5439):509-517 1999.
[3] Stephen Boyd, Neal Parikh, Eric Chu, Borja Peleato, Jonathan Eckstein, al. Distributed opti- mization and statistical learning via the alternating direction method of multipliers. Foundations and Trends@ in Machine learning, 3(1): 1-22, 2011.
[4] Andrei Buciulea, Samuel Rey, Cristobal Cabrera, and Antonio G Marques. Network reconstruc- tion from graph-stationary signals with hidden variables. In Proceedings of the 53rd Asilomar Conference on Signals, Systems, and Computers, pages 56-60. IEEE, 2019.
[5] Andrei Buciulea, Samuel Rey, and Antonio G Marques. Learning graphs from smooth and graph-stationary signals with hidden variables. IEEE Trunsactions on Signal and Information Processing over Networks, 1-273-287, 2022.
[6] Samir Chowdhury and Tom Needham. Generalized spectral clustering via cromov-wasserstein learning. In Proceedings of the 24th International Conference on Artificial Intelligence and Statistics (AISTATS 202/), pages 712-720. PMLR, 2021.
[7] Jose Vinicius de Miranda Cardoso, Jiaxi Ying, and Daniel Palomar. Graphical models in heavy-tailed markets. Advances in Neural Information Processing Systems, 34:19989 20001, 2021.
[8] Steven Diamond and Stephen Boyd. Cvxpy: A python-embeddeed modeling language for convex optimization. The Journal of Machine Learning Research, 17(1):9009-9133, 2016.
[9] Paul Erdős, Alfréd Rényi, et al. On the evolution of random graphs. Publication of the Mathematical Institute of the Hungarian Academy of Sciences, 5(1):17-60, 1960.
[10] Benjamin Girault. Stationary graph signals using an isometric graph translation. In Proceedings of the 23rd European Signal Processing Conference (EUSIPCO 2015), pages 1516-1520, IEEE, 2015.
[11] Yiran He and Hoi- To Wai. Detecting central nodes from low-rank excited graph signals via structured factor analysis. IEEE Transactions Signal Processing, 70:2416-244 22330, 2022.
[12] Alan Hoffman. On approximate solutions of systems of linear inequalities. Journal of Research of the National Bureau of Standards, 49(4):263-265, 1952.
[13] Alexander Jung, Alfred O Hero III, Alexandru Cristian Mara, Saeed Jahromi, Ayelet Heimowitz, and Yonina C Eldar. Semi-supervised learning in network-structured data via total variation minimization. IEEE Transactions on Signal Processing, 67(24):6266-626,, 2019.
[14] Alexander Jung and Nguyen Tran. Localized lincar regression in networked data. IEEE Signal Processing Letters, 26(7): 1090-1094, 2019.
[15] Vassilis Kalofolias. How to learn graph from smooth signals. In Proceedings of the 19th International Conference on Artificial Intelligence and Statistics (AISTATS 2016), pages 920- 929. PMLR, 2016.
[16] Vassilis Kalofolias, Andreas Loukas, Dorina Thanou, and Pascal Frossard. Learning time varying graphs. In Proceedings 2017 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP 2017), pages 2826-2830 IEEE, 2017.
11

---

[17] Eric D. Kolaczyk. Statistical Analysis of Nenwork Data: Methods and Models. Springer- Verlag, New York, NY, USA, 2009.
[18] D. Christopher Manning, Hinrich Schütze, and Prabhakar Raghavan, Introduction Information Retrieval. Cambridge University Press, 2008.
[19] Hermina Petric Maretic and Pascal Frossard. Graph laplacian mixture model. IEEE Transactions on Signal and Information Processing over Networks, 6:261-270, 2020.
[20] Antonio G Marques, Santiago Segarra, Geert Leus, and Alejandro Ribeiro. Stationary graph processes and spectral estimation. IEEE Trunsactions on Signal Processing, 65(22):591155926, 2017.
[21] Gautier Marti, Frank Nielsen, Mikolaj Binkowski, and Philippe Donnat. A review of two decades of correlations, hierarchies, networks and clustering in financial markets. Progress in Information Geometry, pages 245-274, 2021.
[22] Gonzalo Mateos, Santiago Segarra, Antonio G Marques, and Alejandro Ribeiro. Connecting the dots: Identifying network structure via graph signal processing. IEEE Signal Processing Magazine, 36(3):16 2019.
[23] Madeline Navarro, Yuhao Wang, Antonio G Marques, Caroline Uhler, and Santiago Segarra. Joint inference of multiple graphs from matrix polynomials. Journal of Machine Learning Research, 23(76):1-35, 2022.
[24] Bastien Pasdeloup, Vincent Gripon, Grégoire Mercier, Dominique Pastor, and Michael G Rabbat. Characterization and inference of graph diffusion processes from observations of stationary signals. IEEE transactions Signal and Information Processing over Networks, 4(3):481-496 2017.
[25] Nathanaël Perraudin and Pierre Vandergheynst. Stationary signal processing on graphs. IEEE Trunsuctions on Signal Processing, 65(13):3462-3477 2017.
[26] Raksha Ramakrishna, Hoi-To Wai, and Anna Scaglione. A user guide to low-pass graph signal processing and its applications: Tools and applications. IEEE Signal Processing Magazine, 37(6):74 85, 2020.
[27] Samuel Rey, Andrei Buciulea, Madeline Navarro, Santiago Segarra, and Antonio G Marques. Joint inference of multiple graphs with hidden variables from stationary graph signals. In Proceedings 2022 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP 2022). pages 5817-5821. IEEE, 2022.
[28] Tyrrell Rockafellar and Roger J-B Wets. Variational analysis, volume 317. Springer Science & Business Media, 2009.
[29] Johannes O Royset and Roger JB Wets. An Optimization Primer. Springer, 2021.
[30] Santiago Segarra, Sundeep Prabhakar Chepuri, Antonio G Marques, and Geert Leus. Statistical graph signal processing: Stationarity and spectral estimation. Cooperative and Graph Signal Processing, pages 325-347 2018.
[31] Santiago Segarra, Antonio G Marques, Gonzalo Mateos, and Alejandro Ribeiro. Network topology inference from spectral templates. IEEE Transactions Signal and Information Processing over Networks, 3(3:4677-483, 2017.
[32] Santiago Segarra, Yuhao Wang, Caroline Uhler, and Antonio Marques. Joint inference of networks from stationary graph signals. In Proceedings of the 51st Asilomar Conference Signals, Systems, and Computers, pages 975-979. IEEE, 2017.
[33] Rasoul Shafipour, Abolfazl Hashemi, Gonzalo Matcos, and Haris Vikalo. Online topology inference from streaming stationary graph signals. In 2019 IEEE Data Science Workshop (DSW), pages 140-144. IEEE, 2019.
[34] Rasoul Shafipour and Gonzalo Mateos. Online topology inference from streaming stationary graph signals with partial connectivity information. Algorithms, 13(9):228, 2020.
12

---

[35] Rasoul Shafipour, Santiago Segarra, Antonio G Marques, and Gonzalo Mateos. Network topology inference from non-stationary graph signals. In 2017 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), pages 5870-5874. IEEE, 2017.
[36] Stephen M Smith, Karla L Miller, Gholamreza Khiiii-Koossiidi,, Matthew Webster, Christian F Beckmann, Thomas E Nichols, Joseph D Ramsey, and Mark W Woolrich. Network modelling methods for fmnri. Neuroimage, 54(2):875-89) 1.
[37] Oliver Stegle, Sarah A Teichmann, and John C Marioni. Computational and analytical challenges in single-cell transcriptomics. Nature Reviews Genetics, 16(3):133-145, 2015.
[38] Yuichi Tanaka and Yonina C Eldar. Generalized sampling on graphs with subspace and smoothness priors. IEEE Transactions on Signal Processing, 68:2272-2286, 2020.
[39] Yuichi Tanaka, Yonina C Fldar, Antonio Ortega, and Gene Cheung. Sampling signals on graphs: From theory to applications. IEEE Signal Processing Magazine, 37(6):14-33) 2020.
[40] Roman Vershynin. How close is the sample covariance matrix to the actual covariance matrix? Journal of Theoretical Probability, 25(3):655 686, 2012.
[41] Roman Vershynin. High-dimensional probability: An introduction with applications in data science, volume 47. Cambridge University Press, 2018
[42] Xiaolu Wang, Chaorui Yao, Haoyu Lei, and Anthony Man-Cho So. An efficient alternating direction method for graph learning from smooth signals. In Proceedings of 2021 IEEE International Conference on Acousties, Speech and Signal Processing (ICASSP 2021), pages -384. IEEE, 2021.
[43] Chenyue Zhang, Yiran He, and Hoi-To Wai. Product graph learning from multi-attribute graph signals with inter-layer coupling. In Proceedings of 2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP 2023), pages 1-5. IEEE, 2023.
[44] Hui Zhang, Ming Yan, and Wotao Yin. One condition for solution uniqueness and robustness of both 11-synthesis and ll -analllssss minimizations. Advances in Computational Mathematics, 3381::138 2016.
[45] Xiaoqun Zhang, Martin Burger, and Stanley Osher. unified primal-dual algorithm framework based on bregman iteration. Journal of Scientific Computing, 46(1):20-46, 2011.
[46] Licheng Zhao, Yiwei Wang, Sandeep Kumar, and Daniel P Palomar. Optimization algorithms for graph laplacian estimation via admm and mm. IEEE Transactions on Signal Processing. 67(16):4231-4444 2019.
[47] Shuai Zheng, Zhenfeng Zhu, Zhizhe Liu, Zhenyu Guo, Yang Liu, Yuchen Yang, and Yao Zhao. Multi-modal graph learning for disease prediction. IEEE Transactions on Medical Imaging. 41(9):2207-2211 2022.
13

---

(2: + x + y S 0, by + 2 0, x+, σ- - 0. A Organization of the Appendix
The appendix includes the missing proofs, detailed discussions of some argument in the main body and more numerical experiments, We organize the appendix as follows:
The proof of infeasibility condition (Theorem 3.2) is provided in Section B.
Explanations on conditions derived Theorem 3.2 are included in Section C.
The proof of properties the proposed model Lr)geT (Proposition 3.4 & 3.6) given in Section D and some additional properties discussed.
The truncated Hausdorff' distance based proof details of Theorem 4.1 and Corollary 4.4 are given in Section E.
Details L-ADMM and convergence analysis are in Section
Additional experiments and discussions on synthetic and real data are included Section G.
B Proof of Theorem 3.2
Since the linear system (4) has no solution, know from Farkas' lemma that the following system has solutions:
In
B A < 0 (m-1)x
(11)
0 (1m-1)mm-)) x(m-1)
(mm22)
B A < 0 (m-1)mm-21 x1'
Let xx E Rm² be solution to (11). Denote 24 = maxfx", 00. x = max xx" 00. Then, there
exists c E (0, 1] such that
B A (x+ x-- + c1 m1 (x4 x=_ xz) (1mmmmm1:: (m--))m-n) S.
Define == -T
m² (x+ x-),
cl J1 + +. and set 8 = c. For all d E [0, o). (x+ x=, ,3
is solution to the following linear system:
B A (++ - x-) + 2 11 1:m 0 (m 1)(sm 2) S 0,
Again, from Farkas' lemma, this implies that the following linear system does not have a solution:
A, Bs f1, 2 0,
A, Bs
t1,, < 0,
$
(12)
1m-1
Dim
$ 1,
where s pm(mm11/2 and E E. Since (12) is equivalent to:
CC..
SCr 000,00 5 6,
(S1); = 1,
(13)
S S
the above argument indicates that (13) does not have solution. Suppose rSpecT has feasible solution S", then
|| S'C S |C,.S' S'C,llr 5..
Hence, S' also a solution (13). However, (13) does not have solution. We can conclude that :SpecT is infeasible in this case.
14

---

C Explanations on Sufficient Conditions in Theorem 3.2
We elaborate more on the infeasibility condition thatA, has full column rank. An application of the condition is Example 3.1. Specifically, we know that in this case,
h12
h22
h12
h22
h11
0
-112
B
and An =
1
h12
0
h11
h22
h12
h12
h12
0
This implies that
h22
h11
A,
h1
h22
Hence, when hu h22. A, has full column rank. This means that when d is small enough (from Example 3.1 we know 8 /2 h11 h22), the model rSpecT is infeasible.
D Proofs of Properties of Lr)Leg
D.1 Proof of Proposition 3.4
Since the constraint set S is cone, it follows that for all y 0, 78 = S. Then, we know that
Opt(C,a)
argmin SES.CS=SC
||S||,,1 - a1 (S)
Y
argmin 1555.C78
|| 1, a1 1og(?S1)
argmin se -C
>/| - a1 (og(
=Y argmin S
|| |||,,
1 log(S1)
= aC,
where the third equality is from the basic calculus rule the logarithm function. Set ? a and then a)) a Opt(C, 1), which completes the proof.
D.2 Proof of Proposition 3.6
The proof will be conducted by constructing feasible solution rLogST Recall that A, & C, C @ I and the matrix B E Rmm xm(m- 1)/2 that maps a non-negative vector to the vectorization of valid adjacency matrix. Let S = min As, 1}-mat( Bs) with s (2mm21)m/2 being non-negative vector, where mat(-) is the matricization operator. Note that
vec C (CSSSSSC.) @ Cn C. I) vec(S) Ancec(
Then, we know that
|C,S - SC,||p - min | Bs
|IA, Bsll2 5.
Thus, the given S is feasible solution for roogS and it completes the proof.
D.3 Properties of optimal solutions and values of (r)LogSPecT
In this section, we further discuss some properties of the optimal solutions/value of proposed models, which are useful for deriving the recovery guarantee. More specifically, we obtain an upper bound on the optimal solutions (which may not be unique) independent of the sample size n and the inaccuracy parameter 6,, Also, lower bound of optimal values follows.
15

---

Proposition D.1. The following statements hold:
For an optimal solution S (resp. St)) to LogSpecI (resp. rLogSpecT with any given sample size n), it follows that
11,1 am and || S am, Vop 0.
òn 2 2amm|Cn Coo then
amm(1 log S S Vn,
where f" (resp. fi)) denotes the optimal value of LogSpecI (resp. rLogSpecT).
For the first statement, let consider the Karush- Karush-Kuhn- Kuhn-Tucker (KKT) conditions of LogSpecT and rLogSpecTT. Since the LogSpecTT is convex problem and Slater's s condition holds, the KKT conditions are necessary and sufficient for the optimality, i.c., there exists (A1, A2) mm (Ns such
that
Vs(1S -a log(S" C.Ar C. + A 0, - = =S
(14)
ES
where (s((S)) (N E Rm supxes (X - S" N) s 0} is the normal cone of S at S", and 11,, is well-defined since || 111.1 (, 11 at 5 N 0. which is differentiable. Taking further calculation gives that
=11 (Vsl log 1));; = (S-1);
Combining this with (14) by taking inner product of both sides with S we obtain that
( his (A, S S'C + (A2, = 0. (S+1),
ES 10 a >
(15)
From the structure of S and the fact that A2 E NN*SS), one has that (A2, S*) = 0. Also, note that S C Hence, the equation (15) can simplified as the desired result:
(S Jy (5:1);
m
m
(S*);
1S (S)) a 1.j
am.
S*'));
The KKT conditions of rLogSpecT indicate that there exist \1 > 0. Ag E (s and Q
0|| C - -S (i.c., the subgradient of the function - C, SC,,IF at S such that
[s([ - a1 log (S, + Q A = A CI,, Cnlll 5n) = 0, S E.
(16)
Moreover, from the definition of the convex subdifferential we know that 0 > C C.C.lFF (Q, St)) Thus, after taking inner product of both sides of the equation (16) with S it follows that:
Σ (Sh): - om + X (Q, S: + S")
1.j
Σ
(5)) - am + A1 S S Cnlle + (A2,
i.j
Σ (8)) + X16n.
1.
which implics that (S, S am 21dn S am. This completes the proof of the first statement.
16

---

For the second statement, we first prove that in and v* are larger than amm11 - log Define the auxiliary function g : R such that g(x) == log x for any z E +. whose minimum attained at a. Since for any S E S,
1
TR
f(S) = Σ9 9
Σ
Sig
where f is the objective in LogSpecT, it follows that
f(S) Σ g(a) am(1 loga). a).
This implies that o, and v* are larger than amm(1 - log a).. Next, we will show 2 S v². Consider any optimal solution S* LogSpecT. We show that is feasible for rLogSpecT.
C,,''* - S" C. = , - S + S'CCs - Sr
2 1.1 C Coo ll S 2aml|C, - 8n.
where the equality comes from Cx S' = C the first inequality comes from the fact that |XY|| < |X|||| < 1 Y |. the second one comes from the first statement and the last
one is due to on 2 2am C, C Hence, S is feasible for rLogSpec'T, which indicates that i I v². The proof is completed.
E
Proof of Theorem 4.1 & Corollary 4.4
E.1 Truncated Hausdorff distance
In this section, introduce advanced technique in optimization that is efficient in analyzing the recovery guarantee robust formulations. Before that, we introduce the concept of truncated Hausdorff distance between sets.
Definition E.1 (Truncated Hausdorff Distance [29, 6.J]). For any 2 0, the truncated Hausdorff distance between two sets and D is defined as:
(C, D) max BBO, p), D), dD ( (O.,
turns out that the distance between the optimum two minimization problems can bounded with the truncated Hausdorff distance of the epigraphs under some conditions. The result is captured the following lemma.
Lemma E.2 ([29, Theorem 6.56)). Let p E [0, oo). Suppose that the sicended-real-vlleee lower semicontinuous functions f.g: R R satisfy
inf inf g E [-p.
argmin f, arggii C B(0.p).
Then, it follows that
inf f inf yl S dp(epi J, epipig) g).²
(17)
Suppose further that 2 2 (epi (l epi then one has
dist(x" erammin f) 1
(18)
where e-argmin is the E-suboptimal solution set of that is defined as e-argmin f = {ar dom f f(x) & inf c), and x; is minimizer of g.
From the above lemma, we know that two optimization problems are close enough (in the sense of truncated Hausdorff distance), then the optimum of them should be close to each other. Hence, in order to apply this result, we need bound the truncated Hausdorff distance in an explicit way. which is solved by the following Kenmochi condition.
For function f R² - I. its epigraph is defined epi ((x, v) y f(x))
17

---

Lemma E.3 (Kenmochi Condition [29, Proposition 6.58)). Let E [0, oo). Then, for f.g: R R with nonempty epigraphs, one has that
inf 9 S max n. Vo |f s B(a,1)
cppe inf 17 0:
inf max(g(x)), + + Vø E R(m,))
where [f < = (x E R²: f(x) S p}.
E.2 Proof of Theorem 4.1
Before presenting the proof, we first introduce the following lemma.
Lemma E.4 (Hoffman's 's Error Bound [12]). Consider the set S S= E R" Ax S b). There
exists C > 0 such that for any a E R", one has
dist(a,, S) C. |1(Ax - b)+
For the sake of brevity, we denote
(S) ||S|| 1,1 a1 (og((S))))))) in SC,ll 6n) (s()): ((S) = ||||| - a1 log(S (S) + 140 al SCx [|)) + (s)).
Hence, the optimization problem LogSpecT (resp. rl.ogSpecT) equivalent to inf ʃ (resp. inf J.).
Now, we aim to use Lemma E.3 to bound d f, epi fn). Let S S (() satisfy
F(S) S and SC = C.
Then, we know that and consequently S is in the domain of Jn. Then, it follows that for any > 0, we have
||CC, - Sl & 2 I - C & 2pl C S 8n,
inf Jn < (S) = f(S)= < max( -p}, VS E [ S p] n B(0,p) p). R(S.))
(19)
Before verifying the reverse side of the Kenmochi condition, we first consider the non-emptiness of [f S p) BB, p). Since
8n 2 Cp - C ll 2 2ami|IC, C.oll:
it follows from Proposition D.I that || & am S p and f* S S p, which implies that [f S /) (BO p) nonempty. Let S,, E [f S p) ()BB( p). Then, one has that
S, ES and IC.S. S, Cn S 6n.
Hence, it follows that
S, S.C S
Call CCSS, - S.Cullr & Cp - Cull + dn-
Also, note that there exists 3 > 0 such that (S, : 2 3 for all E [m] as Jn < p and ||S,||,,1 ol T log(S)) 1) + oc when S, - O. Thus, applying Lemma E.4 to the linear system
5 = S E mmm. SCx Co S S, (S1); 8, Vie [m]}
yields that there exists such that
dist S,, 5) à. C2p| Cn +
Hence, there exists š in the domain of / such that
|| 51 < ē. (2p) C.ll + 5,) and (51), 2 B, Vi E [m].
18

---

Since the function S + ||S||1,1 - al" (og(S1) is locally Lipschitz continuous when (S1); 2 3,
there exists L > 0 such that
f($) || -a1 log(S1) |S a1 log(S,,1) + L|| S11e
(.S.)) L| $11F
s (.,.)) Le (2pl 00 -C + 8,).
Setting a 2 maxx11, L} i, one can obtain that for any S, E fS S pl n (B)
inf f f(S) - f. (S.) t (2pl1C C..ll + 5.)) S max ( f (S.): -p} + + B(Sn.T)
(20)
where 1) 1-== c C2 C. + 5n). Combining inequality (19) and (20), we can conclude that de(epi f, epi fa) S ci (2 CClll + 6n).
(21)
In order to derive the conclusion (i) and (ii), it remains to check the requirements in Lemma F.2. Since 2 om, the first statement of Proposition D.1 shows that the optimal solutions inf f and inf fr lie in (O, p). Since p 2 f and -p I amm( [ log a), the second statement of the proposition shows that inf inf fn [-p, p]. Hence, applying Lemma E.2 completes the proof of the first two statements.
To prove conclusion (iii), first make the following two claims:
(a) Si is a singleton, whose element is denoted by 1"1,
(b) For any [0, o0), there exists a (()) > such that for all S S and S. E St one has that
||S 1 1| S 8(5) VE.
(22)
Granting these and with the help of Theorem 4.1 (ii), we can derive that for all S E S
dist 1, = || - 1'1 < | S2e., 1112 + ||5 1 - --5'11 1 ||2
S m dist SS, S2c.) + S 1- - 101112
I Ōien + 82 /En.
where 21. 82 are positive constants, and S2e. E Stc. satisfies || - S2e, || = dist( Sts.))
(whose existence is guaranteed since St is convex and compact). Hence,
1i - ŌIEn + ē2VEn.
To proceed, it remains to prove the claims. Define an auxiliary function h: : R + R as h(z) =
Em I3- > log zi for each RTT. Consider the following optimization problem:
min h(x)
(23)
I D EE1SIE E R E
S that is feasible for LogSpecT).
For the sake of brevity, denote the e-suboptimal solution set of (23) as H. In the remaining part, we will first show that S 1 H5 and then, by the strict convexity of h, the desired two claims hold.
The first step is to show that the optimal function value of the problem (23) satisfies h" = f Since it is obvious that a S' 1 is feasible for (23), ht < h(()) f(S*) f*. Suppose to the contrary that h²* < f", from the fact that the objective function is coercive and continuous and the feasible set closed, there exists S such that it is feasible for LogSpecT and x= 51, where xt is an optimal solution to (23). Since ht* h((=') = h(S1 = f(S), this contradicts the fact that f($) > f*. Hence, ht = f² Next, we will show that S H Consider any e-suboptimal solution S St i...
M(1 ((S) -- +e.
Hence, S1 E H and it implies that S C Hi On the other hand, for any e-suboptimal solution He there exists S that feasible for LogSpecT such that S1. Thus,
f(S) h(x)) < h" = +
19

---

This implies that S E S: and consequently H S S Hence, H = S
Since h is strictly convex, its optimal solution set Hj is singleton. Then, Sj H6 is singleton, which proves the first claim. For the second claim, we know that for any S. E S there exists 2t E Hg such that
|S - 101 = - - x
(24)
where x' E Hg. The cocrciveness asserts that 2g and x" are bounded. This together with the fact that h is strongly convex on any bounded set, illustrates that there exists μ >> 0 such that
h(xz) 2 h(x*)) + (Vh(x"), Te x²) + 1| - x 113 h(x²) + | x ==12.
(25)
where the second inequality comes from the global optimality x Combining (24) and (25) gives that
|S S01 = ||=, - x=112 v((((xz) h(x*)) S vue.
This completes the proof of the claims.
E.3 Proof of Corollary 4.4
Suppose to the contrary that there exists a sequence (S where the nihh element an optimal solution to rogSSc'' with sample size such that
dis A
From Proposition D.I, we know that {Sa)n is bounded, and consequently, has a convergent subse- quence. Without loss of generality, may assume that the sequence itself is convergent and the limiting point is S", Note that
- S Sn, C₂ Co and on -
This indicates that S is feasible for LogSpecT.. Then, from Theorem 4.1, we know that ((S))= = f", which leads to f(S") f since f(-) = ||- ||,,1 al log(-1) is continuous. Together with the fact that S" is feasible, we conclude that S is optimal solution to LogSpecT. This further implies that dist(S";, - O, which is contradiction.
Hence, S"
E.4 Proof of Lemma 4.7
Recall the generative model (1). Since follows a sub-Gaussian distribution, it can be shown that for every > 0,
F([|0|12> <
|10022
t [[P(())
S v''
for some positive constant v', which means that æ also follows sub Gaussian distribution. Thus, due to the uub-Gaussian property, |, Cooll can be explictly bounded by the following lemma.
Lemma E.5 ([40, Proposition 2. 2)..). Consider sub-Gaussian, identical, independent random vectors 1, 22,..... Fy E I with n > m. Then for all it follows that
P
Σ
a
E[ax T
E
20 2mm-ss
for some constant
Setting e² = log(n)m) Lemma E.5 indicates that with high probability (lower bounded by 1 n 1),
| - C S (v log n
F Derivations of L-ADMM and Convergence Analysis
This section includes the details of L-ADMM for rLogSpecT.
20

---

F11 Proof of Proposition 5.1
Note that the minimization problem (7) is separable for Z and 9, and can be split into two subprob- lems:
min |C - ((h) C + ( /0-kppp21ll- 203(0,6,,
(26)
min -al log q + 2 (K)T (q- ( 1) 2 - 5(h)
(27)
For problem (26), the optimal solution is the projection of C, S(,) - (k)) Cn + A (2) /p onto B(0, which is given by
8n with Ž C,,s((k (kk) C₁ + A A1) 1p
z(k+l)=
min
I,
||
For problem (27), the first-order optimality condition gives
-a1 / 2/ p (*)) 1) 0.
This together with the fact that the objective function convex implies that
q (x++1)
q" + 4a/p1
with
1 (ps 1 (1)
2
p
F.2
Calculation of (T)(
The projection of X to S can be calculated via an optimization problem:
min
S1? - S
|X
s.t.
S 5
S = 0, i= =1,2.... m,
Sú 0, Vi,j.
which is equivalent to
Σ ((Xg - SSy) + (Xyr Sty)
min
icj
Sú 0, Vi <j
S = 0, Vi.
Hence
(Is(X)); max (0, Xij + Xjih,
itJ, i= j.
F.3
Stopping criterion and updating rule of
We follow the procedures in [3] to update p in each iteration. Similarly, we define the primal residual and dual residual follows:
Pres 1) V ( S S(x+1) C + [|q(++1)
1
d(b++ (0) (((tt))) -
(5(++1) S(()) ( 1 ( - 5(1)
The aim of updating is to control the decaying speed of Pres and dres such that their difference is too large. To this end, we update p adaptively following the scheme:
(k)
if p(+tt) > (R++1.
(k++1)
( /2,
(k+ (k+1)
(k))
otherwise.
When Pres and dres are both smaller than the threshold 10-3, we stop the algorithm.
21

---

F.4
Convergence analysis
Define D = Diag( 1 m=:
m E mmxm" Then, D satisfies (ve S1 and | D = m.
Denote
Q -I - -= A An. Then the linearized ADMM update (8) of S can be written
min L(S) + 2 [[vec(()) vec( ( (()))
where |0ilq == x' Qu. Since T > m+I An 11 we know that Q is positively definite. Consequenty, by treating (Z, q) as one variable, we can apply Theorem 4.2 [45] and directly obtain the result.
G More Experiments and Discussions Synthetic and Real Data
G.1 Influence of graph filters
To make a fair comparison between rSpecT and rLogSpeeT. we test rSpecT on BA graphs with the same graph filters and the results are reported in Figure 6. is obvious that rSpecT fails these cases and cannot benefit from the increase in sample size. This is reasonable since SpecT fails on BA graphs as indicated in Figure 3, let alone the approximation formulation rSpecTT.
rSpecT Performance
rLogSpecT Performance
0.7
lowpass-EXP
1.0
lowpass-EN
highpass-EXP
highpass-EXP
0.6
QUA
0.8
QUA
F-measure
0.5
E-measure
0.7
0.4
0.6
0.5
0.3
0.4
0.2
0.3
0.2
10a
10a
10a
10a
10a
10
10a
102
10a
10a
105
Sample Size
Sample Size
Figure 6: Performance of rSpecT on BA graphs. Figure 7: rLogSPecT on ER graphs with 8n 20yloggn/n n/n.
We further test rLogSp on ER graphs with different numbers of signals observed. The parameter 8. is set as 20 /log n/n and the results are reported in Figure 7. The figure shows that for graph filters that are not high-pass, rLeeS can achieve nearly perfect recovery when the sample size large enough. Also, compared with the performance on BA graphs, rLogSpecT works better on ER graphs. This observation better on ER graphs than BA ones. We further notice that the difference between the low-pass graph filter and the high-pass one huge. To check the conjecture that rLogSpecT generally performs better on low-pass graph filters, we choose different graph filters exp(1.S) with / ranging from -2 to and conduct the experiments on graphs. When the graph shifting operator is the adjacency matrix, the positive low-pass parameter / corresponds to low-pas: graph filters and the negative I corresponds to the high-pass ones [26. 11]. We omit the case when contain any graph information (note that exp(OS)) I).
accordance with the conclusion from Figure 3 that LogSpecT performs
since this filter does not
We then repeat the experiments for 50 times and report the average results in Figure 8. The comparison between the performance of low-pass graph filters and high- pass graph filters indicates that the low- pass graph filters generally outperforms the high-pass ones. A closer look the results shows that the performance grows faster when the absolute value of t is smaller. And eventually, the graph filter with smaller absolute value prevails. This observation interesting since Figure 3 indicates that the choice of graph filters has few impacts on the model performance. One explanation is that both low-pass graph filters and high-pass graph filters attenuate some frequencies of the graph and the larger absolute value of t leads to the more loss of information carried by finite signals.
22

---

Inmpass matermeter
loa pass parameter a-1 Iow pass parameter -05
lowpass persmmeter owipase perameter =1
lowpass parameter = 1:5 low posa parameter -
wwraas safommicte" 1.5
0.9
0.8
0.7
Average
0.5
0.3
0.2
107
16
10a
10
10
10
Sample Size
Figure 8: Effect of Low-Pass Parameter: different performance of graph filters exp(tS) with t
ranging from
G.2 Experiments USPS dataset
USPS dataset is handwritten digit dataset. As shown in [25], it is nearly stationary with respect to the nearest neighbour graph. This dataset collects images of handwritten digits. In the experiment, each pixel is viewed as a node and the value on the pixel is view as the graph signal. We follow the approach in [25] construct the 20 nearest neighbour graph, on which the data is verified to be nearly stationary. More specifically, we pick the 1296 digit images. The weights between two nodes are decided by the Gaussian radial basis function. It turns out that the stationarity measure U, ", equals 0.78 in this dataset. Here U is the eigen-matrix of
:= || diia(U C,U))///I the constructed graph and C, is the covariance matrix of the observed graph signals. We view the 20 nearest neighbour graph as the groundtruth and compare the learned graph from rLogSpecT with it. We present the subgraphs consisting of the top 10 nodes with the largest degrees in the two graphs respectively in Figure 9, The F-measure of these two subgraphs is 0.96. This result corroborates the efficiency of our proposed rLogSpecT model.
graph.10 edges
1.00
0.79
0.75
030
3.0
-3.5
-0.25
-150
-0.0-
-0.25
0.02
1.00
-2.00
-0.25
1.00
Figure 9: Graph Learning on USPS dataset. (Left: subgraph of the symmetric 20 nearest neighbour graph. Right: subgraph of the learned graph from rLooSP The F-measure is 0.96)
23


Reviews:
[{"review_id": "2vfqmcaIPk", "reviewer": "Reviewer_SR4h", "summary": "This is an emergency review, regrettably, the paper is outside of my expertise.\n\nReview:\nThis is a theoretical work concerned with graph learning from signals. The work identifies an issue with the feasibility of a previous well known model in the field, SpecT, then introduces a novel model, LogSpecT, for which this issue does not occur.\n\nWhile I regrettably could not comprehend the main content, the topic, general presentation, and references make me believe that this paper is topically not suited for neurIPS. I concede that the CfP of neurIPS is somewhat vague and broad, but based on the common topics in the field, I believe for the conference this is very much a niche topic to which the vast majority of attendees would not be able to connect.\n\nSuperficial observations\n - There isn't any mention of the word \"neural\"\n - From 45 references, the vast majority is in signal processing. A single one is at neurIPS\n - The intro is exceedingly dry to read, with no motivating example, nor any illustration. Especially, what are examples of signals considered for graph learning? What would the stationarity assumption mean in these examples?\n\nI understand that reading a superficial review like this may be frustrating to the authors, but can do no better. I would much suggest the authors consider a specialist forum, where this type of content might be better understood and appreciated.", "questions": "see above", "limitations": "see above", "rating": 3, "confidence": 1, "soundness": 2, "presentation": 2, "contribution": 2, "strengths": "see above", "weaknesses": "see above", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thanks for your comments. Here are some illustrations based on your questions/suggestions.\n\n**Q1:** There isn't any mention of the word \"neural\"\n\n**Answer:** Graph learning from signals studied in this paper has much relevance in the machine learning and neural science community.\nFor example, it has broad applications in bioinformatics [RX$^+$21].\nMore specifically, we consider the task of graph topology inference from observed data. It has been applied to recognize the functional connectivity related to neural activity patterns. The vast references in this field include [DE06, EO09, GX$^+$21], to name just a few. Due to limits of page length, we are only able to focus on the modeling issue and its theoretical guarantees in this paper. \n\n**Q2:** From 45 references, the vast majority is in signal processing. A single one is at NeurIPS\n\n**Answer:** We are working on a type of data that was newly proposed in the signal processing community and it has not gained much attention from the machine learning community. However, the task we were caring about has long been a focus in machine learning. For example, [K16, KP18] studied how to learn a large graph from smooth graph signals, [YCP20, YCP21, VYP22] investigated this task from the probabilistic graphical model.\n\n**Q3:** The intro is exceedingly dry to read, with no motivating example, nor any illustration. Especially, what are examples of signals considered for graph learning? What would the stationarity assumption mean in these examples?\n\n**Answer:** Thank you for your comments in the introduction part. In practice, the meaning of signals and graphs depends on the scenarios and the models. The signals can be understood as the nodal features in the graph learning task. For instance, in the financial market, the companies can be modeled as the nodes and the signals are the stock returns of each company. The graph learning task is then to investigate the relationship between the companies' stock returns. The relevance of the stationary assumption in several real data sets has been revealed in several existing literature.  [G15] showed that a part of the weather data can be explained by stationary graph signals. [PV17] showed that the well-known USPS dataset and the CMUPIE set of cropped faces are close to stationary. [SM$^+$17] compared the performance of learning protein structures with/without the stationary assumption and validated this assumption. [ZHW22] applied the stationary assumption to US Senate Roll Calls and achieved good performance.  In a nutshell, the stationarity assumption is a theoretical property that stems from the stationary time series and has been verified in several real data. This property has to be understood more from a theoretical viewpoint instead of the intuition.\n\n\n[RX$^+$21] Rui, L., et al. Graph signal processing, graph neural network and graph learning on biological data: a systematic review. IEEE Reviews in Biomedical Engineering, 2021.\n\n[DE06] Danielle, B. and ED, B. Small-world brain networks. The neuroscientist, 12(6):512–523, 2006.\n\n[EO09] ED, B. and Olaf, S. Complex brain networks: graph theoretical analysis of structural and functional systems. Nature reviews neuroscience, 10(3):186–198, 2009.\n\n[GX$^+$21] Gao, S., et al. Smooth graph learning for functional connectivity estimation. NeuroImage, 239:118289, 2021.\n\n[K16] Kalofolias, V. How to learn a graph from smooth signals. In Proceedings of the 19th International Conference on Artificial Intelligence and Statistics (AISTATS 2016), pages 920–929. PMLR, 2016.\n\n[KP18] Kalofolias, V. and Perraudin, N. Large scale graph learning from smooth signals. In International Conference on Learning Representations, 2018.\n\n[YCP20] Ying, J., et al. Nonconvex sparse graph learning under Laplacian constrained graphical model. Advances in Neural Information Processing Systems, 33:7101–7113, 2020.\n\n[YCP21] Ying, J., et al. Minimax estimation of Laplacian constrained precision matrices. In International Conference on Artificial Intelligence and Statistics, pages 3736–3744. PMLR, 2021.\n\n[VYP22] Vinícius, J., et al. Learning bipartite graphs: Heavy tails and multiple components. Advances in Neural Information Processing Systems, 35:14044–14057, 2022.\n\n[G15] Girault, B. Stationary graph signals using an isometric graph translation. In Proceedings of the 23rd European Signal Processing Conference (EUSIPCO 2015), pages 1516–1520. IEEE, 2015.\n\n[PV17] Perraudin, N and Vandergheynst, P. Stationary signal processing on graphs. IEEE Transactions on Signal Processing, 65(13):3462–3477, 2017.\n\n[SM$^+$17] Segarra, S., et al. Network topology inference from spectral templates. IEEE Transactions on Signal and Information Processing over Networks, 3(3):467–483, 2017.\n\n[ZHW22] Zhang, C., et al. Product graph learning from multi-attribute graph signals with inter-layer coupling. arXiv preprint arXiv:2211.00909, 2022."}, {"review_id": "YH0LcB6PX5", "reviewer": "Reviewer_xuW7", "summary": "In this paper, the authors improve upon an existing method, rSpecT, to impute graphs from stationary signals. The latter frames the recovery of the (weighted) adjacency matrix as a convex optimization problem with constraints —- namely the commutativity of the covariance matrix and the adjacency matrix, as well as the first row sum to be equal to 1, in order to avoid the trivial 0 solution. Here, the authors show that this yields infeasible solutions in a number of cases. Consequently, the authors propose to replace the row sum constraint by a log barrier. They show that the solution is always feasible (however, it ceases to be unique).  They then show through a set of experiments that their solution outperforms the state of the art approaches in recovering the underlying graph.", "questions": "See weaknesses.", "limitations": "The method could have been more extensively evaluated. ", "rating": 6, "confidence": 1, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "First, a disclaimer: although familiar with GSP, I am not an expert in the graph recovery literature. My comments must thus be taken sparingly.\n\nOverall, the paper is clearly written and motivated.\nThe discussion on the infeasibility of the current state of the art method, rSpecT, was interesting. The solution that the authors bring is therefore well motivated.", "weaknesses": "As highlighted before, I am not an expert in the field and I am not able to put this paper in context of the general literature on the topic.\n\nA few questions however:\n- the authors claim that their method is more computationally efficient with the use of ADMM. Do they have plots showing the running time? ADMM, although fast at each iteration, can be slow to converge.\n- the experiments could have been expanded: in particular, the authors only present ER and BA graphs, but they could have varied the parameters to try and assess whether their method is robust across all sorts of topologies.\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your comments and advice. We hope the clarification can clear your concerns.\n\n**Q1:** The authors claim that their method is more computationally efficient with the use of ADMM. Do they have plots showing the running time? ADMM, although fast at each iteration, can be slow to converge.\n\n**Answer:** For the running efficiency of our proposed L-ADMM, we refer the reviewer to our global response and our response to Reviewer U8JP that talks about the complexity in each iteration. In the experiment, we show that L-ADMM can reach the accuracy of $10^{-6}$ in 3500 iterations for a graph with 100 nodes. Also, we can observe the local linear convergence for L-ADMM in our model, which indicates that our algorithm can run efficiently. We will add the detailed discussion in our paper.\n    \n**Q2:** The experiments could have been expanded: in particular, the authors only present ER and BA graphs, but they could have varied the parameters to try and assess whether their method is robust across all sorts of topologies. \n\n**Answer:** In order to show the robustness LogSpecT, we investigate its performance when the sparsity level of the underlying graphs changes. Specifically, we set the parameter in the generation of ER graphs from 0.1 to 0.5, which represents the probability that an edge exists between two nodes. The number of nodes is set as 50 and the graph filter is chosen as QUA: $H(S) = S^2 + S + I$. The average results of 10 independent replicates are shown in Figure 1 in the global response. This figure demonstrates the robustness of LogSpecT."}, {"review_id": "bCEfKexd6X", "reviewer": "Reviewer_xFRE", "summary": "This paper addresses the problem of learning a graph from stationary graph signals. To this end, the authors introduce two new optimisation problems called LogSpecT and rLogSpecT.  They prove some recovery guarantees on the optimal value of (r)LogSpecT and give a convergence rate when the signals are sub-Gaussian. They also provide an efficient solver based on L-ADMM that admits closed solutions at each substep.  Finally, some empirical results highlight the performance of their methods.", "questions": "1\\ What is the intuition behind the need for a delta greater than a certain value in Theorem 4.1?\n\n2\\ In the GSP literature, there are many more recent methods than the one of Kalofolias for learning the graph. Why choose to compare your methods only to the one of Kalofolias?\n\n3\\ Furthermore, why exclude this method from the comparison in the results of the section 6.3?\n\n4\\ Experiments are made with a maximum of 60 nodes. What is the maximum number of nodes that can be learned in a reasonable time? if possible, can we have some examples?\n\n=====\n\nMore open questions:\n\n5\\ Is it possible to link the optimization problem to a maximum likelihood estimate?\n\n6\\ It seems that the filter h only has an impact on the covariance estimate. Is it true? If the answer is \"no\", would it be a good idea (if possible) to include information about the h filter in the optimization problem?", "limitations": "-", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "1\\ The paper is very well written and easy to follow\n\n2\\ The problem is well-motivated and the solution seems theoretically sound\n\n3\\ The paper is solid and almost all aspects encountered in the graph learning problem are covered (optimization, recovery guarantees, convergence,...)", "weaknesses": "1\\ There are no error bars. They should be added to really show the results of each method and to be able to compare them.\n\n2\\ A table with all the measurements reported could be a good thing. For example, in all the experiments, only the F measure is given, but precision and recall are also important for assessing the quality of the learned graphs.\n\n3\\ An illustration comparing the real graph with the learned graph (for example by plotting the adjacency matrix) could be interesting.\n\n4\\ The complexity and the time until convergence are not given for the methods. An additional experiment on this aspect should be added.\n\n=======\nRemark\n\nthe step size tau is never defined", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thanks for your comments and advice. We hope that the clarification may help to clear your concerns.\n\n**Q1:** - There are no error bars. They should be added to really show the results of each method and to be able to compare them. - A table with all the measurements reported could be a good thing. For example, in all the experiments, only the F measure is given, but precision and recall are also important for assessing the quality of the learned graphs. - An illustration comparing the real graph with the learned graph (for example by plotting the adjacency matrix) could be interesting.\n    \n**Answer:** Thank you for your helpful advice on presenting the results. We showed the error bars in the real data experiments for fair comparison. However, we omit the error bars in the synthetic data experiments as both models achieve nearly perfect recovery or the difference between them is too obvious. We would like to present all the measurements in the appendix, together with an illustration comparing the groudtruth graph with the learned graph.\n    \n**Q2:** - The complexity and the time until convergence are not given for the methods. An additional experiment on this aspect should be added. - The step size tau is never defined.\n    \n**Answer:** For the running time and complexity of our proposed L-ADMM method, we refer the reviewer to our global response and our response to Reviewer U8JP, where we discussed the complexity in detail and added an experiment to present the residuals until convergence. Also, thank you for your kind reminder.  $\\tau$ in the paper is defined as the parameter in the linearization process and can be viewed as a step size in updating $S$. We will add these clarifications to our paper.\n    \n**Q3:** What is the intuition behind the need for a delta greater than a certain value in Theorem 4.1?\n    \n**Answer:** The intuition that $\\delta$ should be larger than a certain value is that $\\delta$ is to reflect our view on the error of estimating the covariance matrix. We can bear with the overestimate of the errors while the underestimate is unacceptable. From a theoretically viewpoint, the $\\delta$ should be set to ensure that the feasible set of rLogSpecT at least contains the optimal solutions to LogSpecT. This part is required in our proof. Although this condition appears to be stringent, our analysis in Corollary 4.8 then shows that if $\\delta$ is set to decrease in a proper rate, the condition will hold with high probability. \n\n**Q4:** In the GSP literature, there are many more recent methods than the one of Kalofolias for learning the graph. Why choose to compare your methods only to the one of Kalofolias?\n    \n**Answer:** We only compared Kalofolias' method in the experiments and omitted other models for the reason that Kalofolias' method is the only graph learning model with the log barrier. Since our model is another graph learning model with the log barrier, it is natural to compare our model with Kalofolias'. We didn't compare our model with the more recent models as they are based on the probabilistic graphical model instead of the stationarity assumption. However, for fair comparison, we present the experiment results of ALPE, NGL-SCAD, NGL-MCP on the Protein and Reddit datasets in the global response.\n    \n**Q5:** Furthermore, why exclude this method from the comparison in the results of the section 6.3?\n    \n**Answer:** We excluded the performance of Kalofolias' method in section 6.3 for the reason that SpecT is the only efficient model for graph learning from stationary signals and Kalofolias' method is designed for smooth graph signals instead of stationary graph signals. However, for fair comparison, we add the experiments of Kalofolias' method on 50-node ER graphs with different sparsity levels in the global response. The results show that this model cannot learn graphs from the stationary signals.\n\n**Q6:** Experiments are made with a maximum of 60 nodes. What is the maximum number of nodes that can be learned in a reasonable time? if possible, can we have some examples?\n\n**Answer:** The experiment results on 100 node graphs are shown in the global response, which will take around 20 seconds to learn a graph. We note that the capability of our proposed models and algorithms depends not only on the devices, but also on the way to implement the algorithms. As L-ADMM is the first algorithm designed for our proposed LogSpecT, we will investigate the more efficient algorithms based on it in the future.\n\n**Q7:** Is it possible to link the optimization problem to a maximum likelihood estimate? It seems that the filter h only has an impact on the covariance estimate. Is it true? If the answer is \"no\", would it be a good idea (if possible) to include information about the $h$ filter in the optimization problem?\n    \n**Answer:** The two additional open questions are quite insightful and interesting. For the first one, it is a well-known open problem in GSP. To understand why LogSpecT performs well is essential for designing more efficient models for graph learning from stationary signals and further pushing beyond the stationarity assumption. It is actually one of our future directions. For the second question, the impact of $h$ may lie in two folds. On one hand, it affects the estimated covariance matrix $C_n$, which may consequently affect the efficiency of our proposed algorithm. On the other hand, the prior knowledge of $h$ may help to enlarge the class of graphs that can be perfectly recovered by the spectral template. To understand the influence of $h$ in a rigorous way deserves further investigation."}, {"review_id": "VkhiU18vjJ", "reviewer": "Reviewer_f5Sh", "summary": "In this paper, the authors considered the problem of learning graphs from stationary signals. In order to overcome the infeasibility issue of an existing method called rSpecT, the authors proposed a novel formulation by introducing the log barrier term to learn graphs without isolated nodes. The feasibility can be guaranteed by the new formulation. Furthermore, the recovery guarantees of the new formulation was established, and an efficient algorithm based on linearized ADMM was proposed to solve the formulation.", "questions": "Please answer questions from Weakness.", "limitations": "Yes", "rating": 8, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "This paper is very interesting, with novel ideas to solve an issue of a well-known method in graph signal processing. The authors proposed a novel formulation to overcome the infeasibility issue of rSpecT and established the theoretical recovery guarantees. Furthermore, an efficient algorithm based on linearized ADMM was designed to solve the formulation. Numerical results demonstrated the stability and superiority of the new method.\n", "weaknesses": "1. The descriptions of the experiments conducted on real networks require further improvement. Specifically, please provide more details regarding the number of observations, the nature of these observations, and what the edges represent within the networks.\n\n2. As observed in Figures 2 and 4, the F-measure does not increase to 1, even when given a sufficiently large number of observations. The authors should provide commentary on this phenomenon to address any potential concerns.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you so much for your positive comments and advice. We hope the following clarification can answer your questions.\n\n**Q1:** The descriptions of the experiments conducted on real networks require further improvement. Specifically, please provide more details regarding the number of observations, the nature of these observations, and what the edges represent within the networks.\n    \n**Answer:** The experiments on real networks are conducted on graphs from the Protein database [BOS$^+$05] and the Reddit database [YV15]. The Protein database is a Bioinformatics dataset, where nodes are secondary structure elements (SSEs) and there is an edge between two nodes if they are neighbors in the amino-acid sequence or in 3D space. The Reddit database is a social network dataset. In this dataset, each graph corresponds to an online discussion thread, where nodes correspond to users, and there is an edge between two nodes if at least one of them responded to another’s comments. The observed graph signals are then generated synthetically according to the generative model (1), with the graph filtering chosen as $H(S) = \\exp{\\frac{1}{2}S}$. This is to make sure that the observed graph signals are stationary. For the experiments of real signals, we refer the reviewer to our global response, where we discuss the experiments on the handwritten digits USPS dataset, which is verified to be nearly stationary in [PV17].\n\n**Q2:** As observed in Figures 2 and 4, the F-measure does not increase to 1, even when given a sufficiently large number of observations. The authors should provide commentary on this phenomenon to address any potential concerns.\n    \n**Answer:** Figure 2 presents the results of rLogSpecT on BA graphs. As rLogSpecT is an approximation of LogSpecT, its performance will approach that of LogSpecT when the sample size is growing. However, on BA graphs, even LogSpecT cannot achieve 1 F-measure. Hence, the approximation model rLogSpecT cannot achieve 1 F-measure, no matter how large the sample size is. The phenomenon observed in Figure 4 is similar. The ideal model cannot achieve perfect recovery on any real networks, let alone the approximation model rLogSpecT. It remains an open and intriguing problem in GSP how broad the class of networks is, on which LogSpecT can achieve perfect recovery. This is also one of our future directions.\n\n[BOS$^+$05] Borgwardt, K., et al. Protein function prediction via graph kernels. Bioinformatics, 21:47–56, 2005.\n\n[YV15] Yanardag, P and Vishwanathan, S. Deep graph kernels. In Proceedings of the 21th ACM SIGKDD international conference on knowledge discovery and data mining, pages 1365–1374, 2015.\n    \n[PV17] Perraudin, N and Vandergheynst, P. Stationary signal processing on graphs. IEEE Transactions on Signal Processing, 65(13):3462–3477, 2017."}, {"review_id": "MgsVnbKSNq", "reviewer": "Reviewer_U8JP", "summary": "This paper studies the problem of graph learning from stationary signals. The authors propose an algorithmic framework, LogSpecT, which addresses the shortcomings of an existing body of work that relies on graph learning using spectral templates.\n\n", "questions": "Could you comment on the computational complexity of LogSpecT? In particular, I am curious whether the operations involving $C_n$ in Section 5 could be computationally prohibitive for high dimensional datasets.", "limitations": "Limitations have been adequately discussed. ", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. The problem setting is graph learning for stationary graphs using spectral templates. This is a problem of interest in the graph signal processing domain and is motivated by practical applications. The major technical contribution of the paper is to solve the problem by introducing log barrier in the objective function. Introduction of log barrier facilitates tackling the graph learning scenarios that were not feasible by existing algorithms in this domain.\n2. The paper is well written and the technical content is easy to follow.\n3. The experiments demonstrate the LogSpecT ourtperforms the previous versions of spectral templates and two other baselines from the literature in the task of graph learning. Moreover, it is demonstrated that previous versions of graph learning algorithms may encounter infeasibility scenarios on real datasets with high likelihood.", "weaknesses": "There may be a minor concern that the technical contributions of the paper are narrowly focused on overcoming the shortcomings of graph learning with spectral templates. ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your comments. We now provide responses to the questions and concerns you have raised.\n\n**Q1:** There may be a minor concern that the technical contributions of the paper are narrowly focused on overcoming the shortcomings of graph learning with spectral templates.\n    \n**Answer:** It appears to be a strict requirement that the graph signals are stationary. This assumption starts to gain popularity as the extension of the stationarity defined on the regular domain to the graph data [G15]. Since then, many GSP techniques have been developed under this assumption, e.g. graph learning under stationarity [SM$^+$], low-passness detection under stationarity [ZHW23] and multi-attribute graph signal processing under stationarity [ZHW22]. Along these works, many real datasets have been verified to be stationary or nearly stationary. For example, [G15] showed that a part of the weather data can be explained by stationary graph signals and [PV17] showed that the well-known USPS dataset and the CMUPIE set of cropped faces are nearly stationary. In a nutshell, the stationarity assumption has been a cornerstone in the GSP community. Moreover, it is quite insightful and interesting to push beyond the stationarity assumption and tackle the case when the graph signals are non-stationary. This is one of our future directions.\n    \n**Q2:** Could you comment on the computational complexity of LogSpecT? In particular, I am curious whether the operations involved in Section 5 could be computationally prohibitive for high dimensional datasets.\n \n**Answer:** We consider the complexity in each iteration roughly when the number of nodes $m$ is larger than the number of observations $n$. In the update of $Z$ and $q$, the complexity is decided by the the calculation of $C_nS^{(k)}$. This can be written as $(1/nX)(X^\\top S^{(k)})$, where $X \\in \\mathbb{R}^{m\\times n}$ is the data matrix. Hence, the complexity of this step is in the order of $nm^2 + nm$. In the update of $S$, the complexity is decided by the calculation of $C^{(k)}$. The calculation of $C_n\\Lambda$ can be conducted similarly, yielding the complexity of $nm^2 + nm$. In the calculation of $S^{(k)}C_n^2$, we may decompose it as $\\frac{1}{n^2}(S^{(k)}X)(X^\\top X)(X^\\top)$. The complexity in this step is in the order of $nm^2+n^2m+nm$. In the calculation of $C_nS^{(k)}C_n$, we may calculate it as $\\frac{1}{n^2}X(X^\\top S^{(k)}X)X^\\top$. The calculation involved here is in the order of $nm^2 + mn^2$. The update of the dual variables can be decomposed as the summation of the components calculated in the update of the primal variables. In summary, the total complexity of each iteration is in the order of $nm^2 + n^2m + nm$. The coefficients can be improved by storing some common components, e.g. $S^{(k)}X$ and so on. This would be quite efficient even in a high dimensional dataset.\n\n[G15] Girault, B. Stationary graph signals using an isometric graph translation. In Proceedings of the\n23rd European Signal Processing Conference (EUSIPCO 2015), pages 1516–1520. IEEE, 2015.\n\n[SM$^+$17] Segarra, S., et al. Network topology inference from spectral templates. IEEE Transactions on Signal and Information Processing over Networks, 3(3):467–483, 2017.\n\n[ZHW22] Zhang, C., et al. Product graph learning from multi-attribute graph signals with inter-layer coupling. arXiv preprint arXiv:2211.00909, 2022.\n\n[ZHW23] Zhang, C., et al. Detecting low pass graph signals via spectral pattern: Sampling complexity and applications. arXiv preprint arXiv:2306.01553, 2023.\n\n[PV17] Perraudin, N and Vandergheynst, P. Stationary signal processing on graphs. IEEE Transactions on Signal Processing, 65(13):3462–3477, 2017."}]

Author rebuttal:
We thank all the reviewers for your insightful comments and helpful advice. In the global response, we explain the additional experiments and results shown in the attached pdf. We hope our clarifications can help to clear up reviewers' concerns. 

**1. Stability of LogSpecT and comparison with Kalofolias' method.**

We change the probability that an edge exists in the generation of 50-node ER graphs and use the QUA graph filter to generate stationary graph signals. LogSpecT and Kalofolias' method are conducted on the signals. We present the average F-measure over 10 independent replicates in Figure 1. Specifically, the parameter in Kalofolias' method is chosen from $\{0.01, 0.1, 1, 10, 100, 1000\} $ that yields the highest F-measure. The results show that LogSpecT performs well and stably regardless of the sparsity level. However, Kalofolias' method fails on the stationary graph signals. 

**2. Comparison with more recent models on real networks and synthetic signals.**

We test ALPE [YCP21], NGL-MCP, and NGL-SCAD [YCP20] on the Protein and Reddit datasets with the synthetic stationary graph signals generated from lowpass-EXP: $H(S) = \exp{S}$. The true covariance matrix is assumed to be available. The mean and standard deviation of the F-measure from these models are collected in Table 1. The results indicate that the tested advanced algorithms based on the probabilistic graphical model perform poorly on these two datasets with stationary graph signals. A possible reason is that the stationary assumption does not fit the probabilistic graphical model.

**3. Efficiency of the L-ADMM algorithm.**

Since there does not literally exist customized algorithm for rLogSpecT, we compare our proposed L-ADMM with CVXPY [SS16]. The solver used in CVXPY is MOSEK. We conduct the experiments on 100-node BA graphs with 10000 low-pass EXP generated stationary graph signals. $\delta$ is set as $10\sqrt{\frac{\log n}{n}}$. For L-ADMM algorithm, the accuracy is set as $10^{-6}$ and the initialization is set as zero. The running time for L-ADMM takes around $20$ seconds while the solver takes over $110$ seconds. Note that the complexity of each iteration in L-ADMM can be further improved with the treatment of matrix multiplication mentioned in our response to Reviewer U8JP. We present the primal residual and the dual residual with respect to the iteration index in Figure 2. The results corroborate the efficiency of our proposed algorithm. Furthermore, the local linear convergence can be observed from the residuals, which can lead to fast convergence to high accuracy. It is one of our future directions to prove this property theoretically.

**4. Experiments on real datasets: the handwritten digits USPS dataset.**

As shown in [PV17], USPS dataset is nearly stationary with respect to the nearest neighbour graph. This dataset collects images of handwritten digits. In the experiment, each pixel is viewed as a node and the value on the pixel is view as the graph signal. We follow the approach in [PV17] to construct the 20 nearest neighbour graph, on which the data is verified to be nearly stationary. More specifically, we pick the 1296 digit 1 images. The weights between two nodes are decided by the Gaussian radial basis function. It turns out that the stationarity measure $s := \Vert\operatorname{diag}(U^\top C_n U)\Vert_2 / \Vert U^\top C_n U\Vert_F$ equals 0.78 in this dataset. Here $U$ is the eigen-matrix of the constructed graph and $C_n$ is the covariance matrix of the observed graph signals. We view the 20 nearest neighbour graph as the groundtruth and compare the learned graph from rLogSpecT with it. We present the subgraphs consisting of the top 10 nodes with the largest numbers of neighbours in the two graphs respectively in Figure 3. The F-measure of these two subgraphs is 0.96. This result corroborates the efficiency of our proposed rLogSpecT model.

[PV17] Perraudin, N and Vandergheynst, P. Stationary signal processing on graphs. IEEE Transactions on Signal Processing, 65(13):3462–3477, 2017.

[YCP20] Ying, J., et al. Nonconvex sparse graph learning under Laplacian constrained graphical model. Advances in Neural Information Processing Systems, 33:7101–7113, 2020.

[YCP21] Ying, J., et al. Minimax estimation of Laplacian constrained precision matrices. In International Conference on Artificial Intelligence and Statistics, pages 3736–3744. PMLR, 2021.

[SS16] Steven, D. and Stephen, B. Cvxpy: A python-embedded modeling language for convex optimization. The Journal of Machine Learning Research, 17(1):2909–2913, 2016.

Meta-review: This paper is concerned with learning graph structures from signals with recovery guarantees. It was overall well received by the reviewers, and the rebuttal allowed to fix minor concerns.
Please incorporate the remaining reviewers' feedback for the camera-ready version.
Decision: Accept (poster)