Dual Averaging Method for Online Graph-structured Sparsity

TL;DR

This paper introduces GraphDA, a novel online learning algorithm that efficiently handles complex graph-structured sparsity constraints, improving classification and interpretability in large-scale, real-time graph data applications.

Contribution

The paper presents the first online algorithm for graph-structured sparsity constraints, projecting gradients and variables onto lower-dimensional subspaces to capture graph structure.

Findings

Outperforms baseline methods in classification accuracy

Effectively captures graph-structured features

Enhances interpretability of models

Abstract

Online learning algorithms update models via one sample per iteration, thus efficient to process large-scale datasets and useful to detect malicious events for social benefits, such as disease outbreak and traffic congestion on the fly. However, existing algorithms for graph-structured models focused on the offline setting and the least square loss, incapable for online setting, while methods designed for online setting cannot be directly applied to the problem of complex (usually non-convex) graph-structured sparsity model. To address these limitations, in this paper we propose a new algorithm for graph-structured sparsity constraint problems under online setting, which we call \textsc{GraphDA}. The key part in \textsc{GraphDA} is to project both averaging gradient (in dual space) and primal variables (in primal space) onto lower dimensional subspaces, thus capturing the graph-structured sparsity effectively. Furthermore, the objective functions assumed here are generally convex so as to handle different losses for online learning settings. To the best of our knowledge, \textsc{GraphDA} is the first online learning algorithm for graph-structure constrained optimization problems. To validate our method, we conduct extensive experiments on both benchmark graph and real-world graph datasets. Our experiment results show that, compared to other baseline methods, \textsc{GraphDA} not only improves classification performance, but also successfully captures graph-structured features more effectively, hence stronger interpretability.

Tables

Table 1. Table 1. Classification performance on Graph01 of Benchmark dataset

Method	${\text{Pre}}_{{𝒘}_t}\pm$std	${\text{Rec}}_{{𝒘}_t}\pm$std	${\text{F1}}_{{𝒘}_t}\pm$std	${\text{AUC}}_{{𝒘}_t,{\overline{𝒘}}_t}$	${\text{Acc}}_{{𝒘}_t,{\overline{𝒘}}_t}$	${\text{Miss}}_{{𝒘}_t,{\overline{𝒘}}_t}$	${\text{NR}}_{{𝒘}_t,{\overline{𝒘}}_t}$
Adam	0.024$\pm$0.00	1.000$\pm$0.00	0.047$\pm$0.00	(0.618, 0.603)	(0.619, 0.603)	(166.35, 173.10)	(100.0%, 100.0%)
${\ell }_1$-RDA	0.267$\pm$0.11	0.863$\pm$0.09	0.389$\pm$0.13	(0.693, 0.672)	(0.694, 0.673)	(155.30, 166.05)	(11.58%, 83.60%)
AdaGrad	0.256$\pm$0.11	0.877$\pm$0.09	0.379$\pm$0.13	(0.696, 0.636)	(0.696, 0.637)	(156.00, 166.00)	(11.33%, 100.0%)
DA-GL	0.176$\pm$0.11	0.967$\pm$0.04	0.283$\pm$0.12	(0.735, 0.666)	(0.735, 0.667)	(142.90, 162.20)	(15.99%, 100.0%)
DA-SGL	0.523$\pm$0.40	0.854$\pm$0.14	0.506$\pm$0.35	(0.699, 0.647)	(0.699, 0.647)	(151.00, 165.50)	(25.54%, 100.0%)
StoIHT	0.057$\pm$0.04	0.150$\pm$0.08	0.072$\pm$0.03	(0.552, 0.523)	(0.553, 0.523)	(194.55, 195.25)	(7.79%, 40.62%)
GraphStoIHT	0.151$\pm$0.12	0.356$\pm$0.16	0.194$\pm$0.12	(0.603, 0.570)	(0.602, 0.570)	(174.65, 181.40)	(7.84%, 22.06%)
DA-IHT	0.507$\pm$0.20	0.744$\pm$0.12	0.566$\pm$0.11	(0.697, 0.666)	(0.697, 0.666)	(155.65, 162.85)	(4.35%, 39.50%)
GraphDA	0.869$\pm$0.13	0.906$\pm$0.04	0.880$\pm$0.08	(0.749, 0.739)	(0.749, 0.739)	(133.45, 136.20)	(2.56%, 32.12%)

Table 2. Table 2. Performance of Graph02

Method	${\text{Pre}}_{{𝒘}_t}$	${\text{Rec}}_{{𝒘}_t}$	${\text{F1}}_{{𝒘}_t}$	${\text{AUC}}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{Acc}}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{Miss}}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{NR}}_{{𝐰}_t,{\overline{𝐰}}_t}$
Adam	0.042	1.000	0.081	(0.697, 0.663)	(0.696, 0.663)	(144, 151)	(100.0%, 100.0%)
${\ell }_1$-RDA	0.371	0.876	0.494	(0.772, 0.732)	(0.772, 0.731)	(127, 140)	(13.31%, 96.47%)
AdaGrad	0.342	0.888	0.470	(0.771, 0.711)	(0.771, 0.711)	(125, 141)	(14.43%, 100.0%)
DA-GL	0.270	0.976	0.415	(0.809, 0.755)	(0.809, 0.755)	(114, 138)	(17.07%, 100.0%)
DA-SGL	0.283	0.948	0.314	(0.777, 0.738)	(0.777, 0.737)	(123, 141)	(45.42%, 100.0%)
StoIHT	0.102	0.217	0.132	(0.586, 0.557)	(0.586, 0.557)	(171, 179)	(9.48%, 45.60%)
GraphStoIHT	0.279	0.355	0.287	(0.669, 0.620)	(0.669, 0.620)	(150, 158)	(7.31%, 19.29%)
DA-IHT	0.679	0.741	0.694	(0.776, 0.733)	(0.776, 0.733)	(132, 141)	(4.86%, 42.86%)
GraphDA	0.855	0.870	0.850	(0.811, 0.799)	(0.811, 0.799)	(106, 107)	(4.55%, 43.89%)

Table 3. Table 3. Performance of Graph03

Method	${\text{Pre}}_{{𝒘}_t}$	${\text{Rec}}_{{𝒘}_t}$	${\text{F1}}_{{𝒘}_t}$	${\text{AUC}}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{Acc}}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{Miss}}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{NR}}_{{𝐰}_t,{\overline{𝐰}}_t}$
Adam	0.084	1.000	0.156	(0.820, 0.789)	(0.820, 0.788)	(104, 116)	(100.0%, 100.0%)
${\ell }_1$-RDA	0.340	0.940	0.488	(0.870, 0.833)	(0.869, 0.833)	( 88, 104)	(26.15%, 99.51%)
AdaGrad	0.318	0.942	0.462	(0.872, 0.825)	(0.872, 0.824)	( 88, 106)	(29.21%, 100.0%)
DA-GL	0.289	0.990	0.443	(0.894, 0.853)	(0.894, 0.853)	( 78, 100)	(32.07%, 100.0%)
DA-SGL	0.166	0.990	0.283	(0.883, 0.829)	(0.883, 0.828)	( 89, 111)	(52.19%, 100.0%)
StoIHT	0.156	0.208	0.175	(0.635, 0.593)	(0.634, 0.592)	(162, 173)	(11.62%, 50.60%)
GraphStoIHT	0.276	0.223	0.217	(0.666, 0.640)	(0.667, 0.640)	(154, 163)	(8.93%, 23.20%)
DA-IHT	0.716	0.782	0.734	(0.865, 0.834)	(0.865, 0.834)	( 95, 103)	(9.59%, 63.20%)
GraphDA	0.856	0.881	0.864	(0.898, 0.885)	(0.897, 0.885)	(72, 80)	(8.82%, 63.14%)

Table 4. Table 4. Performance of Graph04

