Year: 2021

Abstract:
Graph Representation Learning (GRL) methods have impacted fields from chemistry to social science. However, their algorithmic implementations are specialized to specific use-cases e.g. message passing methods are run differently from node embedding ones. Despite their apparent differences, all these methods utilize the graph structure, and therefore, their learning can be approximated with stochastic graph traversals. We propose Graph Traversal via Tensor Functionals (GTTF), a unifying meta-algorithm framework for easing the implementation of diverse graph algorithms and enabling transparent and efficient scaling to large graphs. GTTF is founded upon a data structure (stored as a sparse tensor) and a stochastic graph traversal algorithm (described using tensor operations). The algorithm is a functional that accept two functions, and can be specialized to obtain a variety of GRL models and objectives, simply by changing those two functions. We show for a wide class of methods, our algorithm learns in an unbiased fashion and, in expectation, approximates the learning as if the specialized implementations were run directly. With these capabilities, we scale otherwise non-scalable methods to set state-of-the-art on large graph datasets while being more efficient than existing GRL libraries – with only a handful of lines of code for each method specialization. GTTF and its various GRL implementations are on: https://github.com/isi-usc-edu/gttf

heading: 1 INTRODUCTION
text: Graph representation learning (GRL) has become an invaluable approach for a variety of tasks, such as node classification (e.g., in biological and citation networks; Veličković et al. (2018); Kipf & Welling (2017); Hamilton et al. (2017); Xu et al. (2018)), edge classification (e.g., link prediction for social and protein networks; Perozzi et al. (2014); Grover & Leskovec (2016)), entire graph classification (e.g., for chemistry and drug discovery Gilmer et al. (2017); Chen et al. (2018a)), etc.
In this work, we propose an algorithmic unification of various GRL methods that allows us to re-implement existing GRL methods and introduce new ones, in merely a handful of code lines per method. Our algorithm (abbreviated GTTF, Section 3.2), receives graphs as input, traverses them using efficient tensor1 operations, and invokes specializable functions during the traversal. We show function specializations for recovering popular GRL methods (Section 3.3). Moreover, since GTTF is stochastic, these specializations automatically scale to arbitrarily large graphs, without careful derivation per method. Importantly, such specializations, in expectation, recover unbiased gradient estimates of the objective w.r.t. model parameters.
1To disambiguate: by tensors, we refer to multi-dimensional arrays, as used in Deep Learning literature; and by operations, we refer to routines such as matrix multiplication, advanced indexing, etc
GTTF uses a data structure Â (Compact Adjacency, Section 3.1): a sparse encoding of the adjacency matrix. Node v contains its neighbors in row Â[v] , Âv, notably, in the first degree(v) columns of Â[v]. This encoding allows stochastic graph traversals using standard tensor operations. GTTF is a functional, as it accepts functions ACCUMULATEFN and BIASFN, respectively, to be provided by each GRL specialization to accumulate necessary information for computing the objective, and optionally to parametrize sampling procedure p(v’s neighbors | v). The traversal internally constructs a walk forest as part of the computation graph. Figure 1 depicts the data structure and the computation. From a generalization perspective, GTTF shares similarities with Dropout (Srivastava et al., 2014).
Our contributions are: (i) A stochastic graph traversal algorithm (GTTF) based on tensor operations that inherits the benefits of vectorized computation and libraries such as PyTorch and Tensorflow. (ii) We list specialization functions, allowing GTTF to approximately recover the learning of a broad class of popular GRL methods. (iii) We prove that this learning is unbiased, with controllable variance. Wor this class of methods, (iv) we show that GTTF can scale previously-unscalable GRL algorithms, setting the state-of-the-art on a range of datasets. Finally, (v) we open-source GTTF along with new stochastic traversal versions of several algorithms, to aid practitioners from various fields in applying and designing state-of-the-art GRL methods for large graphs.

heading: 2 RELATED WORK
text: We take a broad standpoint in summarizing related work to motivate our contribution.
Method Fa m
ily
Sc al
e
L ea
rn in
g
Models GCN, GAT MP 7 exact node2vec NE 3 approx
WYS NE 7 exact Stochastic Sampling Methods
SAGE MP 3 approx FastGCN MP 3 approx LADIES MP 3 approx
GraphSAINT MP 3 approx CluterGCN MP 3 heuristic
Software Frameworks PyG Both inherits / reDGL Both implements
Algorithmic Abstraction (ours) GTTF Both 3 approx
Models for GRL have been proposed, including message passing (MP) algorithms, such as Graph Convolutional Network (GCN) (Kipf & Welling, 2017), Graph Attention (GAT) (Veličković et al., 2018); as well as node embedding (NE) algorithms, including node2vec (Grover & Leskovec, 2016), WYS (Abu-El-Haija et al., 2018); among many others (Xu et al., 2018; Wu et al., 2019; Perozzi et al., 2014). The full-batch GCN of Kipf & Welling (2017), which drew recent attention and has motivated many MP algorithms, was not initially scalable to large graphs, as it processes all graph nodes at every training step. To scale MP methods to large graphs, researchers proposed Stochastic Sampling Methods that, at each training step, assemble a batch constituting subgraph(s) of the (large) input graph. Some of these sampling methods yield unbiased gradient estimates (with some variance) including SAGE (Hamilton et al., 2017), FastGCN (Chen et al., 2018b), LADIES (Zou et al., 2019), and GraphSAINT (Zeng et al., 2020). On the other hand, ClusterGCN (Chiang et al., 2019) is a heuristic in the sense that, despite its good performance, it provides no guarantee of unbiased gradient estimates of the full-batch learning. Gilmer et al. (2017) and Chami et al. (2021) generalized many GRL models into Message Passing and Auto-Encoder frameworks.
These frameworks prompt bundling of GRL methods under Software Libraries, like PyG (Fey & Lenssen, 2019) and DGL (Wang et al., 2019), offering consistent interfaces on data formats.
We now position our contribution relative to the above. Unlike generalized message passing (Gilmer et al., 2017), rather than abstracting the model computation, we abstract the learning algorithm. As a result, GTTF can be specialized to recover the learning of MP as well as NE methods. Morever, unlike Software Frameworks, which are re-implementations of many algorithms and therefore inherit the scale and learning of the copied algorithms, we re-write the algorithms themselves, giving them new properties (memory and computation complexity), while maintaining (in expectation) the original algorithm outcomes. Further, while the listed Stochastic Sampling Methods target MP algorithms (such as GCN, GAT, alike), as their initial construction could not scale to large graphs, our learning algorithm applies to a wider class of GRL methods, additionally encapsulating NE methods. Finally, while some NE methods such as node2vec (Grover & Leskovec, 2016) and DeepWalk (Perozzi et al., 2014) are scalable in their original form, their scalability stems from their multi-step process: sample many (short) random walks, save them to desk, and then learn node embeddings using positional embedding methods (e.g., word2vec, Mikolov et al. (2013)) – they are sub-optimal in the sense that their first step (walk sampling) takes considerable time (before training even starts) and also places an artificial limit on the number of training samples (number of simulated walks), whereas our algorithm conducts walks on-the-fly whilst training.

