Graphical Models in Heavy-Tailed Markets
Authors: José Vinícius De Miranda Cardoso
Jiaxi Ying
Daniel P. Palomar
Venue: NeurIPS 2021 Poster (2021)
Keywords: graphs
financial markets
precision matrix
Laplacian

Abstract: Heavy-tailed statistical distributions have long been considered a more realistic statistical model for the data generating process in financial markets in comparison to their Gaussian counterpart. Nonetheless, mathematical nuisances, including nonconvexities, involved in estimating graphs in heavy-tailed settings pose a significant challenge to the practical design of algorithms for graph learning. In this work, we present graph learning estimators based on the Markov random field framework that assume a Student-$t$ data generating process. We design scalable numerical algorithms, via the alternating direction method of multipliers, to learn both connected and $k$-component graphs along with their theoretical convergence guarantees. The proposed methods outperform state-of-the-art benchmarks in an extensive series of practical experiments with publicly available data from the S\&P500 index, foreign exchanges, and cryptocurrencies.

Graphical Models in Heavy-Tailed Markets
José Vinícius de M. Cardoso, Jiaxi Ying, Daniel P. Palomar Department Electronic and Computer Engineering Department of Industrial Engineering and Decision Analytics The Hong Kong University of Science and Technology Clear Water Bay, Hong Kong SAR China vjvdmc, jx.ying) ust hk, palomar@ust.
Abstract
Heavy-tailed statistical distributions have long been considered a more realistic statistical model for the data generating process in financial markets in comparison to their Gaussian counterpart. Nonetheless, mathematical nuisances, including nonconvexities, involved estimating graphs heavy-tailed settings pose a sig- nificant challenge to the practical design of algorithms for graph learning. In this work, we present graph learning estimators based on the Markov random field framework that assume Student-t data generating process. We design scalable numerical algorithms, via the alternating direction method of multipliers, to learn both connected and k-component graphs along with their theoretical convergence guarantees. The proposed methods outperform state-of-the-ar benchmarks an extensive series of practical experiments with publicly available data from the S&P500 index, foreign exchanges, and cryptocurrencies.
Introduction
Graph learning frameworks are often designed based on the assumption that the observed graph signals are Gaussian distributed [1-9]. While such assumption for graphical models has found great success in many practical areas, which includes brain network analysis [10], psychological networks [11], and single-cell sequencing [12], it inherently neglects scenarios where there may cxist outliers the underlying data is naturally heavy-tailed distributed. As consequence, those methods often lack robustness and may not succeed in capturing a meaningful representation of the underlying
graph [13].
Data from financial instruments are well-known examples of such scenarios where heavy-tailedness and skewedness are present [1 [14-19]. In addition, there has been growing interest in methods for estimating graphical models financial markets, which hence demands the development of scalable and robust learning algorithms [20].
Perhaps one of the most prominent applications, clustering financial time-series via graph techniques has been an active research topic [20-24]. Nonetheless, current techniques rely on the assumption that the underlying graph has tree structure, which does bring advantages due to hierarchical clustering properties, but also have been shown to be unstable [25-27] and not suitable when the data Gaussian distributed [28].
Motivated by practical challenging applications in finance, such as clustering of financial instruments and network estimation, investigate the problem of learning graph matrices whose structure follow that of a Laplacian matrix of an undirected weighted graph for which the data generating process is assumed be Student-t distributed. particular, the main contributions of this paper are as follows:
We propose a novel formulation for learning undirected weighted graphs under the as- sumption that the data generating process is Student-I distributed. We solve the underlying
35th Conference on Neural Information Processing Systems (NcurIPS 2021).

---

learning problem via carefully designed numerical algorithm based on the alternating direction method of multipliers (ADMM), along with the establishment of its theoretical convergence guarantees. We note that the proposed algorithm can be easily extended to account for additional linear constraints on the graph weights.
We extend the proposed framework to account for heavy-tails and k-component graphs simultancously, which enables novel method for clustering financial time-series.
We present extensive practical results, with real-world data from the US stock market, foreign exchanges, and cryptocurrencies, that showcase clear advantages of including heavy-tail assumptions into graph learning frameworks when compared to state-of-the-art, Gaussian-based methods.
Notation: Matrices (vectors) are denoted by bold, italic, capital (lowercase) roman letters like X, æ. Vectors are assumed to be column vectors. The (i, element of matrix X |R"xp is denoted as Xij- The i-th element of a vector x: is denoted as x;- The i-th row of X is denoted as æ; E RPx11 Given symmetric matrix A, >; (A) and Xmax (A) denote the i-th smallest and maximum eigenvalue of A, respectively. The Moore-Penrose inverse of A is denoted as A The Frobenius norm of matrix A is denoted as | F /tr (A A). The operator Diag : matrix with the elements of an input vector along its diagonal. The operator diag RPxp - Re extracts the diagonal a square matrix. For æ E RR, ||0/|o0 = max,/2;). (a))* denotes the projection on to the nonnegative orthant, i.e., (x)+ max x)..
R RRxx creates diagonal
2 Background & Related Works
An undirected, weighted graph is denoted as triple g = (V,E, E, W), where V 21, 2,.......... is the node set, & S {fu, o} u, E V. u + ta is the edge set, that is, a subset of the set of all possible unordered pairs of nodes such that (u, E & iff there exists link between nodes and W RPxP is the symmetric weighted adjacency matrix that satisfies W. O, Wij 0 iff {i,j) E and W = O, otherwise. The combinatorial, unnormalized graph Laplacian matrix is defined, as L D - W. where D A Diag( W1) is the degree matrix.
A p-dimensional, real-valued, Gaussian random variable æ, with mean vector IE [x] - and rank- deficient precision matrix L, is said to form a Laplacian constrained Gaussian Markov random field (LGMRF) [9, 29-31] of rank p k, k 2 1, with respect to graph G, when its probability density function is given as
exp p(æ) x v det" (L) exp {
- 1 -x - 1)) -T L(x - m)
(I)
where det (L) is the pseudo-determinant of L, i.e., the product of its positive eigenvalues [32].
Assume we are given observations of x, i.e., X = [x1, .2...., The goal of graph learning algorithms is to learn Laplacian matrix, or equivalently an adjacency matrix, given only the data matrix X, i.e., often without any knowledge of &
x,,] X E R"XP, x; E R
To that end, the penalized maximum likelihood estimator (MLE) of the Laplacian constrained precision matrix of æ, on the basis of the observed data X, is:
minimize L=0
tr (LS) - log det (L) + +(L),
subject to L1 = 0, L.j- Ljs 0,
where S is a similarity matrix, e.g., the sample covariance (or correlation) matrix S x X xx, and is a regularization function to promote certain properties L such as sparsity low-rankness.
Even though Problem (2) convex, provided we assume a convex choice for h, it is not adequate to be solved by disciplined convex programming languages, such as cvxpy [33]. particularly due to scalability issues related to the computation of the term log del." (L) [6, 7]. Indeed, recently, considerable efforts have been directed towards the design of scalable, iterative algorithms based on block coordinate descent [34], majorization-miniiize Problem (2) in efficient fashion, e.g., [6] and [7].
(MM) [35, 36), and ADMM [37] to solve

---

Estimators based on Gaussian assumptions have been proposed for connected graphs [2, 4-7]. Some of their properties, such sparsity, are yet being investigated [9. 31]. The authors in [29] and [38] proposed optimization programs for learning the class of k-component graphs, as such class is an appealing model for clustering tasks due to the spectral properties of the Laplacian matrix. However, major shortcoming in their formulations is the lack of constraints on the degrees of the nodes, which
allows for trivial solutions, i.e., graphs with isolated nodes.
Elliptical losses along with linear structural constraints that retain the positive-definitenss of the estimated covariance matrix have been proposed in the literature [39]. this work, however, we address the case of Laplacian constraints, which lead to sentifi-sem precision matrices, and nonconvex k-component structural constraints.
3 Proposed Formulations & Algorithms
In this section, we propose optimization formulations and an iterative algorithm to learn a Laplacian matrix from heavy-tailed assumptions. With that goal, express the Laplacian matrix via its linear operator, ie., L Las [29], where RP(e-1)/2 is the vectorized form of the upper triangular part of the adjacency matrix, also known as the vector of graph weights. In addition, we use the fact that, for connected graphs, it follows that det" (Cw) = det(Lw + J). 111 [6].
In order to address the inherent heavy-tailed nature of financial market data [40], we consider the Student-t distribution under the improper Markov random field assumption [41] with Laplacian structural constraints, that is, we assume the data generating process to be modeled as multivariate zero-mean Student-t distribution, whose probability density function can be written as
V V+
p(æ) & (0)
(1
x
Ox
2
4ン 2,
(3)
where @ is a sosidiie- inverse scatter matrix modeled as combinatorial graph Laplacian matrix and v is the number of degrees of freedom, which measures the rate of decay of the tails.
This results in a robustified version of the MLE for connected graph learning, i.e.,
p+v V
n1 log 1 x. Lwwx;
minimize #20,020
log det (O + J),
v
(4)
subject
Lwv,
dun
d.
where 0 [RP(--1)/2 + IRP is the degree operator defined as dw A diag (Cw). The constraint
dw = d enables the learning of additional graph structures such as regular graphs and it crucial for k-component graphs, as discussed in Section 3.2.
From theoretical perspective, the Student-/ model naturally yields sparse graphs. Comparing the objective function in Problem (4) to that of Problem (2), we note that the Student-f contains log(-) term in place of linear term of the graph weights. The usage of log function to promote sparsity is closely related to the iteratively reweighted -nomr as an approximation for the formorm problem [42]. Problem (4) is, in general, nonconvex due to the summation of log terms and hence is challenging be considered directly. Hence, we design iterative algorithm based on the ADMM framework.
3.1 ADMM Solution
The partial augmented Lagrangian function of Problem (4) is given as
p+ log 1 x._ Lwx;
O, w. Y. v)
log det (0 J)
P d) + 1 low - all + (Y, 0 Lw) L f Lw Lwil
+ (.. dw
(5)
2
where Y and y are the dual variables associated with the equality constraints θ = Lw and dw = d. respectively. Note that we deal with the constraints w > 0 and θ > directly, hence there are no dual variables associated with them.
3

---

Given w' and Y the subproblem for @ can be written as
H+1 arg min log det (e + JJ (0, - l 0
Lw 11:
(6)
@=
whose closed-form solution is given by Lemma 1.
Lemma The global minimizer of problem (6) is 143, 441
1 U (T+V 4pI U
+ Given O Y and y', the subproblem for w can be formulated as
2p
where UTUT is the eigenvalue decomposition of p (Lw' + J) Y
P fw" (" (w. L" (V' pQl+1) - (y) pdd- pd)
minimize w20
2 2u
P+V 71
log 1+ x; Lwx,
where d" and L are the adjoint operators of the degree and Laplacian operators, respectively,
In general, subproblem (8) is nonconvex due to the concave nature of the logarithm function. Hence, we resort to the MM method [36] to find stationary point of subproblem (8). We proceed by
constructing a global upper-bound of the objective function of (8) at point w E IRF [RP(P-1)/2
as
g(n, w') g(wr, w²) + (w w',
+
w w 1
(9)
where f is the objective function in the minimization in (8), its gradient is given as V v=l(ww)
at b¹, where
a' = L g - p -(w))
(10)
s = d" (v' - -
pdd-ow))) dw'))
(11)
(p+ + v a
is weighted sample covariance matrix, and l pdmax (0'D+ L). and the maximum eigenvalue of d-0 L is given by Lemma 2, whose proof is presented in the Supplementary Material.
where 5: DE
(w) (,, 1) + +
Lemma The maximum eigenvalue of the matrix d + L L is given as
Amax (d 0 + L C) 222p 1).
(12)
The vector of graph weights w can then be updated by minimizing the function g constructed in (9), which is tantamount solving the following nonnegative, quadratic-constrained, strictly convex problem:
arg min p(2p 1) ||w - w 112 (w, a² + b²).
(13)
w20
whose unique solution can be readily obtained via its KKT optimality conditions and is given as
w a² + b² 2p(2p- - 1) that is, projected gradient descent step with learning rate (2p(2p 1))- Thus, we iterate (14) in order to obtain a stationary point, 1+1 of Problem (8). In practice we observe that few iterations are sufficient to retrieve (+1
(14)
The dual variables Y and y are updated via gradient ascent steps. Algorithm summarizes the implementation find stationary point of Problem (4). We present Theorem proved in the Supplementary Material, which establishes the convergence of Algorithm 1.
Theorem 3 Algorithm 1 converges subsequently for any sufficiently large that is, the sequence {(O', w², Y Y',y") y generated by Algorithm has least one limit point, and each limit point is stationary point of Problem (4).

---

Algorithm 1: Student-t Graph Learning
Data:
Data matrix E RP*x, initial estimate of the graph weights w", desired degree vector d,
penalty parameter p 0, degrees of freedom 2, convergence tolerance >> 0
Result: Graph Laplacian estimation: Lw"
initialize Y 0, y
2 1 1-0 0
3 while |r'| > or || > e do
> update @ via (7)
> update w'+1 by iterating (14)
> update y Y (0/+1 1w
update y yy p (Dw d)
b compute residual r 0l+1 Lw
> compute residual g'+l = dw
1 < 1+1
11 end
3.2
An extension to k-component graphs
The graph learning formulation proposed in (4) is applicable to learn connected graphs. Learning graphs with components, it assumed to be known, poses a considerably higher challenge, as the dimension of the nullspace of the Laplacian matrix Lw is equal to the number of components of the k,
graph [45]. One way to achieve this requirement is by imposing the constraint rank (Cw) which is nonconvex and nondifferentiable, in the maximum likelihood problem generated by (3). We instead relax this rank constraint by noting that via Fan's 's theorem [46]. we have
A (Lw)= minimize tr (VT LWV).
(15)
VERUXM V
Thus, by using the right hand side of (15) as a regularization term, we are able formulate the following optimization problem learn Student-t k-component graph:
p+ V
log 1 x: Lwxx;
minimize w20.0-0.V
log det (e) + (tr((LVVVV)
(16)
subject to
Lw, rank( = p- k. dw d, V I, V E RPXA
The partial augmented Lagrangian function of Problem (16) can be expressed as
L,(O, w. V.Y
+
n Σ log 1 z. LWx,
log det" (0) rjt (CwVVT)
y dw
d) dw
d + YY, Q Lw)
+ p |
Lw (17)
2
Given w and Yd, the subproblem for can be written as
arg min - log det (0) + (0,Y))
0 Lw
(18)
2
rank =-p-k
which is nearly the same as (6). Its solution is obtained as
Lu (T T + 4pI UT.
r+
+
(19)
2p
except that now UTUT is the eigenvalue decomposition of pLw' Y with having the largest p - k cigenvalues along its diagonal and U (p contains the corresponding eigenvectors.
The subproblem to obtain =11 is virtually the same as in (8) except for the additional linear term Vin(LuV.'') Hence, its update is also projected gradient descent step, alike (14) where
á 2 L 5)) + V'V'T Y - ( Lwi))
(20)
5

---

Given w we have the following subproblem for V:
minimize VER .VT tr (VT Lw
(21)
whose closed-form solution is given by the à cigenvectors associated with the hs smallest cigenvalues of Lwt+1 [47, 48]. Algorithm summarizes the implementation to find stationary point of Problem (16), and its convergence is established through Theorem 4, whose proof is presented in the Supplementary Material.
Theorem 4 Algorithm 2 converges subsequently for any sufficiently large that is, the sequence {(O', w', Vi, Y! } generated by Algorithm has at least one limit point, and each limit point is stationary point of Problem (16).
Algorithm 2: k-component Student-t graph learning
Data: Data matrix X RaXp, initial estimate of the graph weights w", number of graph components desired degree vector d, degrees of freedom v, rank hyperparameter 7) > 0, penalty parameter p >> 0, tolerance >> 0
Result: Laplacian estimation: Lw"
initialize Y = 0, y
2 1 0-0
3 while |r'|| > or || 11 > e do
> update e via (19)
b update w/ as in (14) with a² given in (20)
> update V as in (21)
b update y = Y + (0/+1 1
update y/+! y p (dw't1 d)
compute residual nl 1 Lw
> compute residual s dw
10
11
1 + 1+1
12 end
4 Experiments
To evaluate the performance of the proposed graph learning algorithms, we perform experiments using historical daily price time series data, available in Yahoo! Finance from financial instruments in three scenarios: (i) stocks belonging to the S&P500 index, (ii) foreign exchange markets, and (iii) eryptocurrencies. We start by constructing the log-returns data matrix, i.e., a matrix X E IRexp, where n is the number of log-return observations and is the number of instruments, as
Xij log P.,j log P
(22)
where Pr,j is the closing price of the j-th instrument at the i-th day.
Benchmarks: We compare our proposed algorithms with state-of-the-art, Gaussian distribution-bast methods for connected graphs, namely GLE [7] and NGL [9], which use f--norm and minimax concave penalty regularizations, respectively: and CLR [38] and SGL [29] that consider k-component graphs. For a fair comparison among algorithms, we set the degree vector d equal to for the proposed algorithms, i.e., we do not consider any prior information on the degree of nodes. In our ADMM algorithms, we set the penalty parameter to p 1 and the hyperparameter 1j (16) adaptively increased until the rank constraint is satisfied, For GLE and NGL, we use grid search on the sparsity hyperparameter such that the resulting graph yields the highest modularity value. The graph weights in Algorithm I and are initialized using the same procedure in [8].
Our goal with the experiments that follow to verify whether the heavy-tail assumption provides an improved version of the learned graph, which is evaluated based on the modularity of the estimated
Modularity measures the strength of separability of graph into groups [49].

---

graph and the graph visualization. In addition, for the task of clustering stocks, we analyze whether the learned graphs agree with industry standards of sector classification set by the Global Industry Classification Standard (GICS) [50, 51].
4.1 Communities in S&0 Stocks
In this experiment, we consider S&P500 stocks belonging to three sectors, namely, Communication Services (red), Utilities (blue), and Real Estate (green), totalling = 82 stocks, during the time horizon from Jan. 3rd 2014 to Dec. 29th 2017, resulting in n 1006 observations. In order to obtain descriptive insights on this dataset, we measure its degree of heavy-tailedness and annualized volatility'. The former is obtained by fitting the degrees of freedom of Student-t distribution the matrix of log-returns, whereby we obtain // z 5.5 and a R 21%. This scenario can be considered as having a moderate amount of heavy-tailedness.
Figure/ depicts the learned connected graphs on the aforementioned time periods. It can be readily noticed that the graph learned with the Student-t. distribution (Figure Ic) is sparser than those learned with the Gaussian assumption (Figure la and 1b), which results from the fact that the Gaussian distribution is more sensitive to outliers. Moreover, the Student-l graph presents higher degree of interpretability as measured by its modularity value. In addition, larger number of inter-sector connections, indicated by gray-colored edges, which are often spurious from a practical perspective, are present in the graphs learned by NGL and GLE. Sparsity regularization provides a means to remove edges between nodes in the presence of data with outliers and possibly increasing the modularity of the resulting graph. However, they bring the additional task of tunning hyperparameters, which applications. A cleaner graph, without the need for
is often repetitive and impractical for postprocessing or additional regularization, is obtained directly by using the Student-t assumption.
(a) GLE, modularity 0.31.
(b) NGL, modularity 0.49. (c) Algorithm modularity 0.54.
Figure 1: Learned graphs of S&P500 stocks.
Figure 2 illustrates the learned 3-component graphs during the time from Jan. 3rd 2014 until Dec. 29th 2017. We can notice that SGL (Figure 2a) and CLR (Figure 2b) are unable to separate the stocks in way that agrees with their sector information as given by GICS. In addition, the high number of spurious connections (gray-colored edges) are uncharacteristic of the actual expected behavior in stock markets. Figure 2c displays the graph learned by the proposed Algorithm 2, where it presents not only a higher modularity value, but also sparser, more plausible representation of an actual network of stocks with three sectors.
4.1.1 Impact of COVID-19 in the US Stock Market
We collect price data of 85 stocks belonging to three sectors, namely, Communication Services (red), Utilities (blue), and Real Estate (green) from Apr. 22nd 2019 to Dec. 30th 2020, resulting in 429 observations. The degrees of freedom and annualized volatility during this period, which includes the financial crisis caused by the COVID-19 pandemic, were v 2.89 and a ze 41%, indicating a strongly heavy-tailed scenario.
"More often than not, stocks have impacts on multiple sectors, egg. the evident case of technology companics, whose influence affect the prices stocks not only in their own sector, but spans across multiple sectors. Therefore, GICS or other sector classification systems cannot considered as absolute ground-truth labels.
'The annualized volatility is computed as a 1/52 Σ Ti. where a, is the daily sample standard deviation of the i-th stock.
7

---

2015
(a) SGL. modularity = 0.29.
(b) CLR, modularity = 0.33. (c) Algorithm 2, modularity = 0.56.
Figure 2: Learned 3-component graphs of S&P500 stocks.
In Figure 3, we observe that the modularity value of the Caussian-based graphs are severely reduced as compared to the results presented in Figure I, while the opposite occurs with the Student-t graph. That may be explained by the increase in variability and outliers in the data during the financial crisis caused by the COVID-19 pandemic, which is better modelled by a Student-t distribution. the learned graph in Figure 3c clearly displays expected behaviors such the strong correlation between pairs of companies including {Zoom Video Communications Inc., Netflix Inc.), ), and (Twitter Inc., Facebook Inc.), which are not as evident Figures 3a and 3b.
addition,
rebalry
(a) GLE, modularity 0.23.
(b) NGL, modularity 0.46.
(c) Algorithm modularity 0.56.
Figure 3: Learned graphs of S&P500 stocks during the COVID-19 pandemic.
4.2 Communities in Foreign Exchange Markets
We query foreign exchange data from the most traded currencies between the period from Jan. 2nd 2019 to Dec. 31st 2020, totalling log-returns of the currencies prices with respect to the United States Dollar. Similar to the previous experiment, compute the degrees of freedom of Student-t distribution fitted to the log- returessss data matrix, whereby we obtain =s 4.6, which represents a scenario with considerable amount of heavy-tailedness Unlike in the experiment involving S&P500 stocks, there
522 observations. The data matrix is composed by the
classification standard for currencies, hence we rely on a community detection algorithm [52] in order to create classes within the learned graph. In particular, the algorithm in [52] takes input the learned Laplacian matrix of the graph and outputs a membership assignment that maximizes the modularity of the graph.
Figure 4 displays the Iearned graphs. As can be observed, the Student-t graph (Figure 4c) is sparser, more interpretable, and has higher modularity value than that of the Gaussian- based graphs (Figure 4a and 4b). In addition, the expected correlation between currencies of locations geographically close to each other, e.g., (Hong Kong SAR. China), (Taiwan, South Korea}, and (Poland, Czech Republic] are significantly more evident for the Student-t graph
Figure 5 depicts the learned 9-component graphs of currencies during the time window from Jan. 2nd 2019 to Dec. s1 2020. It can observed that the graph learned by the proposed Algorithm (Figure 5c) presents finer structure and higher modularity value than those learned by SGL (Figure 5a) and CLR (Figure 5b). addition, the learned graph in Figure presents more reasonable clusters such as [New Zealand, Australia] and (Poland, Czech Republic, Hungary), which not

---

(a) GLE, modularity 0.34.
(b) NGL, modularity 0.46. (c) Algorithm modularity 0.58.
Figure 4: Learned connected graphs of currencies.
separated in the Gaussian-based graphs. More critically, we can observe that SGL and CLR allow the existence of isolated nodes in the learned graphs. In our proposed algorithm, we avoid such solutions by imposing linear constraints on the degree of the graph.
(a) SGL, modularity 0.62.
CLR, modularity 0.79.
Algorithm modularity 0.84.
Figure 5: Learned 9-component graphs of currencies.
4.3 Communities in Cryptocurrencies
We query daily prices of the p 41 most traded cryptocurrencies during the period starting from Aug. 1st 2017 to Dec. 1st 2020, which amounts n 1218 observations. The degrees of freedom during this time frame was measured as y Rs 3, which is tantamount to strong heavy-tail scenario.
Figure 6 shows the learned graphs by GLE, NGL, and our proposed Algorithm 1. While the graphs in Figures 6a and 6b present a small modularity value and contain a large number of edges, which impairs interpretability. the resulting proposed graph in Figure 6c reveals refined representation of the interactions between pairs of cryptocurrencies, which is possibly more aligned with the actual market scenario. As example, the link between Bitcoin (BTC) and Litecoin (LTC), Bitcoin spinoff established in 2011 1, is substantially more evident in Figure 6c.
poa
STEEM
BDGE
GNT
FMD
VEM
MD0
215
STREM
EIC 250
350
(a) GLE, modularity 0.19.
(b) NGL, modularity 0.40.
(c) Algorithm modularity 0.52.
Figure 6: Learned connected graphs of cryptocurrencies.

---

Figure 7 shows the learned 7-component graphs of cryptocureencies during the aforementioned time window. As in the previous experiments with foreign exchange data, SGL shows isolated nodes in the learned graph (Figure 7a) and CLR contains large number of spurious connections in the main cluster (Figure 7b), whereas the graph learned via Algorithm (Figure 7c) has the largest modularity value. Interestingly, while all three methods agree to cluster (Dogecoin (DOGE), Verge (XVG), Siacoin (SC), DigiByteCoin (DGB)). which may be related to their similar initial release dates, only our proposed algorithm clusters together the coins that mainly focus privacy and anonymity features, i.e., {Monero (XMR), Zcash (ZEC), DASH].
MANG BAT
ENTT
VARVES
MED
AND
MOKA
2C1
ETC
PTC
VAOTA
VUME KIP
STEEN
SGL, modularity 0.36
(b) CLR, modularity 0.66.
(c) Algorithm modularity 0.79.
Figure 7: Learned 7-component graphs of cryptocurrencies.
4.4 Effect of Algorithm Initialization
Initialization is critical phase any iterative algorithm, especially when dealing with nonconvex optimization problems involving nonconvex constraints as it is the case of Algorithm 2. To verify the stability of Algorithm 2, we initialize it as fully connected graph with equal weights, i.e., w x 1. Figure 8 illustrates that there is no significant difference from the estimated graphs with uniform initial point and the ones reported the previous section. This result provides evidence that Algorithm 2 is robust against inaccurate initializations.
LOCE DGA
(a) modularity = 0.56.
(b) modularity = 0.80.
(c) modularity = 0.78.
Figure 8: Learned graphs by Algorithm 2 with an uniform initial point.
5 Conclusions
Heavy-tails are prevalent in time-series of financial markets. Yet, they have been little explored in the context of Laplacian graphical models. In this paper, we have proposed optimization programs to learn graphical models with Laplacian constraints assuming that the data generating process is Student-t distributed. The formulations follow a maximum likelihood approach of Markov random field, for which we designed ADMM algorithms that converge to a stationary point of the resulting nonconvex problems. The proposed algorithms showed significant gains, measured via the modularity values the estimated graphs, when compared to state-of-the-arr counterparts in real-world scenarios that involved data from the US stock market, foreign exchange markets, and cryptocurrencies.
10

---

Acknowledgments
We would like to thank the anonymous reviewers for their helpful comments. This work was supported by the Hong Kong GRF 16207820 research grant.
References
[1] J. Friedman, T. Hastie, and R. Tibshirani. Sparse inverse covariance estimation with the graphical lasso. Biostatistics, 9:432 41, 2008.
[2] B. M. Lake Jnd B. Tenenbaum. Discovering structure by learning sparse graph. In Proceedings of the 33rd Annual Cognitive Science Conference, 2010.
[3] M. Slawski and M. Hein. Estimation of positive definite m-matrices and structure learning for attractive Gaussian Markov random fields. Linear Algebra and its Applications, 473:145-179, 2015.
[4] X. Dong. Thanou, Frossard, and Vandergheynst. Learning Laplacian matrix in smooth graph signal representations. IEEE Transactions on Signal Processing, 64(23):6160-6175 2016.
[5] V. Kalofolias. How learn graph from smooth signals. Proceedings of the 19th In- ternational Conference on Artificial Intelligence and Statistics, volume 51, pages 920-929, 2016.
[6] H. Egilmez, E. Pavez, and A. Ortega. Graph learning from data under Laplacian and structural constraints. IEEE Journal of Selected Topics in Signal Processing, 11(6):22--841, 2017.
[7] L. Zhao, Wang, S. Kumar, and D. P. Palomar. Optimization algorithms for graph Laplacian estimation via ADMM and MM. IEEE Transactions on Signal Processing, 67(16):4231-22) 2019.
[8] Kumar, J. Ying. V. M. Cardoso, and D. P. Palomar. unified framework for structured graph learning via spectral constraints. Journal of Machine Learning Research, 21:1-60, 2020.
[9] J. Ying, J. V. de M. Cardoso, and D. P. Palomar. Nonconvex sparse graph learning under laplacian-structureed graphical model. In Advances Neural Information Processing Systems (NeurIPS), 2020.
[10] Stephen M. Smith, Karla L. Miller, Gholamreza Salimi-Khorshie Khorshidi, Matthew Webster, Christian F. Beckmann, Thomas E. Nichols, Joseph D. Ramsey, and Mark W. Woolrich. Network modelling methods for fmri. Neurolmage, 54(2):875 -
891, 2011.
[11] S. Epskamp, D. Borsboom, and E. Fried. Estimating psychological networks and their accuracy: A tutorial paper. Behavior Research Methods, 50: 195-12 2018.
[12] O. Stegle, S. A. Teichmann, and J. C. Marioni. Computational and analytical challenges in single-cell transcriptomics. Nature Reviews Genetics, 16:133-145, 2015.
[13] Michael Finegold and Mathias Drton. Robust graphical modeling with I-distributions. In Conference on Uncertainty in Artificial Intelligence, 2009.
[14] C. Gourieroux and A. Monfort. Time Series and Dynamic Models. Themes in Modern Econometrics. Cambridge University Press, 1997.
[15] R. Cont. Empirical properties of asset returns: stylized facts and statistical issues. Quantitative Finance, 1:223-236, 2001.
[16] R. S. Tsay. Analysis of Financial Time Series. Wiley, 3rd edition, 2010.
[17] A. C. Harvey. Dynamic models for volatility and heavy tails: with applications to financial and economic time series. Cambridge University Press, 2013.
[18) Y. Feng and D. Palomar. A signal processing perspective on financial engineering. Foundations and Trends in Signal Processing. 9: 9:1-231, 2015.
[19] Nassim Dehouche. Scale matters: The daily, weekly and monthly volatility and predictability of bitcoin, gold, and the s&p 500, 2021.
[20] G. Marti, Nielsen, M. Binkowski, and P. Donnat. A review of two decades of correlations, hierarchies, networks and clustering in financial markets. In arXiv: 1703.00485 2017.

---

[21] R. N. Mantegna. Hierarchical structure in financial markets. The European Physical Journal B, 11(1):193-19) 197, 1999.
[22] C. Dose and S. Cincotti. Clustering of financial time series with application to index and enhanced index tracking portfolio. Physica A: Statistical Mechanics and its Applications, 355(1):145 151, 2005.
[23] G. Marti, S. Andler, F. Nielsen, and P. Donnat. Clustering financial time series: How long is enough? In Proceedings of the Twenty-Fifth International Joint Conference on Artificial Intelligence, 2016.
[24] G. Marti, F. Niclsen, P. Donnat, and S. Andler. On clustering financial time series: need for distances between dependent random variables. Computational Information Geometry, pages 149-174, 2017.
[25] G. Carlsson and F. Mémoli. Characterization, stability and convergence of hierarchical clustering methods. Journal of Machine Learning Research, 11:1425- 1470, 2010.
[26] V. Lemieux, and their effect on portfolio formation and risk analysis. In Proceedings of the International Workshop on Data Science for Macro-Modeling, page 1-6, 2014.
Rahmdel, R. Walker, B. L. W. Wong, and M. Flood. Clustering techniques
[27] G. Marti, P. Very, P. Donnat, and F. Nielsen. A proposal of methodological framework with experimental guidelines to investigate clustering stability on financial time series. In 2015 IEEE 14th International Conference on Machine Learning and Applications (ICMLA), pages 32-37, 2015.
[28] P. Donnat, G. Marti, and P. Very. Toward a generic representation of random variables for machine learning. Pattern Recognition Letters, 70:24-31, 2016.
[29] S. Kumar, J. Ying, J. V. de M. Cardoso, and D. P. Palomar. Structured graph learning via Lapla- cian spectral constraints. Advances in Neural Information Processing Systems (NeurIPS), 2019.
[30] J. Ying, J. V. de M. Cardoso, and D. P. Palomar. Does the f--norm learn sparse graph under Laplacian constrained graphical models? arXiv e-prints: 2006. June 2020.
[31] Ying, V. M. Cardoso, and D. P. Palomar. Minimax estimation of Laplacian constrained precision matrices. In 24th International Conference on Artificial Intelligence and Statistics (AISTATS' 2021.
[32] O. Knill. Cauchy- Binet for pseudo deeuddddetrrmsnntss. Linear Algebra and its Applications, 459:522- 547, 2014.
[33] S. Diamond and S. Boyd. CVXPY: A Python-embedded modeling language for convex optimization. Journal of Machine Learning Research, 17(83):1-5, 2016.
[34] S. J. Wright. Coordinate descent algorithms. Mathematical Programming, 151:3-34, 2015.
[35] David Hunter and Kenneth Lange. A tutorial on mm algorithms. The American Statistician, 58(1):30-37, 2004.
[36] Y. Sun, P. Babu, and D. P. Palomar. Majorization- algorithms in signal processing, communications, and machine learning. IEEE Transactions on Signal Processing, 65(3):794- 816. 2017.
[37] S. Boyd, N. Parikh, E. Chu, B. Peleato, and J. Eckstein. Distributed optimization and statistical learning via the alternating direction method multipliers. Foundations and Trends Machine Learning, 3(1): 2011 1.
[38] F. Nie, X. Wang, M. I. Jordan, and H. Huang. The constrained Laplacian rank algorithm for graph-based clustering. In Proceedings of the Thirtieth AAAI Conference on Artificial Intelligence, AAAI'16, pages 1969-19766 2016.
[39] Y. Wald, N. Noy, G. Elidan, and A. Wiesel. Globally optimal learning for structured elliptical losses. In Advances in Neural Information Processing Systems (NeurIPS' 19), 2019.
[40] S. I. Resnick. Heavy-Tail -Til Phenomena: Probabilistic and Statistical Modeling. Springer- Verlag New York, 2007.
[41] H. Rue and L.. Held. Gaussian Markov Random Fields: Theory And Applications. Chapman Hall/CRC, 2005.
12

---

[42] Emmanuel J. Candès, Michael B. Wakin, and Stephen P. Boyd. Enhancing sparsity by reweighted ε minimization. Journal of Fourier Analysis and Applications, 14(5):877-905. 2008.
[43] D. M. Witten and R. Tibshirani. Covariance-regularizze regression and classification for high dimensional problems. Journal of the Royal Statistical Society. Series B (Statistical Methodology), 71(3):615-636 2009.
[44] Danaher, P. Wang. and D. M. Witten. The joint graphical lasso for inverse covariance estimation across multiple classes. Journal of the Royal Statistical Society Series B, 76(2):373- 397, 2014.
[45] R. K. Chung. Spectral Graph Theory, volume 92. CBMS Regional Conference Series in Mathematics, 1997.
[46] K. Fan. On theorem of Weyl concerning cigenvalues of linear transformations Proceedings of the National Academy of Sciences, 35(1 1):652665,, 1949.
[47] R. A. Horn and C. R. Johnson. Matrix Analysis. Cambridge University Press, 1985.
[48] P.-.. Absil, R. Mahony, and R. Sepulchre. Optimization Algorithms on Matrix Manifolds. Princeton University Press, Princeton, NJ, 2007.
[49] M. E. J. Newman. Modularity and community structure in networks. Proceedings of the National Academy of Sciences of the United States of America, 103, 2006.
[50] Standard & Poor's. Global Industry Classification Standard (GICS). Tech Report, 2006.
[51] Morgan Stanley Capital International and S&P Dow Jones. Revisions to the global industry classification standard (GICS) structure, 2018.
[52] A. Clauset, M. E. J. Newman, and C. Moore. Finding community structure in very large networks. Physical Review E. 70, Dec 2004.
13


Reviews:
[{"review_id": "gcEJYtXO8w6", "reviewer": "Reviewer_C54j", "summary": "The paper proposes a new estimator for learning graphical models under the t-distribution assumption for data generation process. For the new estimator, an ADMM-based algorithm is designed to solve the associated optimization problem in order to find the graph structures. The new estimator is evaluated extensively on real-world datasets to demonstrate its effectiveness. The paper also establishes the convergence of the proposed algorithm.\n", "questions": "", "limitations": "", "rating": 7, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Originality:\n\nThe proposed estimator is not entirely new in the literature, given that t-distribution has been used in graphical models before (see [1]). Also, it is a direct result of well-known maximum-likelihood principle, so deriving such estimator should not be difficult. For the ADMM-type learning algorithm, I acknowledge that there are numbers of technical details to be figured out in order to have the exact algorithm, but at the high level this seems to be a standard application of ADMM recipe. So it seems to me that the novelty of this paper is somewhat limited.\n\n\nQuality:\n\nIn terms of technical quality, though the ADMM algorithm seems straightforward, it is nice to see that the paper proves its convergence to stationary points, which requires some efforts and skills. Regarding the algorithm, one thing that is not clear to me is the reason why the trace regularization term is needed in Eq. (16). Conceptually if we drop the regularization term, the new formulation is still equivalent to problem (16) (in the sense that they would end up with the same set of minimizers if globally optimality could be achieved), as the constraints $rank(\\Theta) = p - k$ and $\\Theta = Lw$ are imposed. I understand that due to the nonconvexity, ADMM algorithm could end up with different solutions for the two formulations (with or without the regularization term), but it is not obvious to me that the regularized one could work better or there is a need to introduce it. This regularization seems to add some novelty to this paper, so I hope that either theoretical or empirical comparison could be provided to support it (e.g., faster convergence or better quality of solutions), otherwise the additional computational cost seems not well-justified.  \n\nFor the practical perspective, I like the extensive evaluation of the proposed algorithm on different datasets, showing its usefulness and advantage, but it would be better to extend the comparison beyond the algorithms that just work for Gaussian distributions. Below are a few questions I have for the experiments.\n\n1. For S&P500 stocks, is there a particular reason for choosing the three GICS sectors (Communication Services, Utilities and Real Estate), as opposed to other sectors?\n2. Have you considered using some type of residual returns from econometrics literature instead of log returns? If you use GICS sectors as the ground truth for clustering, residual returns may be preferred, as sometimes the co-movement of stock returns is mainly driven by some econometrics factors, which may connect stocks across different sectors.  \n\n\nClarity:\n\nOverall the presentation of the paper is clear. The description of the problem setup and the derivation of the algorithm are easy to follow. The figures shown in the experiments give an nice qualitatively demonstration of the superiority of the algorithm.\n\nSignificance:\n\nThe proposed algorithm seems to be readily applicable to a few real-world problems, as demonstrated in the experiment section. So it should be a good addition to this area for handling heavy-tailed distributions.\n\n[1] Finegold, Michael A., and Mathias Drton. \"Robust graphical modeling with t-distributions.\" In UAI 09, pp. 169-176. 2009. \n\n\n\n********** post-rebuttal changes **********\n\nAfter reading author response (especially on the concern about algorithmic novelty), I would like to increase my rating from 6 to 7.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Yes.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "5"}, "extra_text": {}}, {"review_id": "fpupL80il_v", "reviewer": "Reviewer_M9op", "summary": "The paper proposes a framework for learning a graphical model when the data generating distribution is heavy-tailed. More specifically, the author(s) propose a graph learning estimator under Markov Random Field framework when the data is generated from a Student's t-distribution. A scalable algorithm is proposed to handle connected and K-component graphs. The methodology is particularly applicable to financial time series data. The algorithm is shown to be superior in performance compared to existing methods in both simulated and real data applications.", "questions": "", "limitations": "", "rating": 7, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The paper proposes a novel algorithm based on the alternating direction method of multipliers (ADMM) to a pertinent problem in the graphical model viz. estimating graph structure in the heavy-tailed distribution setting. The problem has a number of applications in financial time series especially clustering the financial time series. The algorithm looks technically correct and the numerical performances are good. However, I have following queries for the author(s) mainly centred on the framework and the methodology\n\n-I understand that the lack of constraints on the degree of the nodes may lead to graphs with isolated nodes however, I am not sure whether in real world applications one will always have the number d available. The information on the graph whether it is sparse or dense may be known but the exact number d may not be known. Any other soft constraints to consider here when d is not known?\n\n-Both theorem 3 and theorem 4 discusses the convergence of the respective algorithms for large $\\rho$. But, can author(s) shed some light about any connection with the number of iterations required to achieve the convergence?\n\n-I am not sure whether in the K-connected components case the number k would be always known. One may has to do an estimation of K at the first stage. Am I missing something?\n\n-Do you need any further regularization on Lw in equation (4) to achieve sparsity?\n\n-How does one choose the rank hyperparameter $\\eta$?\n\n\n\n******** Change after author rebuttal ************\n\nI am happy with the technical justification author(s) provided about the degree of the network and the convergence of the algorithm. Further minor technical details about the regularization and the connected components have been appropriately answered. I have now changed my score to \"acceptance\".\n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "I didn't see any negative societal impact of their work.", "review_text": "", "extra_scores": {"ethical_concerns": "None.", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4 hours"}, "extra_text": {}}, {"review_id": "K5-yqx6NC1n", "reviewer": "Reviewer_7Jw4", "summary": "- This paper proposes to use the student-t distribution for data generation in financial markets, which can capture the effect of outliers better than the commonly used Gaussian distribution. \n\n- In particular, the data is modeled using a Laplacian Markov random field with student-t distribution, and a practical ADMM based algorithm to learn the underlying graph structure is proposed.\n\n- The algorithms are designed under two settings namely, connected graphs and k-component graphs, and convergence is proven. \n\n- The proposed approach is shown to outperform the state-of-the-art Gaussian assumption-based approaches with real-world data from stock markets, foreign exchanges, and cryptocurrencies. \n\n", "questions": "", "limitations": "", "rating": 8, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "- The paper addresses an important problem of designing graph learning algorithms under heavy-tailed data distributions which is relevant in the financial domain but challenging due to its non-convex nature. \n\n- The convergence of proposed algorithms is analytically proven and simulations suggest better performance compared to existing algorithms. \n\n- The paper is written well and easy to understand. The mathematical details are clear. ", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": " A section introducing the preliminary concepts such as MRF and the motivation to use them in financial data would be helpful.", "review_text": "", "extra_scores": {"ethical_concerns": "Nil", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "12"}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper formulates a graphical model learning problem by maximizing the likelihood of a multivariate student-t distribution subject to constraints on the Laplacian structured parameter matrix. The formulation is non convex. An explicit ADMM iteration is derived to arrive at a stationary point. 

The primary strength of the paper is in its multiple experiments demonstrating the utility in the context of financial time-series of the modeling assumptions over the more commonly studied Gaussian graphical model setup. 

The main weakness of the paper is that it does not contain a great deal of conceptual depth or technical novelty: the ideas, and their implementation, appear to be relatively straightforward. 

Overall, the paper seems likely to have impact by spurring further work on estimation of graphical models of various types. 
Decision: Accept (Poster)