Method	${\text{Pre}}_{{𝒘}_t}$	${\text{Rec}}_{{𝒘}_t}$	${\text{F1}}_{{𝒘}_t}$	${\text{AUC}}_{{𝐰}_t,{\overline{𝐰}}_t}$	$Ac{c}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{Miss}}_{{𝐰}_t,{\overline{𝐰}}_t}$	${\text{NR}}_{{𝐰}_t,{\overline{𝐰}}_t}$
Adam	0.121	1.000	0.216	(0.884, 0.858)	(0.884, 0.858)	( 77, 90)	(100.0%, 100.0%)
${\ell }_1$-RDA	0.361	0.961	0.513	(0.917, 0.896)	(0.917, 0.896)	( 66, 79)	(36.19%, 99.21%)
AdaGrad	0.376	0.961	0.528	(0.919, 0.889)	(0.919, 0.889)	( 67, 81)	(35.43%, 100.0%)
DA-GL	0.476	0.994	0.640	(0.942, 0.918)	(0.941, 0.918)	( 54, 73)	(26.03%, 100.0%)
DA-SGL	0.238	0.988	0.379	(0.931, 0.894)	(0.931, 0.894)	( 65, 85)	(53.31%, 100.0%)
StoIHT	0.207	0.203	0.204	(0.689, 0.639)	(0.689, 0.639)	(148, 160)	(11.86%, 47.88%)
GraphStoIHT	0.439	0.245	0.299	(0.743, 0.699)	(0.743, 0.699)	(131, 143)	(7.77%, 19.96%)
DA-IHT	0.780	0.801	0.788	(0.919, 0.898)	(0.919, 0.899)	( 74, 82)	(12.51%, 72.72%)
GraphDA	0.931	0.865	0.895	(0.939, 0.925)	(0.939, 0.925)	(56, 61)	(11.30%, 72.80%)

Equations

$$R(T,\mathcal{M}(\mathbb{M})):=\sum _{t=1}^T{f}_t({\mathbf{w}}_t,\{{\mathbf{x}}_t,y_t\})-\underset{\mathbf{w}\in \mathcal{M}(\mathbb{M})}{\min }\sum _{t=1}^T{f}_t(\mathbf{w}),$$

$$\mathcal{M}(\mathbb{M}):=\{\mathbf{w}|\text{supp}(\mathbf{w})\subseteq S\text{ for some }S\in \mathbb{M}\},$$

$$\mathbb{M}=$$

$${\mathbb{V}}_{\mathcal{F}}=F,\gamma (\mathcal{F})=g,\text{ and }\mathbf{c}(\mathcal{F})\le B\},$$

$${\mathbf{w}}_{t+1}=\mathrm{P}({\mathbf{w}}_t-{\eta }_t\nabla f_t({\mathbf{w}}_t),\mathcal{M}(\mathbb{M})),$$

$$\mathrm{P}(\mathbf{w},\mathcal{M}(\mathbb{M}))=\underset{\mathbf{x}\in \mathcal{M}(\mathbb{M})}{\mathrm{argmin}}\Vert \mathbf{w}-\mathbf{x}{\Vert }^2.$$

$${\mathbf{w}}_{t+1}$$

$$\underset{S\in \mathbb{M}}{\max }\Vert P(-\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma },S){\Vert }_2^2$$

$$\underset{S\in \mathbb{M}}{\min }\Vert -\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma }-P(-\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma },S){\Vert }_2^2,$$

$${\mathbf{w}}_{t+1}$$

$$\displaystyle=\operatorname*{arg\,min}_{{\bm{w}}\in{\mathcal{M}}(\mathbb{M})}\Big{\{}\frac{\sqrt{t}}{2\gamma}\|\bar{\bm{s}}_{t+1}\|_{2}^{2}+\langle\bar{\bm{s}}_{t+1},{\bm{w}}\rangle+\frac{\gamma}{2\sqrt{t}}\|{\bm{w}}\|_{2}^{2}\Big{\}}$$

$$\displaystyle=\operatorname*{arg\,min}_{{\bm{w}}\in{\mathcal{M}}(\mathbb{M})}\frac{\gamma}{2\sqrt{t}}\Big{\|}{\bm{w}}-\Big{(}-\frac{\sqrt{t}}{\gamma}\bar{\bm{s}}_{t+1}\Big{)}\Big{\|}_{2}^{2}$$

$$\displaystyle=\operatorname*{arg\,min}_{{\bm{w}}\in{\mathcal{M}}(\mathbb{M})}\Big{\|}{\bm{w}}-\Big{(}-\frac{\sqrt{t}}{\gamma}\bar{\bm{s}}_{t+1}\Big{)}\Big{\|}_{2}^{2},$$

$$\underset{\mathbf{w}\in \mathcal{M}(\mathbb{M})}{\min }\Vert -\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma }-\mathbf{w}{\Vert }_2^2\iff \underset{S\in \mathbb{M}}{\min }\Vert -\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma }-\mathrm{P}(-\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma },S){\Vert }_2^2.$$

$$\Vert \mathbf{x}{\Vert }_2^2-\Vert \mathrm{P}(\mathbf{x},S){\Vert }_2^2=\Vert \mathbf{x}-\mathrm{P}(\mathbf{x},S){\Vert }_2^2.$$

$$\displaystyle\min_{S\in\mathbb{M}}\Big{\{}\|-\frac{\sqrt{t}\bar{\bm{s}}_{t+1}}{\gamma}\|_{2}^{2}-\|{\rm P}({-\frac{\sqrt{t}\bar{\bm{s}}_{t+1}}{\gamma}},S)\|_{2}^{2}\Big{\}}$$

$$=\underset{S\in \mathbb{M}}{\min }\Vert -\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma }-\mathrm{P}(-\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma },S){\Vert }_2^2.$$

$$\Vert -\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma }{\Vert }_2^2+\underset{S\in \mathbb{M}}{\max }\Vert \mathrm{P}(-\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma },S){\Vert }_2^2$$

$$=\underset{S\in \mathbb{M}}{\min }\Vert -\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma }-\mathrm{P}(-\frac{\sqrt{t}{\bar{\mathbf{s}}}_{t+1}}{\gamma },S){\Vert }_2^2.$$

$$\Vert \mathrm{P}(\mathbf{w},H){\Vert }^2\ge c_H\cdot \underset{S\in \mathbb{M}}{\max }\Vert \mathrm{P}(\mathbf{w},S){\Vert }^2,$$

$$\Vert \mathrm{P}(\mathbf{w},T)-\mathbf{w}{\Vert }^2\le c_T\cdot \underset{S\in \mathbb{M}}{\min }\Vert \mathbf{w}-\mathrm{P}(\mathbf{w},S){\Vert }^2,$$

$$\Vert \mathrm{P}(\mathbf{w},H){\Vert }_2^2=\underset{S\in \mathbb{M}}{\max }\Vert \mathrm{P}(\mathbf{w},S){\Vert }_2^2,$$

$$\Vert \mathrm{P}(\mathbf{w},T)-\mathbf{w}{\Vert }_2^2=\underset{S\in \mathbb{M}}{\min }\Vert \mathbf{w}-\mathrm{P}(\mathbf{w},S){\Vert }_2^2.$$

$$|w_{{\tau }_1}|\ge |w_{{\tau }_2}|\ge \dots \ge |w_{{\tau }_s}|\ge \dots \ge |w_{{\tau }_p}|.$$

$$\Vert \mathrm{P}(\mathbf{w},H){\Vert }_2^2\ge \Vert \mathrm{P}(\mathbf{w},S){\Vert }_2^2.$$

$$\Vert \mathrm{P}(\mathbf{w},H){\Vert }_2^2\le \underset{S\in \mathbb{M}}{\max }\Vert \mathrm{P}(\mathbf{w},S){\Vert }_2^2.$$

$$R(T,\mathcal{M}(\mathbb{M}))\le 2DL\sqrt{T}.$$

$$\|{\bm{w}}_{T+1}-{\bm{w}}^{*}\|_{2}^{2}\leq 2\Big{(}D^{2}+\frac{L^{2}}{\gamma^{2}}-\frac{1}{\gamma\sqrt{T}}R(T,\mathcal{M}(\mathbb{M}))\Big{)}.$$

$$Pr{e}_{{\mathbf{w}}_t}$$

$${\text{F1}}_{{\mathbf{w}}_t}$$

$$f_t({\mathbf{w}}_t,\{{\mathbf{x}}_t,y_t\})=\log (1+\exp (-y_t\cdot {\mathbf{w}}_t^{\top }{\mathbf{x}}_t)).$$