heading: 3 GRAPH TRAVERSAL VIA TENSOR FUNCTIONALS (GTTF)
text: At its core, GTTF is a stochastic algorithm that recursively conducts graph traversals to build representations of the graph. We describe the data structure and traversal algorithm below, using the following notation. G = (V,E) is an unweighted graph with n = |V | nodes and m = |E| edges, described as a sparse adjacency matrix A ∈ {0, 1}n×n. Without loss of generality, let the nodes be zero-based numbered i.e. V = {0, . . . , n− 1}. We denote the out-degree vector δ ∈ Zn – it can be calculated by summing over rows of A as δu = ∑ v∈V A[u, v]. We assume δu > 0 for all u ∈ V : pre-processing can add self-connections to orphan nodes. B denotes a batch of nodes.

heading: 3.1 DATA STRUCTURE
text: Internally, GTTF relies on a reformulation of the adjacency matrix, which we term CompactAdj (for "Compact Adjacency", Figure 1c). It consists of two tensors:
1. δ ∈ Zn, a dense out-degree vector (figure 1c, right) 2. Â ∈ Zn×n, a sparse edge-list matrix in which the row u contains left-aligned δu non-zero
values. The consecutive entries {Â[u, 0], Â[u, 1], . . . , Â[u, δu − 1]} contain IDs of nodes receiving an edge from node u. The remaining |V | − δu are left unset, therefore, Â only occupies O(m) memory when stored as a sparse matrix (Figure 1c, left).
CompactAdj allows us to concisely describe stochastic traversals using standard tensor operations. To uniformly sample a neighbor to node u ∈ V , one can draw r ∼ U [0..(δu − 1)], then get the neighbor ID with Â[u, r]. In vectorized form, given node batch B and access to continuous U [0, 1), we sample neighbors for each node in B as: R ∼ U [0, 1)b, where b = |B|, then B′ = Â[B, bR ◦ δ[B]c] is a b-sized vector, with B′u containing a neighbor of Bu, floor operation b.c is applied element-wise, and ◦ is Hadamard product.

heading: 3.2 STOCHASTIC TRAVERSAL FUNCTIONAL ALGORITHM
text: Our traversal algorithm starts from a batch of nodes. It expands from each into a tree, resulting in a walk forest rooted at the nodes in the batch, as depicted in Figure 1d. In particular, given a node batch B, the algorithm instantiates |B| seed walkers, placing one at every node in B. Iteratively, each walker first replicates itself a fanout (f ) number of times. Each replica then samples and transitions to a neighbor. This process repeats a depth (h) number of times. Therefore, each seed walker becomes the ancestor of a f -ary tree with height h. Setting f = 1 recovers traditional random walk. In practice, we provide flexibility by allowing a custom fanout value per depth.
Algorithm 1: Stochastic Traverse Functional, parametrized by ACCUMULATEFN and BIASFN. input: u (current node); T ← [] (path leading to u, starts empty); F (list of fanouts);
ACCUMULATEFN (function: with side-effects and no return. It is model-specific and records information for computing model and/or objective, see text); BIASFN ← U (function mapping u to distribution on u’s neighbors, defaults to uniform)
1 def Traverse(T , u, F , ACCUMULATEFN, BIASFN): 2 if F.size() = 0 then return # Base case. Traversed up-to requested depth 3 f ← F.pop() # fanout duplication factor (i.e. breadth) at this depth. 4 sample_bias ← BIASFN(T, u) 5 if sample_bias.sum() = 0 then return # Special case. No sampling from zero mass 6 sample_bias ← sample_bias / sample_bias.sum() # valid distribution 7 K ← Sample(Â [u, :δu]] ,sample_bias, f) # Sample f nodes from u’s neighbors 8 for k ← 0 to f − 1 do 9 Tnext ← concatenate(T, [u])
10 ACCUMULATEFN(Tnext,K[k], f) 11 Traverse(Tnext,K[k], f,ACCUMULATEFN,BIASFN) # Recursion 12 13 def Sample(N , W , f ): 14 C ← tf.cumsum(W ) # Cumulative sum. Last entry must = 1. 15 coin_flips← tf.random.uniform((f, ), 0, 1) 16 indices← tf.searchsorted(C, coin_flips) 17 return N [indices]
Functional Traverse is listed in Algorithm 1. It accepts: a batch of nodes2; a list of fanout values F (e.g. to F = [3, 5] samples 3 neighbors per u ∈ B, then 5 neighbors for each of those); and more notably, two functions: ACCUMULATEFN and BIASFN. These functions will be called by the functional on every node visited along the traversal, and will be passed relevant information (e.g. the path taken from root seed node). Custom settings of these functions allow recovering wide classes of graph learning methods. At a high-level, our functional can be used in the following manner:
1. Construct model & initialize parameters (e.g. to random). Define ACCUMULATEFN and BIASFN. 2. Repeat (many rounds):
i. Reset accumulation information (from previous round) and then sample batch B ⊂ V . ii. Invoke Traverse on (B, ACCUMULATEFN, BIASFN), which invokes the FN’s, allowing the first
to accumulate information sufficient for running the model and estimating an objective. iii. Use accumulated information to: run model, estimate objective, apply learning rule (e.g. SGD).
ACCUMULATEFN is a function that is used to track necessary information for computing the model and the objective function. For instance, an implementation of DeepWalk (Perozzi et al., 2014) on top of GTTF, specializes ACCUMULATEFN to measure an estimate of the sampled softmax likelihood of nodes’ positional distribution, modeled as a dot-prodct of node embeddings. On the other hand, GCN (Kipf & Welling, 2017) on top of GTTF uses it to accumulate a sampled adjacency matrix, which it passes to the underlying model (e.g. 2-layer GCN) as if this were the full adjacency matrix.
BIASFN is a function that customizes the sampling procedure for the stochastic transitions. If provided, it must yield a probability distribution over nodes, given the current node and the path that lead to it. If not provided, it defaults to U , transitioning to any neighbor with equal probability. It can be defined to read edge weights, if they denote importance, or more intricately, used to parameterize a second order Markov Chain (Grover & Leskovec, 2016), or use neighborhood attention to guide sampling (Veličković et al., 2018), as discussed in the Appendix.
2Our pseudo-code displays the traversal starting from one node rather than a batch only for clarity, as our actual implementation is vectorized e.g. u would be a vector of nodes, T would be a 2D matrix with each row containing transition path preceeding the corresponding entry in u, ... etc. Refer to Appendix and code.
3.3 SOME SPECIALIZATIONS OF ACCUMULATEFN & BIASFN

heading: 3.3.1 MESSAGE PASSING: GRAPH CONVOLUTIONAL VARIANTS
text: These methods, including (Kipf & Welling, 2017; Hamilton et al., 2017; Wu et al., 2019; Abu-ElHaija et al., 2019; Xu et al., 2018) can be approximated by by initializing Ã to an empty sparse n× n matrix, then invoking Traverse (Algorithm 1) with u = B; F to list of fanouts with size h; Thus ACCUMULATEFN and BIASFN become:
def ROOTEDADJACC(T, u, f): Ã[u, T−1]← 1; (1)
def NOREVISITBIAS(T, u): return 1[Ã[u].sum() = 0] ~1δu δu ; (2)
where ~1n is an n-dimensional all-ones vector, and negative indexing T−k is the kth last entry of T . If a node has been visited through the stochastic traversal, then it already has fanout number of neighbors and NOREVISITBIAS ensures it does not get revisited for efficiency, per line 5 of Algorithm 1. Afterwards, the accumulated stochastic Ã will be fed3 into the underlying model e.g. for a 2-layer GCN of Kipf & Welling (2017):
GCN(Ã,X;W1,W2) = softmax( ◦ A ReLu( ◦ AXW1)W2); (3)
with ◦ A = D′1/2D̃′−1Ã′D′−1/2; D′ = diag(δ′); δ′ = ~1>nA ′; renorm trick︷ ︸︸ ︷ Ã′ = In×n + Ã
Lastly, h should be set to the receptive field required by the model for obtaining output dL-dimensional features at the labeled node batch. In particular, to the number of GC layers multiplied by the number of hops each layers access. E.g. hops=1 for GCN but customizable for MixHop and SimpleGCN.

heading: 3.3.2 NODE EMBEDDINGS
text: Given a batch of nodes B ⊆ V , DeepWalk4 can be implemented in GTTF by first initializing loss L to the contrastive term estimating the partition function of log-softmax:
L ← ∑ u∈B log E v∼Pn(V ) [exp(〈Zu, Zv〉)] , (4)
where 〈., .〉 is dot-product notation, Z ∈ Rn×d is the trainable embedding matrix with Zi ∈ Rd is d-dimensional embedding for node u ∈ V . In our experiments, we estimate the expectation by taking 5 samples and we set the negative distribution Pn(V = v) ∝ δ 3 4 v , following Mikolov et al. (2013).
The functional is invoked with no BIASFN and ACCUMULATEFN =
def DEEPWALKACC(T, u, f): L ← L− 〈 Zi,
CT∑ k=1 η [T−k] ( C−k+1 C ) ZT−k
〉 ; η[u] ← η [T−1]
f ; (5)
where hyperparameter C indicates maximum window size (inherited from word2vec, Mikolov et al., 2013), in the summation on k does not access invalid entries of T as CT , min(C, T .size), the scalar fraction ( C−k+1 C ) is inherited from context sampling of word2vec (Section 3.1 in Levy et al., 2015), and rederived for graph context by Abu-El-Haija et al. (2018), and η[u] stores a scalar per node on the traversal Walk Forest, which defaults to 1 for non-initialized entries, and is used as a correction term. DeepWalk conducts random walks (visualized as a straight line) whereas our walk tree has a branching factor of f . Setting fanout f = 1 recovers DeepWalk’s simulation, though we found f > 1 outperforms within fewer iterations e.g. f = 5, within 1 epoch, outperforms DeepWalk’s published implementation. Learning can be performed using the accumulated L as: Z ← Z − ∇ZL;