$$(x_i{)}_{v_i}\sim \{\begin{array}{ll}\mathcal{N}(\mu ,1)& v_i\in F\\ \mathcal{N}(0,1)& v_i\ne F,\end{array}$$

$$f_t({\mathbf{w}}_t,\{{\mathbf{x}}_t,y_t\})=(y_t-⟨{\mathbf{w}}_t,{\mathbf{x}}_t⟩{)}^2.$$

$$y_t={\mathbf{x}}_t^{\top }{\mathbf{w}}^{\ast },$$

Full text

online learning; dual averaging; graph-structured sparsity

††journalyear: 2019††copyright: acmcopyright††conference: The 25th ACM SIGKDD Conference on Knowledge Discovery and Data Mining; August 4–8, 2019; Anchorage, AK, USA††booktitle: The 25th ACM SIGKDD Conference on Knowledge Discovery and Data Mining (KDD ’19), August 4–8, 2019, Anchorage, AK, USA††price: $15.00††doi: 10.1145/3292500.3330915††isbn: 978-1-4503-6201-6/19/08

1. Introduction

As a new paradigm in machine learning, convex online learning algorithms have received enormous attention (ying2006online; Duchi et al., 2011; Xiao, 2010; Hazan et al., 2016; Shalev-Shwartz et al., 2012; Kingma and Ba, 2014; ying2017unregularized). These algorithms update learning models sequentially by using one training sample at each iteration, which makes them applicable to large-scale datasets on the fly and still enjoy non-regret property. For better interpretability and less computational complexity in high dimension data, many online learning algorithms (Xiao, 2010; Yang et al., 2010; Langford et al., 2009; Duchi et al., 2011) exploit $\ell_{1}$ norm or $\ell_{1}/\ell_{2}$ mixed norm to achieve sparse solution (Xiao, 2010; Yang et al., 2010; Jacob et al., 2009). However, these sparsity-inducing models cannot characterize more complex (usually non-convex) graph-structured sparsity constraint, hence, unable to use some important priors such as graph data.

Graph-structured sparsity models have significant real-world applications, for example, social events (Rozenshtein et al., 2014), disease outbreaks (Qian et al., 2014), computer viruses (Draief et al., 2006), and gene networks (Chuang et al., 2007). These applications all contain graph structure information and the data samples are usually collected on the fly, i.e., the training samples have been received and processed one by one. Unfortunately, most of the graph-structured (non-convex) methods (Hegde et al., 2015b, 2016; Chen and Zhou, 2016; Aksoylar et al., 2017) are batch learning-based, which cannot be applied to the online setting. The past few years have seen a surge of convex online learning algorithms, such as online projected gradient descent (Zinkevich, 2003), AdaGrad (Duchi et al., 2011), Adam (Kingma and Ba, 2014), $\ell_{1}$-RDA (Xiao, 2010), FOBOS (Duchi and Singer, 2009), and many others (e.g. (Hazan et al., 2016; Shalev-Shwartz et al., 2012)). However, they cannot be used to tackle online graph-structured sparsity problems due to the limitation of sparsity-inducing norms.

In recent years, machine learning community (Lafond et al., 2015; Hazan et al., 2017; Chen et al., 2018; Gonen and Hazan, 2018; Gao et al., 2018; Yang et al., 2018) have made promising progress on online non-convex optimization with regards to algorithms and local-regret bounds. Nonetheless, these algorithms cannot deal with graph-structured sparsity constraint problems due to the following two limitations: 1) The existing non-convexity assumption is only on the loss functions subject to a convex constraint; 2) Most of these proposed algorithms are based on online projected gradient descent (PGD), and cannot explore the structure information, hardly workable for graph-structured sparsity constraint. To the best of our knowledge, there is no existing work to tackle the combinatorial non-convexity constraint problems under online setting.

In this paper, we aim to design an approximated online learning algorithm that can capture graph-structured information effectively and efficiently. To address this new and challenging question, the potential algorithm has to meet two crucial requirements: 1) graph-structured: The algorithm should effectively capture the latent graph-structured information such as trees, clusters, connected subgraphs; 2) online: The algorithm should be efficiently applicable to online setting where training samples can only be processed one by one. Our assumption on the problem has a non-convex constraint but with a convex objective, which will sustain higher applicability in the practice of our setting. Inspired by the success of dual-averaging (Nesterov, 2009; Xiao, 2010), we propose the Graph Dual Averaging Algorithm, namely, GraphDA. The key part in GraphDA is to keep track of both averaging gradient via dual variables in dual space and primal variables in primal space. We then use two approximated projections to project both primal variables and dual variables onto low dimension subspaces at each iteration. We conduct extensive experiments to demonstrate that by projecting both primal and dual variables, GraphDA captures the graph-structured sparsity effectively. Overall, our contributions are as follows:

$\bullet$ We propose a dual averaging-based algorithm to solve graph-structured sparsity constraint problems under online setting. To the best of our knowledge, it is a first attempt to establish an online learning algorithm for the graph-structured sparsity model.

$\bullet$ We prove the minimization problem occurring at each dual averaging step, which can be formulated as two equivalent optimization problems: minimization problem in primal space and maximization problem in dual space. The two optimization problems can then be solved approximately by adopting two popular projections. Furthermore, we provide two exact projection algorithms for the non-graph data.

$\bullet$ We conduct extensive experiments on both synthetic and real-world graphs. The experimental results demonstrate that GraphDA can successfully capture the latent graph-structure during online learning process. The learned model generated by our algorithm not only achieves higher classification accuracy but also stronger interpretability compared with the state-of-the-art algorithms.

The rest of the paper is organized as follows: Related work is teased out in Section 2. Section 3 gives the notations and problem definition. In Section 4, we present our main idea and algorithms. We report and discuss the experiment results in comparison with other baseline methods in Section 5. A short conclusion ensues in Section 6. Due to space limit, the detailed experimental setup and partial experimental results are supplied in Appendix. Our source code including baseline methods and datasets are accessible at: https://github.com/baojianzhou/graph-da.

2. Related Work

In line with the focus of the present work, we categorize highly related researches into three sub-topics for the sake of clarity.

**Online learning with sparsity. ** Online learning algorithms (Zinkevich, 2003; Bottou and Cun, 2004; Ying and Pontil, 2008; Hazan et al., 2016; Shalev-Shwartz et al., 2012) try to solve classification or regression problems that can be employed in a fully incremental fashion. A natural way to solve online learning problem is to use stochastic gradient descent by using one sample at a time. However, this type of methods usually cannot produce any sparse solution. The gradient of only one sample has such a large variance that renders its projection unreliable. To capture the model sparsity, $\ell_{1}$ norm-based (Bottou, 1998; Duchi and Singer, 2009; Duchi et al., 2011; Langford et al., 2009; Xiao, 2010) and $\ell_{1}/\ell_{2}$ mixed norm-based (Yang et al., 2010) are used under online learning setting; the dual-averaging (Xiao, 2010) adds a convex regularization, namely $\ell_{1}$-RDA to learn a sparsity model. Based on the dual averaging work, online group lasso and overlapping group lasso are proposed in (Yang et al., 2010), which provides us a sparse solution. However, the solution cannot produce methods directly applicable to graph-structured data. For example, as pointed out by (Xiao, 2010), the levels of sparsity proposed in (Duchi and Singer, 2009; Langford et al., 2009) are not satisfactory compared with their batch counterparts.

**Model-based sparsity. ** Different from $\ell_{1}$-regularization (Tibshirani, 1996) or $\ell_{1}$-ball constraint-based method (Duchi et al., 2008), model-based sparsity are non-convex (Baraniuk et al., 2010; Hegde et al., 2015a, b, 2016). Using non-convex such as $\ell_{0}$ sparsity based methods (Yuan et al., 2014; Bahmani et al., 2013; Zhou et al., 2018; Nguyen et al., 2017) becomes popular, where the objective function is assumed to be convex with a sparsity constraint. To capture graph-structured sparsity constraint such as trees and connected graphs, a series of work  (Baraniuk et al., 2010; Hegde et al., 2014b, 2016, 2015b) has proposed to use structured sparsity model $\mathbb{M}$ to define allowed supports $\mathbb{M}=\{S_{1},S_{2},\ldots,S_{k}\}$. These complex models are non-convex, and gradient descent-based algorithms involve a projection operator which is usually NP-hard. (Hegde et al., 2015b, 2016, a) use two approximated projections (head and tail) without sacrificing too much precision. However, the above research work cannot be directly applied to online setting and the objective function considered is not general loss.

**Online non-convex optimization. ** The basic assumption in recent progress on online non-convex optimization (Lafond et al., 2015; Hazan et al., 2017; Chen et al., 2018; Gonen and Hazan, 2018; Gao et al., 2018; Yang et al., 2018) is that the objective considered is non-convex. Local-regret bound has been explored in these studies, most of which are based on projected gradient descent methods, for example, online projected gradient descent (Hazan et al., 2017) and online normalized gradient descent (Gao et al., 2018). However, these online non-convex algorithms cannot deal with our problem setting where there exists a combinatorial non-convex structure.

3. Preliminaries

We begin by introducing some basic mathematical terms and notations, and then define our problem.

3.1. Notations

An index set is defined as $[p]=\{1,\ldots,p\}$. The bolded lower-case letters, e.g., ${\bm{w}},{\bm{x}}\in\mathbb{R}^{p}$, denote column vectors where their $i$-th entries are $w_{i}$, $x_{i}$. The $\ell_{2}$-norm of ${\bm{w}}$ is denoted as $\|{\bm{w}}\|_{2}$. The inner product of ${\bm{x}}$ and ${\bm{y}}$ on $\mathbb{R}^{p}$ is defined as $\langle{\bm{x}},{\bm{y}}\rangle=x_{1}y_{1}+\cdots+x_{p}y_{p}$. Given a differentiable function $f({\bm{w}}):\mathbb{R}^{p}\rightarrow\mathbb{R}$, the gradient at ${\bm{w}}$ is denoted as $\nabla f({\bm{w}})$. The support set of ${\bm{w}}$, i.e., $\text{supp}({\bm{w}}):=\{i|w_{i}\neq 0\}$, is defined as a subset of indices which index non-zero entries. If $|\text{supp}({\bm{w}})|\leq s$, ${\bm{w}}$ is called an $s$ sparse vector. The upper-case letters, e.g., $\Omega$, denote a subset of $[p]$ and its complement is $\Omega^{c}=[p]\backslash\Omega$. The restricted vector of ${\bm{w}}$ on $\Omega$ is denoted as ${\bm{w}}_{\Omega}\in\mathbb{R}^{p}$, where $({\bm{w}}_{\Omega})_{i}=w_{i}$ if $i\in\Omega$; otherwise 0. We define the undirected graph as $\mathbb{G}(\mathbb{V},\mathbb{E})$, where $\mathbb{V}=[p]$ is the set of nodes and $\mathbb{E}$ is the set of edges such that $\mathbb{E}\subseteq\{(u,v)|u\in\mathbb{V},v\in\mathbb{V}\}$. The upper-case letters, e.g., $H,T,S$, stand for subsets of $[p]:=\{1,2,\ldots,p\}$. Given the standard basis $\{{\bm{e}}_{i}:1\leq i\leq p\}$ of $\mathbb{R}^{p}$, we also use $H,T,S$ to represent subspaces. For example, the subspace $S$ is the subspace spanned by $S$, i.e., $\text{span}\{{\bm{e}}_{i}:i\in S\}$. We will clarify the difference only if confusion occurs.

3.2. Problem Definition

Here, we study an online non-convex optimization problem, which is to minimize the regret as defined in the following:

[TABLE]

where each $f_{t}({\bm{w}}_{t},\{{\bm{x}}_{t},y_{t}\})$ is the loss that a learner predicts an answer for the question ${\bm{x}}_{t}$ after receiving the correct answer $y_{t}$, and $\min_{{\bm{w}}\in\mathcal{M}(\mathbb{M})}\sum_{t=1}^{T}f_{t}({\bm{w}})$ is the minimum loss that the learner can potentially get. To simplify, we assume $f_{t}$ is convex differentiable. For example, we can use the least square loss $f_{t}({\bm{w}}_{t},\{{\bm{x}}_{t},y_{t}\})=({\bm{w}}_{t}^{\top}{\bm{x}}_{t}-y_{t})^{2}$ for the online linear regression problem and logistic loss $f_{t}({\bm{w}}_{t},\{{\bm{x}}_{t},y_{t}\})$ $=\log(1+\exp(-y_{t}\cdot{\bm{w}}^{\top}{\bm{x}}_{t}))$ for online binary classification problem where $y_{t}\in\{\pm 1\}$. The goal of the learner is to minimize the regret $R(T,\mathcal{M}(\mathbb{M}))$. Different from the online convex optimization setting in (Shalev-Shwartz et al., 2012; Hazan et al., 2016), $\mathcal{M}(\mathbb{M})\subseteq\mathbb{R}^{p}$ is a generally non-convex set. To capture more complex graph-structured information, in a series of seminal work (Baraniuk et al., 2010; Hegde et al., 2015a, b), a structured sparsity model $\mathcal{M}(\mathbb{M})$ is proposed as follows:

[TABLE]

where $\mathbb{M}=\{S_{1},S_{2},\ldots,S_{k}\}$ is the collection of allowed structure supports with $S_{i}\in[p]$. Basically, $\mathcal{M}(\mathbb{M})$ is the union of $k$ subspaces. Each subspace is uniquely identified by $S_{i}$. Definition (2) is so general that it captures a broad spectrum of graph-structured sparsity models such as trees (Hegde et al., 2014b), connected subgraphs (Hegde et al., 2015b; Chen and Zhou, 2016). We mainly focus on the Weighted Graph Model(WGM) proposed in (Hegde et al., 2015b).

Definition 0 (Weighted Graph Model (Hegde

et al., 2015b)).

Given an underlying graph $\mathbb{G}=(\mathbb{V},\mathbb{E},{\bm{c}})$ defined on the coefficients of the unknown vector ${\bm{w}}$, where $\mathbb{V}=[p]$, $\mathbb{E}\subseteq\mathbb{V}\times\mathbb{V}$ and associated cost vector ${\bm{c}}$ on edges, then the weighted graph model $(\mathbb{G},s,g,B)$-WGM can be defined as the following set of supports:

[TABLE]

where $B$ is the budget on cost of edges in forest $\mathcal{F}$, $\gamma(\mathcal{F})$ is the number of connected component in forest $\mathcal{F}$ denoted as $g$, and $s$ is the sparsity. To clarify, forest $\mathcal{F}$ is the subgraph induced by its nodes set $F$, i.e. $\mathcal{F}:=\mathbb{G}(F,\mathbb{E}^{\prime})$, where $\mathbb{E}^{\prime}=\{(u,v):u\in F,v\in F,(u,v)\in\mathbb{E}\}$. ${\bm{c}}(F)$ is the total edge costs in forest $\mathcal{F}$.

$(\mathbb{G},s,g,B)$-WGM captures a broad range of graph structures such as groups, clusters, trees, and subgraphs. A toy example is given in Figure 1 where we define a graph with 6 nodes in $\mathbb{V}$, 7 edges in $\mathbb{E}$, and the cost of all edges is set to $1$. Let $\mathbb{V}$ be associated with a vector ${\bm{w}}\in\mathbb{R}^{6}$. Suppose we are interested in connected subgraphs111A connected subgraph is the subgraph which has only 1 connected component. with at most 3 nodes, to capture these subgraphs, $\mathbb{M}$ can be defined as $\mathbb{M}=\{S_{i}|\mathbb{G}(S_{i},\mathbb{E}_{i})\text{ is connected },|S_{i}|\leq 3\}$. By letting the budget $B=3$, and the sparsity parameter $s=3$, we can clearly use $(\mathbb{G},3,1,3)$-WGM to represent this $\mathbb{M}$. Figure 1 shows three subgraphs formed by red nodes and edges. The subgraph induced by $\{w_{1},w_{2},w_{5}\}$ on the left is in $\mathbb{M}$. However, the subgraph induced by $\{w_{2},w_{5},w_{6}\}$ in the middle is not in $\mathbb{M}$ because of the non-connectivity. The subgraph formed by $\{w_{3},w_{4},w_{5},w_{6}\}$ on the right is not in $\mathbb{M}$ either, as it violates the sparsity constraint, i.e., $s\leq 3$.

After defining the structure-sparsity model $\mathcal{M}(\mathbb{M})$, we explore how to design an efficient and effective algorithm to minimize the regret under model constraint. An intuitive way to do this is to use online projected gradient descent (Zinkevich, 2003) where the algorithm needs to solve the following projection at iteration $t$:

[TABLE]

where $\eta_{t}$ is the learning rate and ${\rm P}$ is the projection operator onto $\mathcal{M}(\mathbb{M})$, i.e., ${\rm P}(\cdot,\mathcal{M}(\mathbb{M})):$ $\mathbb{R}^{p}\rightarrow\mathbb{R}^{p}$ is defined as

[TABLE]