heading: 4 THEORETICAL ANALYSIS
text: Due to space limitations, we include the full proofs of all propositions in the appendix. 3Before feeding the batch to model, in practice, we find nodes not reached by traversal and remove their corresponding rows (and also columns) from X (and A). 4We present more methods in the Appendix.
4.1 ESTIMATING kTH POWER OF TRANSITION MATRIX
We show that it is possible with GTTF to accumulate an estimate of transition T matrix to power k. Let Ω denote the walk forest generated by GTTF, Ω(u, k, i) as the ith node in the vector of nodes at depth k of the walk tree rooted at u ∈ B, and tu,v,ki as the indicator random variable 1[Ω(u, k, i) = v]. Let the estimate of the kth power of the transition matrix be denoted T̂ k. Entry T̂ ku,v should be an unbiased estimate of T ku,v for u ∈ B, with controllable variance. We define:
T̂ ku,v = ∑fk i=1 t u,v,k i
fk (6)
The fraction in Equation 6 counts the number of times the walker starting at u visits v in Ω, divided by the total number of nodes visited at the kth step from u.
Proposition 1. (UNBIASEDTK) T̂ ku,v as defined in Equation 6, is an unbiased estimator of T ku,v
Proposition 2. (VARIANCETK) Variance of our estimate is upper-bounded: Var[T̂ ku,v] ≤ 1
4fk
Naive computation of kth powers of the transition matrix can be efficiently computed via repeated sparse matrix-vector multiplication. Specifically, each column of T k can be computed in O(mk), where m is the number of edges in the graph. Thus, computing T k in its entirety can be accomplished in O(nmk). However, this can still become prohibitively expensive if the graph grows beyond a certain size. GTTF on the other hand can estimate T k in time complexity independent of the size of the graph, (Prop. 8), with low variance. Transition matrix powers are useful for many GRL methods. (Qiu et al., 2018)

heading: 4.2 UNBIASED LEARNING
text: As a consequence of Propositions 1 and 2, GTTF enables unbiased learning with variance control for classes of node embedding methods, and provides a convergence guarantee for graph convolution models under certain simplifying assumptions.
We start by analyzing node embedding methods. Specifically, we cover two general types. The first is based on matrix factorization of the power-series of transition matrix. and the second is based on cross-entropy objectives, e.g., like DeepWalk (Perozzi et al., 2014), node2vec (Grover & Leskovec, 2016), These two are shown in Proposations 3 and 4
Proposition 3. (UNBIASEDTFACTORIZATION) Suppose L = 12 ||LR− ∑ k ckT k||2F , i.e. factorization objective that can be optimized by gradient descent by calculating∇L,RL, where ck’s are scalar coefficients. Let its estimate L̂ = 12 ||LR− ∑ k ckT̂ k||2F , where T̂ is obtained by GTTF according to Equation 6. Then E[∇L,RL̂] = ∇L,RL.
Proposition 4. (UNBIASEDLEARNNE) Learning node embeddings Z ∈ Rn×d with objective function L, decomposable as L(Z) = ∑ u∈V L1(Z, u) − ∑ u,v∈V ∑ k L2(T k, u, v)L3(Z, u, v), where L2 is linear over T k, then using T̂ k yields an unbiased estimate of∇ZL.
Generally, L1 (and L3) score the similarity between disconnected (and connected) nodes u and v. The above form of L covers a family of contrastive learning objectives that use cross-entropy loss and assume a logistic or (sampled-)softmax distributions. We provide, in the Appendix, the decompositions for the objectives of DeepWalk (Perozzi et al., 2014), node2vec (Grover & Leskovec, 2016) and WYS (Abu-El-Haija et al., 2018).
Proposition 5. (UNBIASEDMP) Given input activations, H(l−1), graph conv layer (l) can use rooted adjacency Ã accumulated by ROOTEDADJACC (1), to provide unbiased pre-activation output,
i.e. E [ ◦ AkH(l−1)W (l) ] = ( D ′−1/2A′D ′−1/2 )k H(l−1)W (l), with A′ and D′ defined in (3).
Proposition 6. (UNBIASEDLEARNMP) If objective to a graph convolution model is convex and Lipschitz continous, with minimizer θ∗, then utilizing GTTF for graph convolution converges to θ∗.

heading: 4.3 COMPLEXITY ANALYSIS
text: Proposition 7. STORAGE complexity of GTTF is O(m+ n). Proposition 8. TIME complexity of GTTF is O(bfh) for batch size b, fanout f , and depth h.
Proposition 8 implies the speed of computation is irrespective of graph size. Methods implemented in GTTF inherit this advantage. For instance, the node embedding algorithm WYS (Abu-El-Haija et al., 2018) is O(n3), however, we apply its GTTF implementation on large graphs.

heading: 5 EXPERIMENTS
text: We conduct experiments on 10 different graph datasets, listed in in Table 1. We experimentally demonstrate the following. (1) Re-implementing baseline method using GTTF maintains performance. (2) Previously-unscalable methods, can be made scalable when implemented in GTTF. (3) GTTF achieves good empirical performance when compared to other sampling-based approaches handdesigned for Message Passing. (4) GTTF consumes less memory and trains faster than other popular Software Frameworks for GRL. To replicate our experimental results, for each cell of the table in our code repository, we provide one shell script to produce the metric, except when we indicate that the metric is copied from another paper. Unless otherwise stated, we used fanout factor of 3 for GTTF implementations. Learning rates and model hyperparameters are included in the Appendix.

heading: 5.1 NODE EMBEDDINGS FOR LINK PREDICTION
text: In link prediction tasks, a graph is partially obstructed by hiding a portion of its edges. The task is to recover the hidden edges. We follow a popular approach to tackle this task: first learn node embedding Z ∈ Rn×d from the observed graph, then predict the link between nodes u and v with score∝ Z>u Zv . We use two ranking metrics for evaluations: ROC-AUC, which is a ranking objective: how well do methods rank the hidden edges above randomly-sampled negative edges and Mean Rank.
We re-implement Node Embedding methods, DeepWalk (Perozzi et al., 2014) and WYS (Abu-ElHaija et al., 2018), into GTTF (abbreviated F ). Table 2 summarizes link prediction test performance. LiveJournal and Reddit are large datasets, where original implementation of WYS is unable to scale to. However, scalable F(WYS) sets new state-of-the-art on these datasets. For PPI and HepTh datasets, we copy accuracy numbers for DeepWalk and WYS from (Abu-El-Haija et al., 2018). For LiveJournal, we copy accuracy numbers for DeepWalk and PBG from (Lerer et al., 2019) – note that a well-engineered approach (PBG, (Lerer et al., 2019)), using a mapreduce-like framework, is under-performing compared to F(WYS), which is a few lines specialization of GTTF.

heading: 5.2 MESSAGE PASSING FOR NODE CLASSIFICATION
text: We implement in GTTF the message passing models: GCN (Kipf & Welling, 2017), GraphSAGE (Hamilton et al., 2017), MixHop (Abu-El-Haija et al., 2019), SimpleGCN (Wu et al., 2019), as their computation is straight-forward. For GAT (Veličković et al., 2018) and GCNII (Chen et al., 2020), as they are more intricate, we download the authors’ codes, and wrap them as-is with our functional.
We show that we are able to run these models in Table 3 (left and middle), and that GTTF implementations matches the baselines performance. For the left table, we copy numbers from the published papers. However, we update GAT to work with TensorFlow 2.0 and we use our updated code (GAT*).

heading: 5.3 EXPERIMENTS COMPARING AGAINST SAMPLING METHODS FOR MESSAGE PASSING
text: We now compare models trained with GTTF (where samples are walk forests) against sampling methods that are especially designed for Message Passing algorithms (GraphSAINT and ClusterGCN), especially since their sampling strategies do not match ours.
Table 3 (right) shows test performance on node classification accuracy on a large dataset: Products. We calculate the accuracy forF (SAGE), but copy from (Hu et al., 2020) the accuracy for the baselines: GraphSAINT (Zeng et al., 2020) and ClusterGCN (Chiang et al., 2019) (both are message passing methods); and also node2vec (Grover & Leskovec, 2016) (node embedding method).

heading: Dataset Split # Nodes # Edges # Classes Nodes Edges Tasks

heading: 5.4 RUNTIME AND MEMORY COMPARISON AGAINST OPTIMIZED SOFTWARE FRAMEWORKS
text: In addition to the accuracy metrics discussed above, we also care about computational performance. We compare against software frameworks DGL (Wang et al., 2019) and PyG (Fey & Lenssen, 2019). These software frameworks offer implementations of many methods. Table 4 summarizes the following. First (left), we show time-per-epoch on large graphs of their implementation of
GraphSAGE, compared with GTTF’s, where we make all hyper parameters to be the same (of model architecture, and number of neighbors at message passing layers). Second (middle), we run their GCN implementation on small datasets (Cora, Citeseer, Pubmed) to show peak memory usage. The run times between GTTF, PyG and DGL are similar for these datasets. The comparison can be found in the Appendix. While the aforementioned two comparisons are on popular message passing methods, the third (right) chart shows a popular node embedding method: node2vec’s link prediction test ROC-AUC in relation to its training runtime.

heading: 6 CONCLUSION
text: We present a new algorithm, Graph Traversal via Tensor Functionals (GTTF) that can be specialized to re-implement the algorithms of various Graph Representation Learning methods. The specialization takes little effort per method, making it straight-forward to port existing methods or introduce new ones. Methods implemented in GTTF run efficiently as GTTF uses tensor operations to traverse graphs. In addition, the traversal is stochastic and therefore automatically makes the implementations scalable to large graphs. We theoretically show that the learning outcome due to the stochastic traversal is in expectation equivalent to the baseline when the graph is observed at-once, for popular GRL methods we analyze. Our thorough experimental evaluation confirms that methods implemented in GTTF maintain their empirical performance, and can be trained faster and using less memory even compared to software frameworks that have been thoroughly optimized.

heading: 7 ACKNOWLEDGEMENTS
text: We acknowledge support from the Defense Advanced Research Projects Agency (DARPA) under award FA8750-17-C-0106.

heading: APPENDIX

heading: A HYPERPARAMETERS
text: For the general link prediction tasks we used a |B| = |V |, C = 5, f = 3, 10 negative samples per edge, Adam optimizer with a learning rate of 0.5, multiplied by a factor of 0.2, every 50 steps, for 200 total iterations. The differences are listed below.
The Reddit dataset was trained using a starting learning rate of 2.0, decaying 50% every 10 iterations.
The LiveJournal task was trained using a fixed learning rate of 0.001, |B| = 5000, f = 2, and 50 negative samples per edge.
For the node classifications tasks:
For F(SimpleGCN) on Amazon, we use f = [15, 15], a batch size of 1024, and a learning rate of 0.02, decaying by a factor 0f 0.2 after 2 and 6 epochs for a total of 25 epochs. On Reddit, it is the same except f = [25, 25]
For F(SAGE) on Amazon we use f = [20, 10], a two layer model, a batch size of 256, and fixed learning rates of 0.001 and 0.002 respectively. On reddit we use f = [25, 20], a fixed learning rate of 0.001, hidden dimension of 256 and a batch size of 1024. On the Products dataset, we used f = [15, 10], a fixed learning learning rate of 0.001 and a batch size of 1024, a hidden dimension of 256 and a fixed learning rate of 0.003.
For GAT (baseline), we follow the authors code and hyperparameters: for Cora and Citeseer, we use Adam with learning rate of 0.005, L2 regularization of 0.0005, 8 attention heads on the first layer and 1 attention head on the output layer. For Pubmed, we use Adam with learning rate of 0.01, L2 regularization of 0.01, 8 attention heads on the first layer and 8 attention heads on the output layer. For F(GAT), we use the same aforementioned hyperparameters, a fanout of 3 and traversal depth of 2 (to cover two layers) i.e. F = [3, 3]. For F(GCN), we use the authors’ recommended hyperparameters. Learning rate of 0.005, 0.001 L2 regularization, and F = [3, 3], for all datasets. For both methods, we apply “patience” and stop the training if validation loss does not improve for 100 consecutive epochs, reporting the test accuracy at the best validation loss. For F(MixHop), we wrap the authors’ script and use their hyperparameters. For F(GCNII), we use F = [5, 5, 5, 5, 5, 5], as their models are deep (64 layers for Cora). Otherwise, we inherit their network hyperparameters (latent dimensions, number of layers, dropout factor, and their introduced coefficients), as they have tuned them per dataset, but we change the learning rate to 0.005 (half of what they use) and we extend the patience from 100 to 1000, and extend the maximum number of epochs from 1500 to 5000 – this is because we are presenting a subgraph at each epoch, and therefore we intuitively want to slow down the learning per epoch, which is similar to the practice when someone applies Dropout to a neural networks. We re-run their shell scripts, with their code modified to use the Rooted Adjacency rather than the real adjacency, which is sampled at every epoch.
MLP was trained with 1 layer and a learning rate of 0.01.