However, there are two essential drawbacks of using (3): First, the projection in (3) only uses single gradient $\nabla f_{t}({\bm{w}}_{t})$ which is too noisy (large variance) to capture the graph-structured information at each iteration; Second, the training samples coming later are less important than these coming earlier due to the decay of learning rate $\eta_{t}$. Recall that $\eta_{t}$ needs to decay asymptotically to $\mathcal{O}(1/\sqrt{t})$ in order to achieve a non-regret bound. Fortunately, inspired by (Nesterov, 2009; Xiao, 2010), the above two weaknesses can be successfully overcome by using dual averaging. The main idea is to keep tracking both primal vectors (corresponding to ${\bm{w}}_{t}$ in primal space) and dual variables (corresponding to gradients, $\nabla f_{t}({\bm{w}}_{t})$ in dual space222Notice that we use $\ell_{2}$-norm, i.e. $\|\cdot\|_{2}$, which is defined in the Euclidean space $X=\mathbb{R}^{p}$. By definition, the dual norm of $\ell_{2}$ is identical to itself, i.e., $\|\cdot\|_{*}=\|\cdot\|_{2}$. Also, recall that the dual space $X^{*}$ of the Euclidean space is also identical with each other $X^{*}=X=\mathbb{R}^{p}$.) at each iteration. In the next section, we focus on designing an efficient algorithm by using the idea of dual averaging to capture graph-structured sparsity under online setting.

4. Algorithm: GraphDA

We try to develop a dual averaging-based method to minimize the regret (1). At each iteration, the method updates ${\bm{w}}_{t}$ by using the following minimization step:

[TABLE]

where $\beta_{t}$ is to control the learning rate implicitly and ${\bm{g}}_{i}$ is a subgradient in $\partial f_{i}({\bm{w}},\{{\bm{x}}_{i},y_{i}\})=\{{\bm{g}}:f_{i}({\bm{z}},\{{\bm{x}}_{i},y_{i}\})\geq f_{i}({\bm{w}},\{{\bm{x}}_{i},y_{i}\})+\langle{\bm{g}},{\bm{z}}-{\bm{w}}\rangle,\forall{\bm{z}}\in\mathcal{M}(\mathbb{R})\}$333As we assume $f_{i}$ is convex differentiable, then we have $\partial f_{i}({\bm{w}},\{{\bm{x}}_{i},y_{i}\})=\{\nabla f_{i}({\bm{w}},\{{\bm{x}}_{i},y_{i}\})\}$, i.e. ${\bm{g}}_{i}=\nabla f_{i}({\bm{w}},\{{\bm{x}}_{i},y_{i}\}).$. Different from the convexity explored in (Nesterov, 2009) and (Xiao, 2010), the problem considered here is generally non-convex, which makes it NP-hard to solve. Initially, the solution of the primal is set to zero, i.e., ${\bm{w}}_{0}={\bm{0}}$.444There are two advantages: 1. $\bm{w}_{0}=\bm{0}$ is trivially in the $\mathcal{M}(\mathbb{M})$; 2. $\bm{w}_{0}=\bm{0}\in\operatorname*{arg\,min}_{\bm{w}}\|\bm{w}\|_{2}^{2}$ under convex setting (Xiao, 2010) so that sublinear regret can obtain. Then at each iteration, it computes a subgradient ${\bm{g}}_{t}$ based on current data sample $\{{\bm{x}}_{t},y_{t}\}$ and then updates ${\bm{w}}_{t}$ by using (5) averaging gradient from the dual space. The algorithm terminates after receiving $T$ samples and returns the model ${\bm{w}}_{T}$ or $\bar{{\bm{w}}}_{T}=1/T\sum_{t=0}^{T}{\bm{w}}_{t}$ depending on needs. The dual averaging step (5) has two advantages: 1) The gradient information of training samples coming later will not decay when new samples are coming; 2) The averaging gradient can be accumulated during the learning process; hence we can use it to capture graph-structure information more effectively than online PGD-based methods.

Due to the NP-hardness to compute (5), it is impractical to directly use (5). Thus, we have to treat this minimization step more carefully for $\mathcal{M}(\mathbb{M})$. The minimization step (5) has the following equivalent projection problems, specified in Theorem 1.

Theorem 1.

Assume $\beta_{t}=\gamma\sqrt{t}$, where $\gamma>0$ and denote $\bar{\bm{s}}_{t+1}=\frac{1}{t+1}\sum_{i=0}^{t}{\bm{g}}_{i}$. The minimization step of (5) can be expressed as the following two equivalent optimization problems:

[TABLE]

where $P({\bm{s}},S)$ is the projection operator that projects ${\bm{s}}$ onto the subspace spanned by $S$.

Proof.

The original minimization problem in (5) can be equivalently expressed as

[TABLE]

where the second equality follows by adding a constant to the minimization objective and  (8) follows by multiplying $2\sqrt{t}/\gamma$ on the third equation. Hence, (5) is equivalent to the minimization of (8). Clearly, (8) is essentially the projection ${\rm P}(-(\sqrt{t}{\bm{s}}_{t+1})/\gamma,\mathcal{M}(\mathbb{M}))$ defined in (4). To further explore (8), notice that one needs to solve the following equivalent minimization problem:

[TABLE]

Here, for any ${\bm{x}}$, ${\rm P}({\bm{x}},S)$ is an orthogonal projection operator that projects ${\bm{x}}$ onto subspace spanned by $S$. By the projection theorem, for any ${\bm{x}}$, it always has the following property:

[TABLE]

Replacing ${\bm{x}}$ by $-\sqrt{t}\bar{\bm{s}}_{t+1}/\gamma$ and adding minimization to both sides with respect to subspace $S$, we obtain:

[TABLE]

By moving the minimization into the negative term, we obtain

[TABLE]

We prove the theorem. ∎

The above theorem leads to a key insight that the NP-hard problem (8) can be solved either by maximizing $\|{\rm P}(\sqrt{t}\bar{\bm{s}}_{t+1}/{\gamma},S)\|_{2}^{2}$ or by minimizing $\|\sqrt{t}\bar{\bm{s}}_{t+1}/{\gamma}-{\rm P}(\sqrt{t}\bar{\bm{s}}_{t+1}/{\gamma},S)\|_{2}^{2}$ over $S$. Inspired by (Hegde et al., 2016, 2015a, 2015b), instead of solving these two problems exactly, we apply two approximated algorithms provided in (Hegde et al., 2015b) to solve the problem approximately. We present the following two assumptions:

Assumption 1 (Head Projection (Hegde

et al., 2016)).

Let $\mathbb{M}$ and $\mathbb{M}_{H}$ be the predefined subspace models. Given any ${\bm{w}}$, there exists a $(c_{H},\mathbb{M},\mathbb{M}_{H})$ Head-Projection which is to find a subspace $H\in\mathbb{M}_{H}$ such that

[TABLE]

where $0<c_{H}\leq 1$. We denote ${\rm P}({\bm{w}},H)$ as ${\rm P}({\bm{w}},\mathbb{M},\mathbb{M}_{H})$.

Assumption 2 (Tail Projection (Hegde

et al., 2016)).

Let $\mathbb{M}$ and $\mathbb{M}_{T}$ be the predefined subspace models. Given any ${\bm{w}}$, there exists a $(c_{T},\mathbb{M},\mathbb{M}_{T})$ Tail-Projection which is to find a subspace $T\in\mathbb{M}_{T}$ such that

[TABLE]

where $c_{T}\geq 1$. We denote ${\rm P}({\bm{w}},T)$ as ${\rm P}({\bm{w}},\mathbb{M},\mathbb{M}_{T})$.

To minimize the regret $R(T,\mathcal{M}(\mathbb{M}))$, we propose the approximated algorithm, presented in Algorithm 1 below. Initially, the primal vector ${\bm{w}}_{0}$ and dual vector $\bar{\bm{s}}_{0}$ are all set to ${\bm{0}}$. At each iteration, it works as the following four steps:

$\bullet$ Step 1: The learner receives a question ${\bm{x}}_{t}$ and makes a prediction based on ${\bm{x}}_{t}$ and ${\bm{w}}_{t}$. After suffering a loss $f_{t}({\bm{w}}_{t},\{{\bm{x}}_{t},y_{t}\})$, it computes the gradient ${\bm{g}}_{t}$ in Line 4;

$\bullet$ Step 2: In Line 5, the current gradient ${\bm{g}}_{t}$ has been accumulated into $\bar{\bm{s}}_{t+1}$, which is ready for the next head projection555Pseudo-code of these two projections are provided in Appendix A for completeness.;

$\bullet$ Step 3: The head projection inputs the accumulated gradient $\bar{\bm{s}}_{t+1}$ and outputs the vector ${\bm{b}}_{t+1}$ so that $\text{supp}({\bm{b}}_{t+1})\in\mathbb{M}_{H}$;