heading: B PROOFS

heading: B.1 PROOF OF PROPOSITION 1
text: Proof. E[T̂ ku,v] = E [ fk∑ i=1 tu,v,ki
fk
] = fk∑ i=1 E[tu,v,ki ]
fk =
fk∑ i=1 P [tu,v,ki = 1]
fk =
fk∑ i=1 T ku,v
fk = T ku,v

heading: B.2 PROOF OF PROPOSITION 2
text: Proof. Var[T̂ ku,v] = ∑fk i=1 V ar[t u,v,k i ]
f2k = fkT ku,v(1− T ku,v) f2k = T ku,v(1− T ku,v) fk
Since 0 ≤ T ku,v ≤ 1, then T ku,v(1−T ku,v) is maximized with T ku,v = 1
2 . Hence V ar[T̂ ku,v] ≤
1
4fk

heading: B.3 PROOF OF PROPOSITION 3
text: Proof. Consider a d-dimensional factorization of ∑ k ckT k, where ck’s are scalar coefficients:
L = 1 2 ∣∣∣∣∣ ∣∣∣∣∣LR−∑
k
ckT k ∣∣∣∣∣ ∣∣∣∣∣ 2
F
, (7)
parametrized by L,R> ∈ Rn×d. The gradients of L w.r.t. parameters are:
∇LL = ( LR> −
∑ k
ckT k ) R> and ∇RL = L> ( LR> −
∑ k
ckT k ) . (8)
Given estimate objective L (replacing T̂ with using GTTF-estimated T̂ ):
L̂ = 1 2 ∣∣∣∣∣ ∣∣∣∣∣LR−∑
k
ckT̂ k ∣∣∣∣∣ ∣∣∣∣∣ 2
F
. (9)
It follows that: E [ ∇LL̂ ] = E [( LR> −
∑ k
ckT̂ k ) R> ]
= E [( LR> −
∑ k
ckT̂ k )] R> Scaling property of expectation
= ( LR> −
∑ k
ckE [ T̂ k ]) R> Linearity of expectation
= ( LR> −
∑ k
ckT k ) R> Proof of Proposition 1
= ∇LL The above steps can similarly be used to show E [ ∇RL̂ ] = ∇RL

heading: B.4 PROOF OF PROPOSITION 4
text: Proof. We want to show that E[∇ZL(T̂ k, Z)] = ∇ZL(T k, Z). Since the terms of L1 are unaffected by T̂ , they are excluded w.l.g. from L in the proof.
E[∇ZL(T̂ k, Z)] = E [ ∇Z(− ∑ u,v∈V ∑ k∈{1..C} L2(T̂ k, u, v)L3(Z, u, v)) ]
(by linearity of expectation) = −∇Z ∑ u,v∈V ∑ k∈{1..C} L2(E[T̂ k], u, v)L3(Z, u, v)
(by Prop 1) = −∇Z ∑ u,v∈V ∑ k∈{1..C} L2(T k, u, v)L3(Z, u, v) = ∇ZL(T k, Z)
The following table gives the decomposition for DeepWalk, node2vec, and Watch Your Step. Node2vec also introduces a biased sampling procedure based on hyperparameters (they name p and q) instead of uniform transition probabilities. We can equivalently bias the transitions in GTTF to match node2vec’s. This would then show up as a change in T̂ k in the objective. This effect can also be included in the objective by multiplying 〈Zu, Zv〉 by the probability of such a transition in L3. In this format, the p and q variables appear in the objective and can be included in the optimization. For WYS, Qk are also trainable parameters.
For methods in which the transition distribution is not uniform, such as node2vec, there are two options for incorporating this distribution in the loss. The obvious choice is to sample from a biased transition matrix, Tu,v = W̃u,v , where W̃ is the transition weights. Alternatively, the transition bias can be used as a weight on the objective itself. This approach is still unbiased as
E v∼W̃u [L(v, u)] = ∑ v∈V Pv∼W̃u [v]L(v, u) = ∑ v∈V W̃u,vL(v, u)

heading: B.5 PROOF OF PROPOSITION 5
text: Proof. Let Ã be the neighborhood patch returned by GTTF, and let .̃ indicate a measurement based on the sampled graph, Ã, such as the degree vector, δ̃, or diagonal degree matrix, D̃. For the remainder of this proof, let all notation for adjacency matrices, A or Ã, and diagonal degree matrices, D or D̃, and degree vector, δ, refer to the corresponding measure on the graph with self loops e.g. A← A+ In×n. We now show that the expectation of the layer output is unbiased.
E [ ◦ AkH(l−1)W (l) ] = [ E [ ◦ Ak ] H(l−1)W (l) ] implies that E [ ◦ AkH(l−1)W (l) ] is unbiased if
E [ ◦ Ak ] = ( D−1/2AD−1/2 )k .
E [ ◦ Ak ] = E [ D1/2 ( D̃−1Ã )k D−1/2 ] = D1/2E [( D̃−1Ã )k] D−1/2
Let Pu,v,k be the set of all walks {p = (u, v1, ..., vk−1, v)|vi ∈ V }, and let p∃Ã indicate that the path p exists in the graph given by Ã. Let tu,v,k be the transition probability from u to v in k steps, and let tp be the probability of a random walker traversing the graph along path p.
E [( D̃−1Ã )k u,v ] = E [ T̃ ku,v ] = Pr [ t̃u,v,k = 1 ] = ∑ p∈Pu,v,k Pr[p∃Ã]Pr[t̃p = 1|p∃Ã]
= ∑
p∈Pu,v,k k∏ i=1 1[A[pi, pi+1] = 1] f + 1 δ[pi] k∏ i=1 (f + 1)−1 = ∑ p∈Pu,v,k k∏ i=1 1[A[pi, pi+1] = 1]δ[pi] −1
= ∑
p∈Pu,v,k Pr[tp = 1] = Pr[tu,v,k] = (T )ku,v =
( D−1A )k u,v
Thus, E [ ◦ Ak ] = ( D−1/2AD−1/2 )k and E [ ◦ AkH(l−1)W (l) ] = ( D−1/2AD−1/2 )k H(l−1)W (l)
For writing, we assumed nodes have degree, δu ≥ f , though the proof still holds if that is not the case as the probability of an outgoing edge being present from u becomes 1 and the transition probability becomes δ−1u i.e. the same as no estimate at all.

heading: B.6 PROOF OF PROPOSITION 6
text: GTTF can be seen as a way of applying dropout (Srivastava et al., 2014), and our proof is contingent on the convergence of dropout, which is shown in Baldi & Sadowski (2014). Our dropout is on the adjacency, rather than the features. Denote the output of a graph convolution network5 with H:
H = GCNX(A;W ) = T XW
We restrict our analysis to GCNs with linear activations. We are interested in quantifying the change of H as A changes, and therefore the fixed (always visible) features X is placed on the subscript. Let Ã denote adjacency accumulated by GTTF’s ROOTEDADJACC (Eq. 1).
H̃c = GCNX(Ãc).
Let A = {Ãc}|A|c=1 denote the (countable) set of all adjacency matrices realizable by GTTF. For the analysis, assume the graph is α-regular: the assumption eases the notation though it is not needed. Therefore, degree δu = α for all u ∈ V . Our analysis depends6 on 1|A| ∑ Ã∈A Ã ∝ A. i.e. the average realizable matrix by GTTF is proportional (entry-wise) to the full adjacency. This is can be shown when considering one-row at a time: given node u with δu = α outgoing neighbors, each of its neighbors has the same appearance probability = 1δu . Summing over all combinations ( δu f ) , makes
each edge appear the same frequency = 1δu |A|, noting that |A| evenly divides ( δu f ) for all u ∈ V .
We define a dropout module:
d
A = |A|∑ c zcÃc with z ∼ Categorical
( |A| of them︷ ︸︸ ︷ 1
|A| ,
1
|A| , . . . ,
1
|A|
) , (10)
where zc acts as Multinoulli selector over the elements of A, with one of its entries set to 1 and all others to zero. With this definitions, GCNs can be seen in the droupout framework as: H̃ = GCNX( d
A). Nonetheless, in order to inherit the analysis of (Baldi & Sadowski, 2014, see their equations 140 & 141), we need to satisfy two conditions which their analysis is founded upon:
(i) E[GCNX( d
A)] = GCNX(A): in the usual (feature-wise) dropout, such condition is easily verified. (ii) Backpropagated error signal does not vary too much around around the mean, across all realizations of d A.
Condition (i) is satisfied due to proof of Proposition 5. To analyze the error signal, i.e. the gradient of the error w.r.t. the network, assume loss function L(H), outputs scalar loss, is λ-Lipschitz continuous.
5The following definition averages the node features (uses non-symmetric normalization) and appears in multiple GCN’s including Hamilton et al. (2017).
6If not α-regular, it would be 1|A| ∑ Ã∈A Ã ∝ D −1A
The Liptchitz continuity allows us to bound the difference in error signal between L(H) and L(H̃):
||∇HL(H)−∇HL(H̃)||22 (a) ≤λ (∇HL(H)−∇HL(H̃))>(H − H̃) (11) (b)
≤λ ||∇HL(H)−∇HL(H̃)||2 ||H − H̃||2 (12) w.p.≥1− 1
Q2 ≤ λ ||∇HL(H)−∇HL(H̃)||2 W>X>Q √ V ar[T ]XW (13)
= λQ
2 √ f ||∇HL(H)−∇HL(H̃)||2 ||W ||21 ||X||21 (14)
||∇HL(H)−∇HL(H̃)||2 ≤ λQ
2 √ f ||W ||21 ||X||21 (15)
where (a) is by Lipschitz continuity, (b) is by Cauchy–Schwarz inequality, “w.p.” means with probability and uses Chebyshev’s inequality, with the following equality because the variance of T is shown element-wise in proof for Prop. 2. Finally, we get the last line by dividing both sides over the common term. This shows that one can make the error signal for the different realizations arbitrarily small, for example, by choosing a larger fanout value or putting (convex) norm constraints on W and X e.g. through batchnorm and/or weightnorm. Since we can have ∇HL(H) ≈ ∇HL(H̃1) ≈ ∇HL(H̃2) ≈ · · · ≈ ∇HL(H̃|A|) with high probability, then the analysis of Baldi & Sadowski (2014) applies. Effectively, it can be thought of as an online learning algorithm where the elements of A are the stochastic training examples and analyzed per (Bottou, 1998; 2004), as explained by Baldi & Sadowski (2014) .

heading: B.7 PROOF OF PROPOSITION 7
text: The storage complexity of CompactAdj is O(sizeof(δ) + sizeof(Â)) = O(n+m). Moreover, for extemely large graphs, the adjacncy can be row-wise partitioned across multiple machines and therefore admitting linear scaling. However, we acknolwedge that choosing which rows to partition to which machines can drastically affect the performance. Balanced partitioning is ideal. It is an NP-hard problem, but many approximations have been proposed. Nonetheless, reducing inter-communication, when distributing the data structure across machines, is outside our scope.