$\bullet$ Step 4: The next predictor ${\bm{w}}_{t+1}$ is then updated by using the tail projection, i.e., $\text{supp}({\bm{w}}_{t+1})\in\mathbb{M}_{T}$. The weight $-\sqrt{t}/\gamma$ is to control the learning rate.

The algorithm repeats the above four steps until some stop condition is satisfied. The main difference between our method and the methods in (Nesterov, 2009; Xiao, 2010) lies in that, we, as a first attempt, use two projections (Line 6 and Line 7), to project dual vector $\bar{\bm{s}}_{t+1}$ and primal vector ${\bm{w}}_{t+1}$ onto a graph-structured subspaces $\mathbb{M}_{H}$ and $\mathbb{M}_{T}$ respectively. In dual projection step, most of the irrelevant gradient entries have been effectively set to zero values. In primal tail projection step, we make sure ${\bm{w}}_{t+1}$ has been projected onto $\mathbb{M}_{T}$ so that the constraint of interest is satisfied.

In real applications, graph data is not always available, i.e., $\mathbb{M}$ cannot be explicitly constructed by $\mathbb{G}(\mathbb{V},\mathbb{E})$, so we often have to deal with non-graph data but still with the aim to pursue structure sparsity constraint. To compensate, we provide Dual Averaging Iterative Hard Thresholding, namely DA-IHT, presented in Theorem 13, to handle non-graph data cases.

Theorem 2.

Assume that the graph information is not available or the graph is a complete graph and the budget $B$ is large enough. We can define our model $\mathbb{M}$ such that it includes all possible $s$-sparse subgraphs, i.e., $\mathbb{M}=\{S:|S|\leq s\}$. Then there exists exactly head and tail projection algorithm such that

[TABLE]

and

[TABLE]

Proof.

Since the graph is a complete graph (i.e., all subgraphs are connected.) and the budget constraint $B$ is large enough, any subset $S$ that has $s$ elements belongs to $\mathbb{M}$. In this case, $\mathbb{M}$ contains all $s$-subsets, i.e., $\mathbb{M}=\{S_{i}:|S_{i}|\leq s\}$. By sorting the magnitudes of ${\bm{w}}$ in a descending manner, we have

[TABLE]

Let $H=T=\{\tau_{1},\tau_{2},\ldots,\tau_{s}\}$. For any s-sparse set $S$, by the fact that $|w_{\tau_{1}}|,|w_{\tau_{2}}|,\ldots,|w_{\tau_{s}}|$ are the largest magnitude $s$ entries, we always have

[TABLE]

At the same time, $H\in\mathbb{M}$, then

[TABLE]

Hence, we prove (12). In a similar vein, we can also prove (13). ∎

By Theorem 13, one can implement the two projections in Line 6 and 7 of Algorithm 1 by sorting the magnitudes $\bar{\bm{s}}_{t+1}$ and $-\sqrt{t}/\gamma{\bm{b}}_{t+1}$ respectively, to deal with non-graph data. DA-IHT will be used as a baseline in our experiment to compare with the graph-based method, GraphDA.

**Time Complexity. ** At each iteration of GraphDA, the time complexity of two projections depends on the graph size $p$ and the number of edges $|\mathbb{E}|$. As proved in (Hegde et al., 2015b), two projections have the time complexity $\mathcal{O}(|\mathbb{E}|\log^{3}(p))$. In many real-world applications, the graphs are usually sparse, i.e., $\mathcal{O}(p)$, and then the total complexity of each iteration of GraphDA is $\mathcal{O}(p+p\log^{3}(p))$. Our method is characterized by two merits: 1) The time cost of each iteration is nearly-linear time; 2) At each iteration, it only has $\mathcal{O}(p+|\mathbb{E}|)$ memory cost, where $\mathcal{O}(p)$ stores the averaging gradient and current solution and $\mathcal{O}(|\mathbb{E}|)$ is to save the graph. For DA-IHT, we need to select the top $s$ largest magnitude entries at each iteration. Thus, the time cost of per-iteration is $\mathcal{O}(sp)$ with $\mathcal{O}(p)$ memory cost.

**Regret Discussion. ** Given any online learning algorithm, we are interested in whether the regret $R(T,\mathcal{M}(\mathbb{M}))$ is sub-linear and whether we can bound the estimation error $\|{\bm{w}}_{t}-{\bm{w}}^{*}\|_{2}$. We first assume the primal vectors are $\{{\bm{w}}_{t}\}_{t=0}^{T}$ and the dual gradient sequences $\{{\bm{g}}_{t}\}_{t=0}^{T}$. We then assume that the potential solution ${\bm{w}}$ is always bounded in $D$, i.e., $\|{\bm{w}}\|_{2}\leq D$ and gradients are also bounded, i.e., $\|{\bm{g}}_{t}\|_{2}\leq L$. Then for any $T\geq 1$ and any ${\bm{w}}\in\mathcal{M}(\mathbb{M})$, the regret in (Xiao, 2010) can be bounded as the following:

[TABLE]

Given any optimal solution ${\bm{w}}^{*}\in\operatorname*{arg\,min}_{{\bm{w}}\in\mathcal{M}(\mathbb{M})}\sum_{i=1}^{T}f_{i}({\bm{w}})$ and the solution ${\bm{w}}_{T}$, the estimation error, i.e., $\|{\bm{w}}_{T+1}-{\bm{w}}^{*}\|$ is bounded as the following:

[TABLE]

However, the regret bound (14) and estimation error (15) are under the assumption that the constraint set $\mathcal{M}(\mathbb{M})$ is convex. For GraphDA, an approximated algorithm, it is difficult to establish a sublinear regret bound. The reasons are two-fold: 1) Due to the non-convexity of $\mathcal{M}(\mathbb{M})$, it is possible that GraphDA converges to a local minimal, so the regret will potentially be non-sublinear; 2) The solution of model projection is approximated, making the regret analysis harder. Although recent work (Gao et al., 2018; Hazan et al., 2017) shows that it is still possible to obtain a local-regret bound when the objective function is non-convex, it is different from our case since we assume the objective function convex subject to a non-convex constraint. We leave the theoretical regret bound analysis of GraphDA an open problem.

5. Experiments

To corroborate our algorithm, we conduct extensive experiments, comparing GraphDA with some popular baseline methods. Note DA-IHT derived from Theorem 13 is treated as a baseline method. We aim to answer the following questions:

• •

**Question Q1: ** Can GraphDA achieve better classification performance compared with baseline methods?

• •

**Question Q2: ** Can GraphDA learn an stronger interpretative model through capturing more meaningful graph-structure features compared with baseline methods?

5.1. Datasets and evaluation metrics

**Datasets. ** We use the following three publicly available graph datasets: 1) Benchmark Dataset (Arias-Castro et al., 2011). Four benchmark graphs (Arias-Castro et al., 2011) are shown in Figure 2. The four subgraphs are embedded into $33\times 33$ graphs with 26, 46, 92, and 132 nodes respectively. Each graph has $p=1,089$ nodes and $m=2,112$ edges with unit weight $1.0$. We use the Benchmark dataset to learn an online graph logistic regression model; 2) MNIST Dataset (LeCun, 1998). This popular hand-writing dataset is used to test GraphDA on online graph sparse linear regression. It contains ten classes of handwritten digits from 0 to 9. We randomly choose each digit as our target graph. Each pixel stands for a node. There exists an edge if two nodes are neighbors. We set the weights to 1.0 for edges; 3) KEGG Pathway Dataset (LeNail et al., 2017). The Kyoto Encyclopedia of Genes and Genomes (KEGG) dataset contains 5,372 genes. These genes (nodes) form a connected graph with 78,545 edges. The edge weights stand for correlations between two genes. We use KEGG to detect a related pathway.

**Evaluation metrics. ** We have two categories of metrics to answer Question Q1 and Question Q2 respectively. To measure classification performance of ${\bm{w}}_{t}$ or $\bar{\bm{w}}_{t}$666For the comparison, we also evaluate the averaged decision $\bar{\bm{w}}_{t}$ similar as done in (Xiao, 2010)., we use classification Accuracy(Acc), the Area Under Curve (AUC) (Hanley and McNeil, 1982), and the number of Misclassified samples (Miss). To evaluate feature-level performance (interpretability), we use Precision (Pre), Recall (Rec), F1-score (F1), and Nonzero Ratio (NR). To clarify, given any optimal ${\bm{w}}^{*}\in\mathbb{R}^{p}$ and learned model ${\bm{w}}_{t}$, Pre, Rec, F1, and NR are defined as follows:

[TABLE]

5.2. Baseline methods

We consider the following eight baseline methods: 1) $\ell_{1}$-RDA (Xiao, 2010). We use the enhanced Regularized Dual-Averaging ($\ell_{1}$-RDA) method in Algorithm 2 of (Xiao, 2010); 2) DA-GL (Yang et al., 2010). Online Dual Averaging Group Lasso (DA-GL) is the dual averaging method with group Lasso; 3) DA-SGL (Yang et al., 2010). It also uses dual averaging, but with sparse group Lasso as the regularization; 4) AdaGrad (Duchi et al., 2011). The adaptive gradient with $\ell_{1}$ regularization is different from $\ell_{1}$-RDA (Xiao, 2010). AdaGrad yields a dedicated step size for each feature inversely. In order to capture the sparsity, we use its $\ell_{1}$ norm-based method for comparison; 5) Adam (Kingma and Ba, 2014). Since there is no sparsity regularization in Adam, it generates totally dense models. We use its online version777One can find more details of the online version in Section 4 of (Kingma and Ba, 2014). to compare with these sparse methods; 6) StoIHT (Nguyen et al., 2017). We use this method with block size 1, which can be treated as online learning setting; 7) DA-IHT, derived from Theorem 13 in this paper. We use it to compare with GraphDA, which has graph-structure constraint; 8) GraphStoIHT. We apply the head and tail projection to StoIHT to generate GraphStoIHT.

Online Setting. All methods are completely online, i.e., all learning algorithms receive a single sample per-iteration. Due to the space limit, the parameters of all baseline methods including GraphDA are in Appendix A. All numerical results are averaged from 20 trials. The following three sections report and discuss the experimental results on each dataset to answer Q1 and Q2.

5.3. Results from Benchmark dataset

Given the training dataset $\{{\bm{x}}_{i},y_{i}\}_{i=1}^{t}$, where ${\bm{x}}_{i}\in\mathbb{R}^{p}$ and ${y_{i}}\in\{\pm 1\}$ on the fly, the online graph sparse logistic regression is to minimize the regret $R({\bm{t}},\mathcal{M}(\mathbb{M}))$ where $f_{t}({\bm{w}}_{t})$ is a logistic loss defined as

[TABLE]

We simulate the negative and positive samples as done in (Arias-Castro et al., 2011). $y_{t}=-1$ stands for no signals or events (“business-as-usual”). $y_{t}=+1$ means a certain event happens such as disease outbreak/computer virus hidden in current data sample ${\bm{x}}_{t}$, and feature values in subgraphs are abnormally higher. That is, if $y_{t}=-1$, then $x_{v_{i}}\sim\mathcal{N}(0,1)\ \forall\ v_{i}\in\mathbb{V}$; and if $y_{t}=+1$, then,

[TABLE]

where $F$ stands for the nodes of a specific subgraph showcased in Figure 2. Then each entry $({\bm{w}}^{*})_{i}$ is $\mu$ if $i\in F$; otherwise 0. We first fix $\mu=0.3$ and then generate validating, training, and testing samples, each with 400 samples. All methods stop at $t=400$ after seeing all training samples once. Parameters are tuned on 400 validating samples. We test ${\bm{w}}_{t}$ and $\bar{\bm{w}}_{t}$ on testing samples.

**Classification Performance on fixed $\mu$. ** Table 1 shows that four all three indicators of classification performance, GraphDA scores higher than the other baseline methods. Specifically, it has the highest Acc (0.749, 0.739) and AUC (0.749, 0.739) with respect to ${\bm{w}}_{t}$ and $\bar{{\bm{w}}}_{t}$. The averaged number of misclassified samples (Miss) is lower (133.45, 136.20), than other methods by quite a large margin. Figure 3 further shows that the number of misclassified samples of GraphDA keeps the lowest during the entire online learning course for all four graphs (Arias-Castro et al., 2011).

The sparsity $s$ is an important parameter for GraphDA. We explore how $s$ affects the test error rate888The test error rate is calculated as $1-Acc_{{\bm{w}}_{t}}$ similar as done in (Duchi et al., 2011).. We compare the error rate of GraphDA with that of the non-sparse method Adam. Figure 5 clearly demonstrates that GraphDA has the least test error rate corresponding to the true model sparsity (26, 46, 92, 132 for these four subgraphs). When $s$ reaches the true sparsity ($26,46,92,132$ respectively), the testing error rate of GraphDA is the minimum.

**Classification Performance on different $t$ and $\mu$. ** We explore how different numbers of training sample and different $\mu$ affects the performance of each method similarly done in (Yang et al., 2010). First, we choose $t$ from set $\{100,200,300,\ldots,1000\}$, and tune the model based on classification accuracy. Results in Figure 6 show that when the number of training sample increases, the classification accuracy of all methods are increasing accordingly, but GraphDA enjoys the highest classification accuracy on both ${\bm{w}}_{t}$ and $\bar{\bm{w}}_{t}$. Second, we choose $\mu$ from the set $\{0.1,0.2,\ldots,1.0\}$ and fix $t=400$. As is reported in Figure 7, when $\mu$ is small (a harder classification task), all methods achieve lower accuracy; when $\mu$ is large (an easier task), all methods can obtain very high accuracy except StoIHT and GraphIHT. Again, Acc of GraphDA is the highest.

**Model interpretability. ** Table 1 shows that three out of the four indicators of feature-level performance for GraphDA score higher than for the other baseline methods. To be more specific, our method has the highest F1, 0.880, exceeding other methods by a large margin, which means that graph-structured information does help improve its performance. It also testifies that the head/tail projection during the online learning process does help capture more meaningful features than others. The nonzero ratio (NR) of ${\bm{w}}_{t}$ and ${\bar{{\bm{w}}}_{t}}$ is the least. The learned ${\bm{w}}_{t}$ in Figure 4 shows GraphDA successfully captures these subgraphs in ${\bm{w}}_{t}$, which are the closest to true models in terms of shapes and values (white colored pixels). Adam learns a totally dense model, and hence has worse performance. DA-IHT and $\ell_{1}$-RDA obtain very similar performance, probably because both of them use the dual averaging techniques. The results of StoIHT and GraphStoIHT testify that the online PGD-based methods hardly learn an interpretable model. In brief, our algorithm exploits the least number of features to learn the best model among all of the methods.

5.4. Results from MNIST dataset

The goal of online graph sparse linear regression is to minimize the regret where each $f_{t}({\bm{w}}_{t})$ is the least square loss defined as

[TABLE]

On this dataset, we use the least square loss as the objective function. The experiment is to compare the feature-level F1 score of different algorithms. We generate 1,400 data samples by using the following linear relation:

[TABLE]

where ${\bm{x}}_{t}\in\mathcal{N}({\bm{0}},{\bm{I}})$. We use three different strategies to obtain ${\bm{w}}^{*}$. The first one is to directly use the sparse images and then normalize them to the range $[0.0,1.0]$, which we call Normalized Model. The second is to generate ${\bm{w}}^{*}$ by letting all non-zeros be $1.0$, which is called Constant Model. The third is to generate the nonzero nodes by using Gaussian distribution $({\bm{w}}^{*})_{i}\sim\mathcal{N}(0,1)$ independently, which is Gaussian Model. Again, our dataset is partitioned into three parts: training, validating and testing samples. We increase the number of training samples $n$ from $\{50,100,\ldots,1000\}$ and then use 200 samples as validating dataset to tune the model. For all the eight online learning algorithms, we pass each training sample once and stop training when all training samples are used. The results shown in Figure 8 are generated from the 200 testing samples.

Adam, StoIHT and GraphStoIHT are excluded from comparison because of their inferior performance. From Figure 8, we can observe that when the training samples increase, the F1 score of all methods is increasing correspondingly. But the F1 score values of GraphDA in Normalized Model, Constant Model, and Gaussian Model are the highest among the six methods.

5.5. Results from KEGG dataset

To demonstrate that GraphDA can capture more meaningful features during online learning process, we test it on a real-world protein-protein iteration (PPI) network in (LeNail et al., 2017)999It was originally provided in KEGG (Kanehisa et al., 2016). This online learning scenario could be realistic since the training samples can be collected on the fly. More Details of the dataset including the data preprocessing are in Appendix B.3. We explore a specific gene pathway, HSA05213, related with endometrial cancer101010Details of pathway HSA05213(50 genes) can be found in https://www.genome.jp/dbget-bin/www_bget?hsa05213. Due to the lack of true labels and ground truth features (genes), we directly use the two data generation strategies in (LeNail et al., 2017), namely Strategy 1 (corresponding to a hard case) and Strategy 2 (corresponding to an easy case). After the data generation, we have 50 ground truth features. The number of positive and negative samples are both 100. The goal is to find out how many endometrial cancer-related genes are learned from different algorithms as done by (LeNail et al., 2017). All algorithms stop when they have seen 200 training samples.

We report the feature F1 score in Figure 9. GraphDA outperforms the other baseline methods in terms of both two strategies, with F1 score about 0.5 for Strategy 1 and about 0.9 for Strategy 2, higher than the rest methods. Interestingly, DA-OL and DA-SOL achieve better results only between 60 and 70 training samples and then become worse between 70 and 100. A possible explanation is that the learned model selected by the tuned parameters is not steady when the number of training samples seen is small. In addition to a better F1 score, another strength of GraphDA and DA-IHT is that the standard deviation of F1 score is smaller than other convex-based methods, including $\ell_{1}$-RDA, DA-GL, and DA-SGL.

In Figure 10, we show the identified genes by different methods. Clearly, GraphDA can find more meaningful genes, indicated by less blue nodes (genes found but not in HSA05213) and more red nodes (genes found and in HSA05213). However, all the other five baseline methods have many isolated nodes (not connected to cancer-related genes).

6. Conclusion and Future Work

In this paper, we propose a dual averaging-based method, GraphDA, for online graph-structured sparsity constraint problems. We prove that the minimization problem in the dual averaging step can be formulated as two equivalent optimization problems. By projecting the dual vector and primal variables onto lower dimensional subspaces, GraphDA can capture graph-structure information more effectively. Experimental evaluation on one benchmark dataset and two real-world graph datasets shows that GraphDA achieves better classification performance and stronger interpretability compared with the baseline methods so as to answer the two questions raised at the beginning of the experiment section. It remains interesting if one can prove that both the exact and approximated projections have non-regret bound under some proper assumption, and if one can explore learning a model under the setting that true features are time evolving (Hou et al., 2017).

7. Acknowledgements

The work of Yiming Ying is supported by the National Science Foundation (NSF) under Grant No #1816227. The work of Baojian Zhou and Feng Chen is supported by the NSF under Grant No #1815696 and #1750911.

Appendix A Reproducibility

A.1. Implementation details

All experiments are tested on a server of Intel Xeon(R) 2.40GHZ E5-2680 with 251GB of RAM. The code is written in Python2.7 and C language with the standard C11. The implementation of the head and tail projection follows the original implementation111111The two projections were originally implemented in C++, which are available at: https://github.com/ludwigschmidt/cluster_approx. We implement them by using C language. Taking the advantage of the continuous memory of arrays in C, our code is faster than original one.. We present the pseudo code in Algorithm 2 below. The two projections are essentially two binary search algorithms. Each iteration of the binary search executes the Prize Collecting Steiner Tree (PCST) algorithm (Johnson et al., 2000) on the target graph. Both projections have two main parameters: a lower bound sparsity $s_{l}$ and an upper bound sparsity $s_{h}$. In all of the experiments, two sparsity parameters have been set to $s_{l}=p/2$ and $s_{h}=s_{l}*(1+\omega)$ for the head projection, where $\omega$ is the tolerance parameter set to $0.1$. For the tail projection, we set $s_{l}=s$ and $s_{h}=s_{l}*(1+\omega)$. The binary search algorithm terminates when it reaches $max\_iter=20$ maximum iterations. Line 7 of Algorithm 2 is the PCST algorithm proposed in (Hegde et al., 2014a). We use a non-root version and Goemans-Williamson pruning method to prune the final forest.

A.2. Parameter tuning

Initial parameters ${\bf w}_{0}$ of all baseline methods and proposed algorithms are zero vectors ${\bf w}_{0}={\bf 0}$, which means we train all methods starting from a zero point. We list all related methods and their corresponding parameter settings below. (1) $\ell_{1}$-RDA is the enhanced version provided in Algorithm 2 of (Xiao, 2010). There are three parameters: The $\ell_{1}$-regularization parameter $\lambda$ is chosen from $\{$0.0001, 0.0005, 0.001, 0.005, 0.01, 0.03, 0.05, 0.1, 0.3, 0.5, 1.0, 3.0, 5.0, 10.0$\}$ which is a superset used in (Xiao, 2010). The parameter $\gamma$ to control the learning rate is chosen from $\{$1.0, 5.0, 10.0, 50.0, 100.0, 500.0, 1000.0, 5000.0, 10000.0 $\}$, and the sparsity-enhancing parameter $\rho$ is chosen from $\{$0.0, 0.00001, 0.000005, 0.0001, 0.0005, 0.001, 0.005, 0.01, 0.05, 0.1, 0.5, 1.0$\}$, where 0.0 is for the basic regularization. All the three parameter sets are supersets used in (Xiao, 2010). (2) Adam. We directly use the parameters $\beta_{1}=0.9,\beta_{2}=0.999,\epsilon=10^{-8}$ provided in (Kingma and Ba, 2014). For the magnitude of steps in parameter space $\alpha$, we choose it from $\{$0.0001, 0.0005, 0.001, 0.005, 0.01, 0.05, 0.1, 0.5 $\}$. (3) DA-GL/SGL have two main parameters, $\lambda$ to control the sparsity and $\gamma$ to control the learning rate. We choose $3\times 3$ grids as groups for Benchmark dataset and choose $2\times 2$ grids for MNIST dataset. (4) DA-SGL has an additional parameter $\gamma_{g}$, which is set to 1.0 for all groups as done in (Yang et al., 2010). For each group $i$, there exists an additional parameter $r_{i}$ for DA-SGL. We set it as default value $r_{i}=1$ as recommended by the authors. (5) AdaGrad has two main parameters, $\lambda$ to control sparsity and $\eta$ to control the learning rate, which is from $\{$0.0001 , 0.0005, 0.001, 0.005, 0.01, 0.05, 0.1, 0.5, 1.0, 5.0, 10.0, 50.0, 100.0, 500.0, 1000.0, 5000.0$\}$. (6) StoIHT has two parameters: sparsity $s$ from $\{$5, 10, 15, 20, 25, 26, 30, 35, 40, 45, 46, 50, 55, 60, 65, 70, 75, 80, 85, 90, 92, 95, 100, 105, 110, 115, 120, 125, 130, 132, 135, 140, 145, 150$\}$, and $\gamma$ to control the learning rate. (7) GraphStoIHT shares the same parameter settings (sparsity $s$ and $\gamma$) as GraphDA. The block size of GraphStoIHT and StoIHT are set to 1. (8) GraphDA has parameters $\gamma$ and $s$.

Appendix B More experimental results

B.1. More results from Benchmark dataset

We present the results of Graph02, Graph03 and Graph04 in Table 2, 3, 4, respectively. Basically, we show the classification performance (Acc, Miss, AUC) and feature-level performance (Pre, Rec, F1, NR). The size of validating and testing dataset are both 400. All results are averaged from 20 trials of experiment.

B.2. More results from MNIST Dataset

We show the results on image id $\{1,2,3,6,7,8,9\}$. To make the task more challenging, these 7 images are the sparsest images (with the digits forming a connected component) selected from MNIST dataset. The sparsity parameter $s$ of DA-IHT and GraphDA is chosen from $\{$30, 32, …, 100$\}$ with step size 2. Figure 11 reports the results.

B.3. More results from KEGG Pathways

This PPI network contains a total of 229 pathways. Each pathway often involves a specific biological function, e.g. metabolism. We restrict our analysis on 225 pathways (by removing 4 empty pathways), which contains 5,374 genes with 78,545 edges. These genes form a connected graph. There exists an edge if two proteins (genes) physically interact with each other (Li et al., 2017). Weights of edges stand for the confidence of these interactions. There are 7,368 genes with null values. We sample these null values from $\mathcal{N}(0,1)$.

Due to the inferior performance of Adam, StoIHT and GraphStoIHT, we exclude them from experimental evaluation. Notice that DA-GL and DA-SGL need groups as priors. However, the groups (pathways) of this PPI network have overlapping features. To remedy this issue, we simply replicate these overlapping features as suggested in (Yang et al., 2010; Jacob et al., 2009), and by doing so, the two baselines are still applicable. For these two non-convex methods, we choose the sparsity parameter from $\{$40, 45, 50, 55, 60 $\}$. We report the averaged results from 20 trials in Figure 12.