heading: B.8 PROOF OF PROPOSITION 8
text: For each step of GTTF, the computational complexity is O(bhf ). This follows trivially from the GTTF functional: each nodes in batch (b of them) builds a tree with depth h and fanout f i.e. with hf tree nodes. This calculation assumes random number generation, ACCUMULATEFN and BIASFN take constant time. The searchsorted function is linear, as it is called on a sorted list: cumulative sum of probabilities.

heading: C ADDITIONAL GTTF IMPLEMENTATIONS

heading: C.1 MESSAGE PASSING IMPLEMENTATIONS

heading: C.1.1 GRAPH ATTENTION NETWORKS (GAT, VELIČKOVIĆ ET AL., 2018)
text: One can implement GAT by following the previous subsection, utilizing ACCUMULATEFN and BIASFN defined in (1) and (2), but just replacing the model (3) by GAT’s:
GAT(Ã,X;A,W1,W2) = softmax((A ◦ ◦ A) ReLu((A ◦ ◦ A)XW1)W2); (16)
where ◦ is hadamard product and A is an n× n matrix placing a positive scalar (an attention value) on each edge, parametrized by multi-headed attention described in (Veličković et al., 2018). However, for some high-degree nodes that put most of the attention weight on a small subset of their neighbors, sampling uniformly (with BIASFN=NOREVISITBIAS) might mostly sample neighbors with entries inAwith value≈ 0, and could require more epochs for convergence. However, our flexible functional
allows us to propose a sample-efficient alternative, that is in expectation, equivalent to the above:
GAT(Ã,X;A,W1,W2) = softmax(( √ A ◦ ◦ A) ReLu(( √ A ◦ ◦ A)XW1)W2); (17)
def GATBIAS(T, u): return NOREVISITBIAS(T, u) ◦ √ A[u, Â[u]]; (18)

heading: C.1.2 DEEP GRAPH INFOMAX (DGI, VELIČKOVIĆ ET AL., 2019)
text: DGI implementation on GTTF can use ACCUMULATEFN=ROOTEDADJACC, defined in (1). To create the positive graph: it can sample some nodes B ⊂ V . It would pass to GTTF’s Traverse B, and utilize the accumulated adjacency Â for running: GCN(Â,XB) and GCN(Â,Xpermute), where the second run randomly permutes the order of nodes in X . Finally, the output of those GCNs can then be fed into a readout function which outputs to a descriminator trying to classify if the readout latent vector correspond to the real, or the permuted features.

heading: C.2 NODE EMBEDDING IMPLEMENTATIONS

heading: C.2.1 NODE2VEC (GROVER & LESKOVEC, 2016)
text: A simple implementation follows from above: N2VACC , DEEPWALKACC; but override BIASFN =
def N2VBIAS(T, u): return p−1[i=T−2]q−1[〈A[T−2],A[u]〉>0]; (19)
where 1 denotes indicator function, p, q > 0 are hyperparameters of node2vec assigning (unnormalized) probabilities for transitioning back to the previous node or to node connected to it. 〈A[T−2], A[u]〉 counts mutual neighbors between considered node u and previous T−2. An alternative implementation is to not override BIASFN but rather fold it into ACCUMULATEFN, as:
def N2VACC(T, u, f): DEEPWALKACC(T, u, f); η[u] ← η[u] × N2VBIAS(T, u); (20) Both alternatives are equivalent in expectation. However, the latter directly exposes the parameters p and q to the objective L: allowing them to be differentiable w.r.t. L and therefore trainable via gradient descent, rather than by grid-search. Nonetheless, parameterizing p & q is beyond our scope.

heading: C.2.2 WATCH YOUR STEP (WYS, ABU-EL-HAIJA ET AL., 2018)
text: First, embedding dictionaries R,L ∈ Rn× d2 can be initialized to random. Then repeatedly over batches B ⊆ V , the loss L can be initialized to estimate the negative part of the objective:
L ← − ∑ u∈B log σ(−Ev∈U(V ) [〈Ru, Lv〉+ 〈Rv, Lu〉]),
Then call GTTF’s traverse passing the following ACCUMULATEFN=
def WYSACC(T, u): if T.size() 6= Q.size(): return; t← T [0]; U ← T [1 :] ∪ [u];
ctx_weighted_L← ∑ j QjLUj ; ctx_weighted_R← ∑ j QjRUj ;
L ← L− log(σ(〈Rt,ctx_weighted_L〉+ 〈Lt,ctx_weighted_R〉));

heading: D MISCELLANEOUS

heading: D.1 SENSITIVITY
text: The following figures show the sensitivity of fanout and walk depth for WYS on the Reddit dataset.

heading: D.2 RUNTIME OF GCN ON CITATION NETWORKS
text: Citation networks (Cora, Citeseer, Pubmed) have been popularized by Kipf & Welling (2017), and for completeness, we report in Figure 4 the runtime versus test accuracy of GCN on these networks. We compare against PyG, which is optimized for GCN. We find that the methods are somewhat comparable in terms of training time.

Review summaries:
The authors propose a "meta-algorithm" for approximating various graph representation learning schemes: generate batches of random trees with fixed fanout (and possibly biased probabilities of selecting different edges), and use them to accumulate information to approximate operations on the graph.  The idea is beautifully simple, and generalizes running independent random walkers, an approach that is used in deriving many related algorithms.  The biasing and accumulation functions are user provided, and the authors show how different choices of these functions can be used to approximate a number of graph representation learning schemes.  The authors also provide a software framework, though it was inaccessible at review due to anonymization.  Experiments show the approach is much more scalable than competing approaches (though, to be fair, some of the competition was not targeting scalability).