State Space Models on Temporal Graphs: A First-Principles Study

Authors: jintang  li
ruofan  wu
xinzhou  jin
boqun  ma
liang  chen
zibin  zheng
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-14 14:00:37

Abstract: Over the past few years, research on deep graph learning has shifted from static graphs to temporal graphs in response to real-world complex systems that exhibit dynamic behaviors. In practice, temporal graphs are formalized as an ordered sequence of static graph snapshots observed at discrete time points. Sequence models such as RNNs or Transformers have long been the predominant backbone networks for modeling such temporal graphs. Yet, despite the promising results, RNNs struggle with long-range dependencies, while transformers are burdened by quadratic computational complexity. Recently, state space models (SSMs), which are framed as discretized representations of an underlying continuous-time linear dynamical system, have garnered substantial attention and achieved breakthrough advancements in independent sequence modeling. In this work, we undertake a principled investigation that extends SSM theory to temporal graphs by integrating structural information into the online approximation objective via the adoption of a Laplacian regularization term. The emergent continuous-time system introduces novel algorithmic challenges, thereby necessitating our development of GraphSSM, a graph state space model for modeling the dynamics of temporal graphs. Extensive experimental results demonstrate the effectiveness of our GraphSSM framework across various temporal graph benchmarks.

Main paper: 
 introduction: 
 As a class of neural networks designed to operate directly on graph-structured data, graph neural networks (GNNs) [21, 43, 24] have achieved remarkable success and have established new state-ofthe-art performance across a broad spectrum of graph-based learning tasks [19] . While significant progress has been made in researching static graphs, many real-world networks, such as social, traffic, and financial networks may exhibit temporal behaviors that carry valuable time information [20, 23] . This gives rise to temporal (dynamic) graphs, wherein the nodes and edges of the graph may undergo constant or periodic changes over time. In applications where temporal graphs arise, modeling and exploiting the dynamic nature of the continuously evolving graph is crucial in representing the underlying data and achieving high predictive performance [22, 44, 41] .Learning over temporal graphs is typically approached as a sequence modeling problem in which graph snapshots form a sequence [34] . This often involves challenges related to long graph sequences and scalability issues [25] . Recurrent neural networks (RNNs) [46, 4, 18] have historically dominated sequence modeling over the last years. However, they have long been plagued by poor capability in modeling long sequences due to rapid forgetting. This hampers their performance in temporal graphs that require a broader context or longer time window to capture relevant dependencies and patterns. Recently, the advancement of Transformers [42] has led to a shift in this paradigm, given their superior performance. Yet, Transformers also struggle with long sequence learning because the computational and memory complexity of self-attention is quadratically dependent on the sequence length. The overwhelming computation and memory requirements/costs associated with Transformers makes them less applicable in practical applications handling long-term sequences [51] .Recently, state space models (SSMs) have emerged as a powerful tool for sequence modeling [11, 13, 29, 40, 6, 10] . The salient characteristic that distinguishes state space models as particularly compelling is their conceptualization of sequential inputs as discrete observations from an underlying process evolving in continuous time, which naturally arises in scenarios such as speech processing [40] and time series analysis [49] . SSMs sustain a latent state throughout an input sequence and formulate state update equations through the discretization of an underlying linear dynamical system (LDS). Owing to their invariant state size, SSMs exhibit an efficient inferential time complexity, akin to that of RNNs. Simultaneously, they overcome the long-range modeling deficiencies inherent to RNNs through meticulous initializations of state matrices which are theoretically shown to achieve an optimal compression of history [11] .Temporal graphs often manifest as discrete snapshots capturing the evolution of an underlying graph that is inherently dynamic and continuous in nature [20] . In this context, the SSM methodology could be appropriated as a foundational primitive for temporal graph modeling. However, SSMs are predominantly architected towards independent sequence modeling. Hence, the task of systematically incorporating time-varying structural information into the SSM framework poses significant challenges. Specifically, it remains unexplored as to whether the foundational methodology of discretized LDS is readily applicable to the domain of temporal graphs.In this work, we advance the SSM methodology to encompass temporal graphs from the first principles. Rather than presupposing the evolution of the underlying temporal graph, we dive into the fundamental problem of online function approximation that underpins the theoretical development of SSMs for sequence modeling [11] . By solving a novel Laplacian regularized online approximation objective, we derive a piecewise dynamical system that compresses historical information of temporal graphs. The piecewise nature of the obtained continuous-time system poses new challenges toward discretization into linear recurrences, thereby motivating our design of GRAPHSSM, a state space framework for temporal graphs. The main contributions of this work are summarized as follows:• We introduce the GHIPPO abstraction, a novel construct predicated on the objective of Laplacian regularized online function approximation. This abstraction can alternatively be conceptualized as a memory compression primitive that simultaneously compresses both the feature dynamics and the evolving topological structure of the underlying temporal graph. The solution to GHIPPO is characterized by a dynamical system that is piecewise linear in node feature inputs. • We introduce GRAPHSSM, a flexible state space framework designed for temporal graphs, which effectively addresses the key algorithmic challenge of unobserved graph mutations that impedes the straightforward discretization of the GHIPPO solution into (linear) recurrences through employing a novel mixed discretization strategy. • Experimental results on six temporal graphs have validated the effectiveness of GRAPHSSM. In particular, GRAPHSSM has the advantages in scaling efficiency compared to existing state-of-thearts, which can generalize to temporal graphs with long-range snapshots.2 Related workWe consider two types of mixing mechanisms: convolution (CONV1D) and Scaled interpolation (INTERP) which we describe below: CONV1D. This is the usual convolution operation along the sequence dimension using shared parameters. We use a kernel size of 2 so that only consecutive representations are mixed.INTERP. This is an input-dependent weighted average strategy followed by an input-dependent scaling, implemented asEQUATIONwhere Z 1 , Z 2 ∈ R N V ×d are node representation matrices corresponding to consecutive snapshots. ρ and ξ are scale and weight functions that map two inputs into positive real numbers of identical shape with Z 1 or Z 2 , defined byEQUATIONwhere W ρ , W ξ ∈ R 2d×d and b ρ , b ξ ∈ R d are learnable parameters therein. 
temporal graph learning: 
 A major branch of temporal graph learning methods consists of snapshot-based methods, which handle discrete-time temporal graphs by learning the temporal dependencies across a sequence of time-stamped graphs. Early works mainly focus on learning node representations by simulating temporal random walks [39] or modeling the triadic closure process [50] on multiple graph snapshots. These methods typically generate piecewise constant representations and may suffer from the staleness problem [20] . In recent years, the most established solution has been switched to combine sequence models (e.g., RNNs [46] and SNNs [38, 8] ) with static GNNs to capture temporal dependencies and correlations between snapshots [47, 34, 39, 25] . To better translate the success achieved on static graphs in both their design and training strategies, recent frameworks such as ROLAND [48] and its variants [53, 17] have been proposed to repurpose static GNNs to temporal graphs. There is another important line of research that focuses on continuous-time temporal graphs, we kindly refer readers to [27] and [20] for comprehensive surveys on this research topic.State space models (SSMs) have historically served as a pivotal tool in fields such as signal processing [30] and time series analysis [3] . In recent advancements, they have also seen active adoption as a layer within neural sequence modeling frameworks [11, 13, 29, 40, 10] . The linear nature of SSMs confers several significant advantages. Key among these is the better-controlled stability that enables effective long-range modeling through careful initializations of state space layer parameters [11, 31] , with the most representative method being HIPPO [11] , a theory-driven framework notable for its optimal memory compression on continuous sequence inputs. Moreover, the computational efficacy of SSMs is notably enhanced through the use of techniques such as convolutions [13, 6] or parallel scans [40] . The promising properties of SSMs also attracts further explorations on graphs [2] .Comparison. The usual paradigms for designing sequence models over graphs involve recurrence (e.g. RNNs [46] ), integrate-and-fire (e.g. SNNs [38, 8] ), or attention (e.g. Transformers [42] ), which each come with tradeoffs [14] . For example, RNNs are a natural recurrence model for sequential modeling that require only constant computation/storage per time step, but are slow to train and suffer from the rapid forgetting issue. This empirically limits their ability to handle long sequences. SNNs share a similar recurrent architecture with RNNs while using 1-bit spikes to transmit temporal information, which would sacrifice expressivity and potentially suffer from optimization difficulties (e.g., the "vanishing gradient problem") [26] . Transformers encode local context via attention mechanism and enjoy fast, parallelizable training, but are not sequential, resulting in more expensive inference and an inherent limitation on the context length. Compared to the aforementioned architectures, SSMs particularly the promising Mamba (S6) model [10] , offer advantages such as fast training and inference, along with fewer parameters and comparable performance. These characteristics make SSMs particularly well-suited for sequence modeling, even (or especially) on extremely long sequences. Comparisons among these architectures are illustrated in table 1 Table 1 : Comparisons of different neural network architectures for sequence modeling.RNNs [46, 4, 18] SNNs [38, 8] Transformers [42] SSMs (S6 [10] ) 3 The GRAPHSSM frameworkSSM Graph Convolution + × K Linear SSM Graph Convolution + Output × K … … Feature Mixing 𝐺0 𝐺𝐿 𝐺1 𝐺𝐿-1 SSM Graph Convolution ∈ ∈ GCN, GAT, SAGE, … S4, S5, S6, … 𝐻 1:𝐿 (𝑘) Figure 1: GRAPHSSM framework.The primary motivation of our framework is the fact that discrete-time temporal graphs are sequential observations of an underlying temporal graph that evolves continuously. Adopting this functional viewpoint, we will first develop a piecewise recurrent memory update scheme in section 3.1 that optimally approximates the underlying continuoustime temporal graph, utilizing a novel extension of the HIPPO abstraction to graph-typed inputs [11] . The proposed framework retains many nice properties of HIPPO while posing the new challenge of unobserved graph mutation when handling discretely-observed observations, which we analyze in section 3.2 and propose a mixing mechanism to improve the recurrent approximation. Finally, we present GRAPHSSM framework in section 3.3. An overview of GRAPHSSM is shown in figure 1 . 
ghippo: hippo on temporal graphs: 
 Setup. We fix a time interval [0, T ]. A temporal graph on [0, T ] is characterized by two processes G and X: For each t ∈ [0, T ], the process G maps t to a graph object G(t) = (V (t), E(t)). We assume the node process V (t) to be fixed over time, i.e., V (t) ≡ V, t ∈ [0, T ] with N V = |V | and discuss the case for varying node processes in appendix B.2. The edge process E(t) is a piecewise-constant process with a finite number M of mutations over [0, T ] that are described via a sequence of events:EQUATIONEach event E m constitutes an interaction between node pair (u m , v m ) at time t m with action a m , the action could be either insertion or deletion. The evolution process is thus depicted as the following:EQUATIONThe process X maps t to a node feature matrix X(t) ∈ R N V ×d with feature dimension d. Throughout this paper, it is often helpful to view G and X as graph-valued and matrix-valued functions. In typical discrete-time temporal graph learning problems, the underlying graph is observed at timestamps τ 1 , . . . , τ L with time gaps∆ l = τ l -τ l-1 , 2 ≤ l ≤ L.The observations thus form a sequence of snapshots {G(τ l ), X(τ l )} 1≤l≤L which are abbreviated as {G 1:L , X 1:L }. Notably, the observation times are usually interleaved with the mutation times, resulting in the majority of mutation times remain unobserved. This situation presents significant challenges in effectively modeling the dynamics of graph evolution, a topic that will be further explored subsequently.The HIPPO abstraction. Algorithmically, the goal of continuous-time dynamic modeling is to design a memory module that optimally compresses all the historical information [11] . Under the context of univariate sequence modeling, the HIPPO framework [11] formalizes the memory compression problem into an online approximation problem in some function space and derives HIPPO operators under specific types of basis functions, among which the HIPPO-LEGS configuration has become the state-of-the-art in state-space sequence modeling paradigms [13, 40] . However, naively extending HIPPO abstraction to graph learning scenarios (via treating node features as inputs) could be deemed inadequate since HIPPO handles distinct inputs independently, without the capability to incorporate the interconnectivity information among various inputs which could potentially enhance the efficiency of memory compression. For illustrative purposes, in instances where input observations are noisy, the exploitation of neighborhood information has the potential to facilitate a denoising step, as evidenced in image processing applications [33] and semi-supervised learning primitives [52, 45] . To systematically utilize the connectivity information, we propose a new approximation paradigm, the Laplacian-regularized online approximation that extends HIPPO to graph modeling frameworks. Formally, we start with the simple setup with d = 1, i.e., each node possesses a scalar feature, and we propose an approximation scheme that simultaneously approximates the history of all the N V inputs up until time t, i.e., {X(s), s ∈ [0, t]} using their corresponding memories at time t, i.e., Z(t) = {z v (t)} v∈V ∈ R N V ×1 according to the following objective at time t:EQUATIONHere α > 0 is a balancing constant, µ t is a time-dependent measure that is supported on the interval [0, t] which controls the importance of various parts of the input domain 3 and L(t) is a normalized Laplacian at time t, which allows definition such as the symmetric normalized Laplacian L sym (s) = I -D(s) -1/2 A(s)D(s) -1/2 where D(s) is a diagonal matrix whose diagonals are node degrees, or random walk normalized Laplacian L rw (s) = I -D(s) -1 A(s). The objective (3) is understood as the ordinary HIPPO approximation objective augmented with a regularization component that encourages the smoothness of memory compression with respect to adjacent nodes. 4The imposition of smoothness constraints commonly emerges as a beneficial relational inductive bias in the context of graph learning [1] . By leveraging the data from adjacent nodes, one can potentially achieve a more effective denoising effect during the process of node memory compression. To specify a suitable approximation subspace for memories Z, we adopt the approach of HIPPO that uses some N -dimensional subspace of polynomials which we denote as P N . Now we define a graph memory projection operator GPROJ t that maps the temporal graph up until time t to a collection of N V polynomials with each one lies in P N , i.e., GPROJ t (G, X) = arg minEQUATIONWe further define a coefficient operator COEF t that maps each polynomial in the collection in (4) to the coefficients of the basis of orthogonal polynomials defined with respect to µ t , the following definition formalizes our extension of HIPPO to continuous-time temporal graphs which we term GHIPPO: Definition 1 (GHIPPO). Given a continuous-time temporal graph (G, X), a time-varying measure family µ t , an N -dimensional subspace of polynomials P N , the GHIPPO operator at time t is the composition of GPROJ t and COEF t that maps the temporal graph and node features to a collection of projection coefficientsU (t) ∈ R N V ×N , or GHIPPO (G, X) = COEF t (GPROJ t (G, X)).The most favorable property of the HIPPO framework on independent inputs is that the outputs of HIPPO operators are characterized via a concise ordinary differential equation (ODE) that takes the form of a linear time-invariant state space model (LTI-SSM). The following theorem states that most of the desirable properties of HIPPO are retained by GHIPPO except for the LTI property: Theorem 1. Let G evolve according to (2) . Taking µ t to be the scaled Legendre measure (LegS) with µ t = 1 t I [0,t] where I [0,t] stands for the indicator function of the interval [0, t], the evolution of the outputs of GHIPPO operator is characterized by M ODEs according to mutation times as follows:EQUATIONwhere A ∈ R N ×N and B ∈ R N ×1 takes the same form as in the HIPPO formulation [11] :EQUATIONEQUATIONAccording to theorem 1, the solution (5) is LTI over each interval [t m , t m+1 ) during which the graph structure remains fixed. This property further extends to a piecewise LTI perspective over the interval [0, T ]. Moreover, we may view the solution (5) as a two-stage procedure that could be intuitively described as diffuse-then-update. Specifically, this procedure entails a sequential execution, wherein an initial diffusion operation is applied to the features of the input nodes, succeeded by an update to the memory of these nodes.Theorem 1 establishes an analogue of HIPPO theory on temporal graphs. It is straightforward to verify that most of the subsequent refinements of HIPPO apply to GHIPPO as well. Among these we will utilize the popular technique of diagonal state spaces [15, 12] that simply sets A as a diagonal matrix with negative diagonal elements5 . To apply the GHIPPO framework to discrete-time temporal graphs, a critical step is to develop a discretized version of (5). However, unlike ordinary HIPPO where we can use standard discretization techniques of ODEs to discretize LTI equations, the GHIPPO ODE contains discontinuities that correspond to mutation times of the underlying temporal graph, which are often not observed given only access to a list of snapshots. This issue of unobserved dynamics complicates the development of a viable discretization scheme for GHIPPO, as is pictorially illustrated in figure 2. To devise a solution to this challenge, we start by analyzing a hypothetical oracle scenario in which all mutations are observable.An oracle discretization. We consider a time range [τ l-1 , τ l ) between the l -1th and the lth snapshot, and assume there are altogether M l mutation events {E l,i } 1≤i≤M l happened during this period. Let G l,0 = G l-1 be the graph snapshot at τ l-1 , the following process describes the structural evolution inside the interval [τ l-1 , τ l ):G l-1 = G l,0 E l,1 -→ G l,1 E l,2 -→ G l,2 -→ • • • -→ G l,M l -1 E l,M l -→ G l,M l = G l (7)Next, we derive a discretization formula under the strategy of zeroth-order-hold (ZOH). We assume that all intermediate mutations are observed, with the node features staying fixed between mutations, i.e., X(t) ≡ X l,i , t ∈ [t l,i-1 , t l,i ). The following theorem characterizes the resulting state evolution: Theorem 2 (Oracle discretization of ( 5)). Assume A is a diagonal matrix with negative diagonals, for any 1 ≤ l ≤ L. Let L l,i be some Laplacian of G l,i , we have the following oracle update rule:EQUATIONwhereU l ∈ R N V ×N denotes the discretized state at step l with U 0 = 0. For each 1 ≤ l ≤ L, 0 ≤ i ≤ M l , Λ i ∈ R N ×Nare non-negative diagonal matrices with values depending only on the mutation times, which satisfy (8) , given all the (unobserved) mutation information, the state update rule is equivalent to applying ZOH to X l which is an element-wise convex combination of all the diffused node features. In practice, among all the components of X, we only have access to X l-1 , X l , G l-1 , G l with the rest left unobserved. Therefore, we propose mixed discretization as an approach to approximate X l . Specifically, we introduce the following mechanisms:M l i=0 Λ i = I. Mixed discretization. According toX (O) l = GNN θ (X l , G l ) ,(ordinary ZOH)EQUATIONEQUATION) which are compositions of inter-node mixing (a consequence of diffusion) and intra-node mixing (mixing node features of consecutive snapshots). For the process of inter-node mixing, we opt to approximate the diffusion operation with a learnable shallow graph neural network (typically a 1-layer GNN) parameterized by θ to alleviate the computation burden and improve flexibility 6 . A detailed discussion considering the relation between certain GNN formulations and the choice of Laplacian is presented in appendix B.1. In the context of intra-node mixing, we introduce a MIX module parameterized by ϕ to merge either consecutive node features (as illustrated in (feature mixing)) or consecutive node representations produced by the GNN model (as illustrated in (representation mixing)). In this paper, we assess two simple MIX instantiations: Convolution with a kernel size of 2 (CONV1D) and a gating mechanism that interpolates between the two inputs (INTERP). We postpone a comprehensive description of the mixing methods to appendix D.1. The resulting discretized system is presented as the following matrix-valued state space model:EQUATIONWhen exact timestamps for snapshots are unavailable, we use the adaptive time step strategy as in [11, 10] that models ∆ a 1-dimensional affine projection of the inputs followed by a non-negative activation like softplus. Finally, we utilize the approximation A -1 e ∆A -I ≈ ∆I for diagonal As, and equip the system with an output Y with a state projection matrix C ∈ R N ×1 .Having established the SSM equation ( 9), we are ready to introduce our main framework GRAPHSSM.In alignment with conventional design paradigms in the SSM literature, we define a depth-K GRAPHSSM model through the sequential composition of K GRAPHSSM blocks, with each block characterized as follows:EQUATIONwhere we useH (k)1:L to denote the concatenation of the hidden representation at depth k of all the snapshots along the sequence dimension and H (0) 1:L are the node features X 1:L . The GRAPHSSM blocks, as outlined in (10) , incorporate an SSM layer that operates on graph snapshot inputs. This is followed by the application of a nonlinear activation σ and the integration of a residual connection which we denote as the addition of a linear projection of inputs with Linear denotes a linear projection layer that ensures dimension compatibility.GRAPHSSM-S4. The architectural formulation of the SSM layer essentially involves the expansion of the one-dimensional recurrence, as specified in (9) , to accommodate general dimensions, i.e., d > 1. This expansion is achieved in a straightforward manner by utilizing an individual SSM for each dimension. Consequently, the emergent SSM layer adopts a Single-Input, Single-Output (SISO) configuration. Such a design is intuitively understood as the graph learning analogue of S4 [13] , which we term GRAPHSSM-S4.GRAPHSSM-S5 and GRAPHSSM-S6. In addition to the SISO implementation, we further introduce two variants within the GRAPHSSM framework. The first alternative represents a Multiple-Input, Multiple-Output (MIMO) extension of (9) , wherein a single SSM system is applied across all dimensions. This variant serves as a graph-informed analogue to the S5 model [40] . The second variant extends the S4 model by facilitating input-controlled time intervals and state matrices (∆, B, and C). This innovation yields a selective state space model, drawing parallels to the latest SSM architectures such as S6 [10] .A detailed exposition of the GRAPHSSM-S4 (resp. GRAPHSSM-S5, GRAPHSSM-S6) layer is provided in algorithm 1 (resp. algorithm 2, algorithm 3) in appendix D.2. The overall end-to-end architecture is briefly illustrated in figure 1 , where we use feature mixing as the mixing mechanism for illustration. Remark 1 (Choice of mixing mechanisms). In the GRAPHSSM architecture, each SSM layer incorporates a mixing mechanism. Based on our empirical investigations, we have observed that employing more sophisticated mixing strategies such as (feature mixing) and (representation mixing), yields benefits predominantly when these are applied exclusively to the lowermost layer. Specifically, this entails utilizing either X (F) l or X (R) l configurations in the initial layer, while defaulting to X (O) l for the layers that follow. An intuitive rationale behind this strategic layer-specific choice will be elucidated in appendix B.3. 
experiments: 
 This section presents our key experimental findings on the temporal node classification task. Also, ablation studies of the key design choices are presented. Due to space limitation, the detailed experimental settings are deferred to appendix F.Node classification performance. The node classification performance of all methods is presented in table 2. It has been observed that graph embedding methods, especially static ones, tend to underperform in most cases. This is expected since these methods are typically trained in an unsupervised manner, solely focusing on exploiting the graph structure. We note that continuous-time methods HTNE and M 2 DNE exhibit poor performance in DBLP-3, Brain, and Reddit even when compared to static methods. This indicates that continuous-time methods are not well-suited for handling discrete-time graphs, particularly in the absence of temporal continuity. As can also be observed from table 2, most temporal graph neural networks demonstrate good performance on DBLP-3 and Brain datasets, where the node labels are largely dominated by node attribute information [47] . However, for datasets like Reddit and DBLP-10, where graph topology information plays a more significant role in classification, the performance has notably degraded. This indicates that the baseline methods struggle to effectively capture the underlying evolving graph structure and exploit it for accurate classification. In contrast, our most performant architecture, GRAPHSSM-S4, exhibits an average improvement of 14% and 2% in Micro-F1 and Macro-F1 scores, respectively, compared to state-of-the-art baselines on the Reddit and DBLP-10 datasets. In addition, GRAPHSSM-S4 is a more preferable choice for long graph sequences, achieving new state-of-the-art performance on the DBLP-10 dataset. SSM architectures. As GRAPHSSM is a general framework that generalizes SSMs to temporal graphs, we conduct experiments on extending GRAPHSSM with different ad-hoc SSMs, including S5 [40] and S6 [10] . The node classification results on four datasets are shown in table 4 . By comparing different variants of GRAPHSSM, we can find that S4 is the best architecture for learning over temporal graph sequences. S5, being a simplified version of S4 with fewer parameters, achieves poor performance on all datasets. Notably, while S6 shows impressive performance in other modalities such as language or images [10, 51] , it is observed that they underperform when applied to graph sequences. This indicates that the selective mechanism may not be a good fit for graph data. Mixing mechanism. We assess the effectiveness of various mixing mechanisms introduced in section 3.2 through a series of experiments conducted using the S4 variant of GRAPHSSM. The analysis spans four distinct configurations: no intra-node mixing (X (O) 1 + X (O)2 ), feature mixing at the first layer (X (F) 1 + X (O)2 ), and representation mixing at either the first (X (R) 1 + X (O) 2 ) or second ( X (O) 1 + X (R)2 ) layers. The findings, presented in table 5 , indicate that the integration of the MIX module at the first layer generally leads to enhanced model performance. An intuitive explanation for this observed phenomenon is elaborated in appendix B.3. Initialization strategy. Recent advancements have highlighted the crucial role of initialization in SSMs [12] , prompting our investigation into the effects of various initialization strategies for the A matrix. Specifically, we explore "hippo", "constant", and "random" initializations, with their comprehensive definitions provided in appendix D.2. The result, as shown in figure 3 exhibits distinct performance variations across different initialization strategies, with HIPPO being typically the dominant one which corroborates our theoretical motivations. 
conclusion: 
 In this work, we introduce a conceptualized GHIPPO abstraction on temporal graphs. Building upon GHIPPO, we propose GRAPHSSM, a theoretically motivated state space framework for modeling temporal graphs derived from a novel memory compression scheme. The proposed framework is computationally efficient and versatile in its design, which is further corroborated by strong empirical performance across various benchmark datasets. We also point out the unobserved graph mutation issue in temporal graphs and propose different mixing mechanisms to ensure temporal continuity across consecutive graph snapshots. Despite the promising results, the applicability of GRAPHSSM is presently confined to discrete-time temporal graphs. A discussion of our framework's current limitations and the scope for future extensions is presented in appendix E.with corresponding penalties written asEQUATIONEQUATIONThe above display reveals the inductive bias of Laplacian regularizations as a promoting closeness in a weighted ℓ 2 metric regarding adjacent nodes' memory compressions, with distinct choices of Laplacians utilizing different weighting schemes. In particular, let α → ∞ in objective 3 then when the Laplacian is chosen as L sym , the solution Z sym (s) must satisfyEQUATIONIt then follows that Z sym compresses all the historical degree profiles over connected components of G. Analogously, when L rw is chosen, it follows that the solution Z rw (s) must satisfyEQUATIONwhich compresses the composition of connected components of G.Diffusion and GNN approximation. We consider approximations of the following diffused node features with respect to some type of Laplacian:EQUATIONThe right-hand side of the preceding display is equivalent to performing infinite rounds of message passing. If we drop most of the higher-order terms, we arrive at models similar to graph neural networks. In particular, we keep only the first order terms, i.e., k = 1, then for the two Laplacians listed above, for each v ∈ V , we have the resulting approximations:h sym v ≈ (1 -α)Bx v + u∈N (v) α √ d u d v Bf u (GCN-Like) h rw v ≈ (1 -α)Bx v + u∈N (v) α d u Bf u . (SAGE(MEAN)-Like)The above display exhibits a similar pattern to the design of graph neural networks with a aggregatethen-combine procedure, with the corresponding aggregation steps mirroring two typical GNN architectures GCN [21] and SAGE with mean pooling [16] . Furthermore, note that the effect of the balancing constant α would be absorbed into the learnable parameters of the GNN. 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [Yes] Justification: The discussion on both potential positive societal impacts and negative societal impacts of the work is provided in appendix A.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: This paper poses no such risks.Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [Yes]Justification: We have cited the original paper that produced the code package or dataset and have explicitly stated the license used for the open-source frameworks.Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. • If this information is not available online, the authors are encouraged to reach out to the asset's creators. 
new assets: 
 Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?Answer: [NA] Justification: This paper does not release new assets.Guidelines:• The answer NA means that the paper does not release new assets.• Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc. • The paper should discuss whether and how consent was obtained from people whose asset is used. • At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file. 
crowdsourcing and research with human subjects: 
 Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?Answer: [NA] Justification: This paper does not involve crowdsourcing nor research with human subjects.Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Including this information in the supplemental material is fine, but if the main contribution of the paper involves human subjects, then as much detail as possible should be included in the main paper. • According to the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector. 
institutional review board (irb) approvals or equivalent for research with human subjects: 
 Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?Answer: [NA] Justification: This paper does not involve crowdsourcing nor research with human subjects.Guidelines:• The answer NA means that the paper does not involve crowdsourcing nor research with human subjects. • Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: The methodology described in section 3.1 applies to temporal graphs with a fixed node set. To extend our approach to accommodate graphs featuring varying node sets, we initially focus on the continuous-time context, subsequently delving into discussions on discretization strategies. Suppose on the time interval T = [0, T ], the node set evolves as depicted in the following sequence:EQUATIONThat is, throughout the interval T , the node set undergoes alterations on R distinct occasions, with associated changes occurring at times t 1 , . . . , t R , respectively. We denote these evolving node sets as V 0 , . . . , V R . To systematically analyze this temporal evolution, we partition the entire interval T into R + 1 segments:EQUATION)According to the formulations in section 3.1, on each T r , we have a well defined GHIPPO operator and the solutions are characterized by theorem 1. With an approximation order of N , we let the resulting projection coefficients beEQUATIONTo address the issue of shape incoherence arising from variations in node sets, we employ a memory alignment procedure. This technique facilitates the mapping from U r (t r+1 -) to U r+1 (t r+1 ), ensuring that the memory associated with each node is aligned according to the following scheme:EQUATIONThe memory alignment procedure (21) retains the continuity of states for nodes that persist over time.For nodes that emerge anew within the graph, it assigns a default initial state, which could either be an all-zero state or an estimation derived a priori from the states of neighboring nodes.Discretizations. Within the established context, Theorem 2 remains applicable on each segment T r . Consequently, our primary concern becomes the treatment of nodes that emerge between consecutive snapshots. Adhering to the ZOH discretization rule, newly emerged nodes lack historical states and therefore do not undergo the MIX strategy, and use their initial state during their first appearance in the recurrent update. This initial state can be set to zero or determined through aggregation from neighboring nodes.The various mixing mechanisms introduced in this paper are designed to facilitate an estimation of a weighted average of unobserved graph representations that occur amidst successive observational time points. Starting with the output generated by the initial SSM block, these outputs inherently encapsulate the information pertaining to the current snapshot, as well as that of its antecedent. Thus, the incorporation of mixing mechanisms at a second-layer may inadvertently result in the assimilation of superfluous information, extending beyond the target scope of back-to-back snapshots. Therefore, confining the deployment of mixing solely to the first SSM layer ensures the strict conservation of temporal locality. We have empirically verified that such an approach yields enhanced performances.In this section we present the proof of theorem 1 and theorem 2. We first present some necessary technical preparations: For any t ∈ [0, T ], let µ t be some finite measure and let H µt denote the Hilbert space induced by the inner product ⟨f, g⟩ µt :=EQUATIONLet P N (t) be the space of polynomials constructed via the restriction of each element in P N to [0, t].We assume the measure family be chosen such that P N (t) ⊂ H µt , ∀t ∈ [0, T ]. For each v ∈ V , we assume that the restriction of x v (viewing as a function on [0, T ]) to [0, t] is an element of H µt . Note that these assumptions are trivially satisfied for the scaled Legendre measure (LegS)µ t = 1 t I [0,t] .Proof of theorem 1. Hereafter we omit the dependence on µ t and write the inner product simply as ⟨•, •⟩ without misunderstandings. Let P 0 , . . . , P N -1 be a set of orthogonal polynomials in P N with ⟨P i , P j ⟩ = 0 for i ̸ = j and the degree of P n is n for each 0 ≤ n ≤ N -1. Then for any f ∈ H µt , the optimal approximation in L 2 (µ t ) distance in P N is given byEQUATIONwhere we define Π to be the projection operator. Now we turn to L t (Z; G, X, µ), viewing x v as a function on [0, t] for any v, we have:EQUATIONEQUATIONEQUATIONEQUATIONThe preceding display suggest that the minimizer of L t (Z; G, X, µ) is the same as the minimizer of L t (Z; G, X, µ). It thus suffices to analyze L t (Z; G, X, µ) which is easier to work with since Π(x v ) ∈ P N , ∀v ∈ V and the solution is a direct application of Laplacian regularization with respect to the integrand at any s ∈ [0, t], yielding:EQUATIONNow let the coefficient matrixQ ∈ R N V ×N be defined as Q v,n = ⟨x v , P n ⟩, ∀v ∈ V, n ∈ [N ],we obtain the GHIPPO operator as:EQUATIONNext we take derivatives to the coefficients. Note that L(t) is discontinuous and we can only apply derivative on intervals where L(t) remains same. First note that if we choose µ t to be the scaled Legendre measure (LegS) with µ t = 1 t I [0,t] , and P n as basic Legengre polynomials [11, Appendix B.1.1], then we have the HIPPO property:EQUATIONwhereEQUATIONFix some 1 ≤ m ≤ M and for t ∈ [t m-1 , t m ) we have:EQUATION)which finishes the proof.Proof of theorem 2. For ease of presentation, we operate on the node level instead of graph level.Recall the unobserved dynamics:EQUATION)Following the assumptions, we can intuitively write the update process as follows:EQUATION)For each 0 ≤ i ≤ M l , let D i := (I + αL l,i ) -1 X l,i B ⊤ . Let d v,i be the v-th row of D i and u v,i be the v-th row of U l,i . We first write the ZOH update corresponding to each step in (36) for every v ∈ V :u v,i = e (ti-ti-1A) u v,i-1 + A -1 e (ti-ti-1A) -I d v,i , for 1 ≤ i ≤ M l u v,l for i = 0 (37)Next we do the recursion from the rightmost to the leftmost according to (8):EQUATION)where we defineΥ = A -1 e (τ l -t M l )A -I u v,M l + M l i=1e (τ l -ti)A A -1 e (ti-ti-1)A -I u v,i-1 (39) in the above display we define t 0 = τ l-1 . Note that A -1 and e Aβ are simultaneouly diagonalizable for any β, therefore the matrix multiplication commutes and we further writeEQUATIONWith some abuse of notation now we let A ∈ R N denote the diagonal vector of the matrix. We provide the following construction:λ i =       e (τ l -t M l )A -I e (τ l -τ l-1 )A -I i = M l e (τ l -ti-1)A -e (τ l -ti)A e (τ l -τ l-1 )A -IEQUATION)Here note that λ i ∈ R N . It is straightforward to verify that:EQUATIONwhere {λ i } 0≤i≤M l are non-negative N -dimensional vectors satisfying M l i=0 λ i = 1 N , with 1 N denoting the all-one vector of dimension N . As the values of λs are independent of v, the proof finishes by combining (38) , (42) and write the above conclusion in matrix form via settingΛ i = diag(λ i ), 0 ≤ i ≤ M l D Algorithm descriptions D.In this section, we elucidate on the methodology of GRAPHSSM through three specific instantiations.For clarity in our explanation, we employ certain notational conventions that might be somewhat different from the main text: the term V refers to the number of vertices in each graph snapshot G l within a sequence of L graph snapshots {G l } 1≤l≤L which we further denote as G 1:L , and D represents the dimensionality of node features. The symbol LINEAR is used to represent a linear projection layer including a bias term, where the dimensions for input and output are typically clear from the context to ensure compatibility. The notation X 1:L denotes the concatenation of L tensors of the same dimensions along their second axis. For operations on tensors of order higher than two, we use the einsum notation, as defined by the einops framework [37] . We present the algorithmic description of our design of SSM layers, namely GRAPHSSM-S4 (resp. GRAPHSSM-S5, GRAPHSSM-S6) in algorithm 1 (resp. algorithm 2, algorithm 3). 7 Subsequently, we adopt the Algorithm 1 GRAPHSSM-S4 layer Input: A sequence of graph (snapshots) G 1:L with each of size V . Node (hidden) feature inputsX 1:L ∈ R V ×L×D .A graph neural network GNN θ parameterized by θ.A mixing mechanism MIX ϕ parameterized by ϕ.A ∈ R D×N , B ∈ R D×N , C ∈ R D×N .A linear layer for adaptive time gaps LINEAR τ . Output:Y 1:L ∈ R V ×L×D 1: # Approximate diffusion via GNN 2: for t = 1 to L do 3: Z l = GNN θ (X l , G l );4:H l = Z l if l = 1 else MIX(Z l , Z l-1 ); 5: end for 6: Initialize state U 0 = 0; # SISO state of shape V × D × N 7: for t = 1 to L do 8: ∆ l = softplus (LINEAR τ (H l )); 9: A = exp (einsum(∆ l , A, "V, DN → V DN ")); 10: B = einsum(∆ l , B, "V, DN → V DN "); 11: U l = U l-1 ⊙ A + einsum(B, H l , "V DN, V D → V DN "); 12:Y l = einsum(U l , C, "V DN, DN → V D"); 13: end for; 14: return Y 1:L ; following neural architecture composed of K blocks, with each block composed of one SSM layer followed by nonlinear activation and a residual connection:EQUATIONwhere H1:L are the node features X 1:L . The SSMLAYER in (45) may be chosen as any of {GRAPHSSM-S4 , GRAPHSSM-S5, GRAPHSSM-S6}. In our implementation of GRAPHSSM-S6, we add an additional layer normalization as the last operation of each block.Input: A sequence of graph (snapshots) G 1:L with each of size V .Node (hidden) feature inputs X 1:L ∈ R V ×L×D . A graph neural network GNN θ parameterized by θ.A mixing mechanism MIX ϕ parameterized by ϕ.A ∈ R N ×1 , B ∈ R D×N , C ∈ R N ×D .A linear layer for adaptive time gaps LINEAR τ . Output:Y 1:L ∈ R V ×L×D 1: # Approximate diffusion via GNN 2: for t = 1 to L do 3: Z l = GNN θ (X l , G l );4:H l = Z l if l = 1 else MIX(Z l , Z l-1 ); 5: end for 6: Initialize state U 0 = 0; # MIMO state of shape V × N 7: for t = 1 to L do 8: ∆ l = softplus (LINEAR τ (H l )); 9: A = exp ∆ l A ⊤ ; 10: B = einsum(∆ l , B, "V, DN → V DN "); 11: U l = U l-1 ⊙ A + einsum(B, H l , "V DN, V D → V N "); 12:Y l = U l C; 13: end for; 14: return Y 1:L ; Initialization strategy. Recent developments in state space modeling have underscored the significance of initializing the state matrices A, B, and C, with the initialization of A frequently emerging as the most critical factor for the performance of the SSM [12] . Building upon the progress made in S4 [13] and S4D [29, 15] , we evaluate three disparate initialization strategies for the matrix A. Note that since A is diagonal, we instead represent A as a N -dimensional vector:EQUATIONS4D-Real (HIPPO) This is the diagonal part of the original HIPPO matrices (6).S4D-Const (Constant) This is the real part of the eigenvalues corresponding to the S4N matrix as defined in [13] , which equals -1 2 . Random This initialization is generated via a negative transform of a random number χ, which we generated using the Glorot initialization method.Additionally, we initialize the B matrices using a constant of all-1 vector, and we initialize C randomly using Glorot.As detailed in section 3.1 and the algorithmic outlines provided, the implementations of GRAPHSSM across all three variants can be stratified into two primary phases: a diffuse-and-mixing step, and a linear recurrence step. The diffuse-and-mixing stage facilitates straightforward parallelization through the employment of methods such as graph batching. The inherent linear characteristic of the recurrence operation permits the utilization of efficient computation strategies, notably the selective scan technique as introduced in [10] . This approach yields a FLOP complexity of O(V LDN ) per SSM layer with work-efficient parallelization, concurrently achieving IO efficiency. Furthermore, note that if we replace the adaptive time gap mechanism into a constant, i.e., we use ∆ l ≡ 1 L , 1 ≤ l ≤ L in line 8 of algorithm 1 and algorithm 2, the resulting linear system is time-invariant and we can use other computational accelerations like convolution [13, 7] and parallel scan [40] .Input: A sequence of graph (snapshots) G 1:L with each of size V .Node (hidden) feature inputs X 1:L ∈ R V ×L×D . A graph neural network GNN θ parameterized by θ. Three graph neural networks for selective state spaces GNN θ B , GNN θ C , GNN ∆ . A mixing mechanism MIX ϕ parameterized by ϕ. State-space parameters A ∈ R D×N .Output:Y 1:L ∈ R V ×L×D 1: # Approximate diffusion via GNN 2: for t = 1 to L do 3: Z l = GNN θ (X l , G l );4:H l = Z l if l = 1 else MIX(Z l , Z l-1 ); 5: end for 6: Initialize state U 0 = 0; # SISO state of shape V × D × N 7: for t = 1 to L do 8: ∆ l = softplus (GNN ∆ (X l , G l ) + b); 9: A = exp (einsum(∆ l , A, "V D, DN → V DN ")); 10: B = einsum(∆ l , GNN θ B (X l , G l ), "V D, V N → V DN "); 11: U l = U l-1 ⊙ A + einsum(B, H l , "V DN, V D → V DN "); 12: C = GNN θ C (X l , G l ); 13:Y l = einsum(U l , C, "V DN, DN → V D"); 14: end for; 15: return Y 1:L ;In this section, we discuss the limitations of the GRAPHSSM framework and propose a few future research directions that might be of interest.In this study, we focus on modeling discrete-time temporal graphs (DTTGs) through the lens of discretizing continuously evolving systems. The continuous-time viewpoint holds promise for encapsulating the modeling of continuous-time temporal graphs (CTTGs), a domain of growing importance in graph learning literature. However, the current GHiPPO framework has its limitations when extending to continuous-time setups. We provide a brief discussion as follows:DTDG, CTDG and GHiPPO Recall that in our formulation of the underlying graph process (2), the node features evolve continuously and the topological relations among nodes allow finite (countable) mutations. In DTTG representations, we do not directly observe the events, but we observe the entire graph at certain time spots resulting in a serious of snapshots. In this spirit, DTTGs have complete latitudinal information, but are lossy regarding longitudinal information. In CTTG representations, we have complete observations of events, but upon each event information, we do not observe the features of the rest of the nodes (that do not participate in those specific events). Therefore, CTTGs have complete longitudinal information, but are lossy regarding latitudinal information. In this regard, we may view DTDG and CTDG as two different lossy observation schemes of the underlying graph process in the GHiPPO abstraction.Handing CTTGs using SSM discretizations is challenging In section 3.2 of our paper (especially theorem 2), we established the discretization scheme upon an ideal, discrete observation (We observe the graph snapshot at each mutation events). We believe that this result might reasonably hints the gap between possible empirical approximations in either DTTG or CTTG scenarios: In DTTGs, we believe approximations using available snapshots are possible since from hindsight, the ideal representation is a convex combination of the snapshot representations at the mutation times. The approximation bias mostly comes from fewer snapshots, and we use mixing strategies to mitigate the biases. However, in CTTG scenarios, we miss the majority of information in each snapshot. Besides, consturcting snapshots from CTDGs is itself a very impractical method. Hence, we regard the modeling of CTDG to be beyond the scope of GraphSSM.The distinguishing algorithmic feature of GHIPPO compared to the conventional HIPPO framework lies in the piecewise nature of the dynamical system it generates. This characteristic leads to the challenge of dealing with unobserved dynamics, a factor that motivated the development of our MIX module. However, it's important to acknowledge that the mixing module serves as an approximation of the actual underlying dynamics, thus representing a limitation within the framework. This acknowledgment raises an intriguing question: might there exist alternative problem formulations capable of yielding a smoother dynamical system that mitigates the issue of discontinuities? One potential pathway could involve adopting smoother versions of the Laplacian or revising the approximation objective specified in (3) towards one that fosters a smooth solution. Such a solution would promote consistency in the dynamics across the complete temporal interval. Implementing these innovations would, however, necessitate the incorporation of more sophisticated technical assumptions and theoretical tools which we left for future explorations. Temporal continuity. As illustrated in figure 2 , our work has highlighted the problem of unobserved graph mutations in learning from discrete-time temporal graphs. The issue of unobserved graph mutations greatly hampers the temporal continuity of such graphs, presenting a significant challenge for learning if not properly addressed. To quantitatively measure the temporal continuity of a temporal graph, we calculate the average proximity between consecutive graph snapshots in the graph sequence. Specifically, we utilize Jaccard distance and Cosine similarity to measure the temporal continuity in terms of graph structure and node features, respectively:EQUATION)Datasets. We focus on the node classification task in discrete-time temporal graphs, which is a straightforward extension of static graphs. The experiments are conducted on six temporal graph benchmarks with different scales and time snapshots, including DBLP-3 [47] , Brain [47] , Reddit [47] , DBLP-10 [25] , arXiv [19] , and Tmall [25] . Dataset statistics are summarized in table 6 including the corresponding temporal continuity. The graph datasets are collected from real-world networks belonging to different domains. It should be noted that in the arXiv dataset, the time information is associated with the nodes rather than the edges. As a result, we split the snapshots of arXiv based on the occurrence of nodes. Each snapshot graph in the dataset shares the same attribute information but not the topology. Therefore, TC feature = 1.000 for arXiv in our experiments.Baselines. We compare GRAPHSSM with the following baselines: (i) static graph embedding methods: DeepWalk [36] , Node2Vec [9] ; (ii) temporal graph embedding methods: HTNE, M 2 DNE, and DynamicTriad [50] ; (iii) discrete-time temporal graph neural networks: MPNN [32] , STAR [47] , tNodeEmbed [39] , EvolveGCN [34] , SpikeNet [25] , and ROLAND [48] . For baselines that are originally designed for static graphs, we accumulate historical information (edges) in the graph snapshot sequence and represent the static graph structure at the last time point. All baselines are carefully tuned to achieve their best results based on the code officially provided by the authors.Implementation details. GRAPHSSM is built on the success of SSMs, where in this work we have derived variants of GRAPHSSM-S4, GRAPHSSM-S5, and GRAPHSSM-S6, under different SSM settings. Our experiments are mainly conducted on the S4 architecture. we employ feature mixing for DBLP-10 and representation mixing for other datasets. The graph convolution networks used to learn the graph structure are SAGE [16] for all datasets, except for arXiv, where GCN [21] is used. We implement our models as well as baselines with PyTorch [35] and PyTorch Geometric [5] , which are open-source software released under BSD-style8 and MIT9 license, respectively. All datasets used throughout experiments are publicly available. All experiments are conducted on an NVIDIA RTX 3090 Ti GPU with 24 GB memory. Code will be made available at https: //github.com/EdisonLeeeee/GraphSSM.Evaluation protocol. We adopt the conventional transductive learning setting, where the graph structure of all snapshots is visible during both training and inference stages. This is analogous to the standard node classification task, but with the additional incorporation of time information to facilitate the learning. For the DBLP-3, Brain, and Reddit datasets, we adopt the 81%/9%/10% train/validation/test splits as suggested in [47] . For the DBLP-10 and Tmall datasets, we follow the experimental settings of previous works [25] , where 80% of the nodes are randomly selected as the training set, and the remaining nodes are used as the test set. Note that stratified sampling is employed to ensure that the class distribution remains consistent across splits. For the arXiv dataset, we use the fixed public splits. We use Micro-F1 and Macro-F1 to evaluate the node classification performance. We report the average performance with standard deviation across 5 runs for each method.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors). • It should be clear whether the error bar is the standard deviation or the standard error of the mean. Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction).Technically, we require Z and X to reside within some appropriately defined Hilbert space. A comprehensive treatment will be provided in appendix C.To be more precise, at any time s ∈ [0, t], the integrand Z(s) ⊤ L(s)Z(s) inside the second term of (3) attains its minimum when Z(s) satisfies certain smoothness criterion which is determined via the choice of graph Laplacian. A more detailed explanation is deferred to appendix B.1.More concretely, diagonal SSMs are defined by diagonal A matrices whose diagonal elements lie on the complex plane with negative real parts[12], yet recent developments have found that complex state matrices are often not necessary[10]. For ease of representation, we only explore real state matrices in this paper.We let the balancing constant α be absorbed into the learnable parameters. Indeed, for GNNs that employ asymmetric aggregation[16], it is plausible to conceptualize the GNN as engaging a form of auto-balancing.In these algorithmic descriptions, we illustrate using the representation mixing mechanism. The case for feature mixing is similarly defined.https://github.com/pytorch/pytorch/blob/master/LICENSEhttps://github.com/pyg-team/pytorch_geometric/blob/master/LICENSE 


Decision: This paper presents a well-motivated and technically sound approach to applying state space models (SSMs) to temporal graph neural networks (TGNNs), supported by strong theoretical foundations and promising experimental results. However, it has some limitations, such as the lack of comparison with recent memory-based TGNNs, which are state-of-the-art, and the omission of key tasks like link prediction in the empirical evaluation. Additionally, the model's applicability is limited to discrete-time graphs, and potential issues with heterophilic graphs warrant further exploration. Despite these drawbacks, the paper offers a valuable and interesting research direction.

Reviews: review: paper_summary: paper_summary: This paper investigates SSM theory to temporal graphs by integrating structural information into the online approximation with laplacian regularization term.
strength_weakness: strength_weakness: 1. Since one advantage of the SSM based methods is its small parameter size, then what's the size of the model compared with other baseline methods? It's better to list the number of parameters to show the efficiency of the proposed method.

2. The transformer based methods have better performance as shown in Table 1, do you include any transformer-based methods for the experimental comparison?
questions: questions: 1. Since one advantage of the SSM based methods is its small parameter size, then what's the size of the model compared with other baseline methods?

2. The transformer based methods have better performance as shown in Table 1, do you include any transformer-based methods for the experimental comparison?

3. What is the time complexity of the GraphSSM?
limitations: limitations: The authors provide the limitation in appendix E.
paper_summary: paper_summary: The paper introduces GRAPHSSM, a novel state space model framework for temporal graphs. GRAPHSSM extends state space models (SSMs) by incorporating structural information through Laplacian regularization to handle the dynamic behaviors of temporal graphs. The framework aims to overcome limitations of recurrent neural networks (RNNs) and Transformers in modeling long-range dependencies and managing computational complexity. The authors propose the GHIPPO abstraction for memory compression and various mixing mechanisms to handle unobserved graph mutations. Extensive experiments on multiple temporal graph benchmarks demonstrate the effectiveness and efficiency of GRAPHSSM.
strength_weakness: strength_weakness: - There already exists an important family of temporal graph neural networks (TGNNs) that use memory mechanism (see [1] and [2] for examples) and resemble a state transition model. The paper seems to neglect/be unaware of this line of related work.

- There are parts of the method/extension that are not clear (see questions below for more details)

- There are a couple of serious limitations regarding the empirical study: 
1) the baselines for the paper are rather outdated. The paper did not compare with the more recent TGNN models (e.g. above-mentioned memory-based temporal graph neural network, which is the SOTA for TGNN) 
2) the performance improvement from the proposed baselines is still marginal
3) the conical task for temporal graph neural should be link prediction which is missing from the current study

[1] Rossi, Emanuele, et al. "Temporal graph networks for deep learning on dynamic graphs." arXiv preprint arXiv:2006.10637 (2020).

[2] Zhang, Yao, et al. "Tiger: Temporal interaction graph embedding with restarts." Proceedings of the ACM Web Conference 2023. 2023.
questions: questions: - Q1: I fail to understand the second part of Eq.(3).  Why does minimizing the second term of Eq.(3) amounts to increasing smoothness?  In addition, I don't think the interpretation/remark given in line 179-180 are correct. If $\alpha$ goes to infinity, the whole objective, Eq.(3), is dominated by the second term. There is a degenerate solution $Z(s) = \mathbf{0}$ (zero vector) that would minimize the objective. Which part of the objective/model prevent this degenerate solution?

- Q2: Theorem 1 is saying that GHIPPO is "nice" because its parameter update can be described by the given ODEs. Can you further explain why this property is nice? does the benefit comes from better efficiency or inference performance? and how?

- Q3: the discretization in the HIPPO paper was for computation purposes. In the case of temporal graph, where event and snapshot are inherently discrete. What is the (optimization) objective for discretization? i.e., what should be considered as a "good discretization" in this case?
limitations: limitations: N.A
paper_summary: paper_summary: The paper proposes an approach for processing discrete-time temporal graphs via an extension of state-space models (SSMs). In this direction, the main contribution of the paper is a generalization of the HIPPO framework for graph structured data (named GHIPPO), which defines how the state of different nodes should be updated through time to compress the evolution of nodes features (a solution capable of dealing with multiple changes that might happen between two snapshots is additionally introduced in the paper). In experimental evaluation, an SSM-based solution implemented with the proposed GHIPPO framework (named GraphSSM) outperforms a variety of approaches on both the moderately sized DBLP-3, Brain, Reddit and DBLP-10 datasets, and the large-scale arXiv and Tmall datasets.
strength_weakness: strength_weakness: Generally I do not have many weaknesses to highlight in the paper. From a technical perspective, I would want to highlight that: 

1) the proposed approach might face some issues with heterophilic graphs (as the smoothness in the node state imposed in equation (3) is not necessarily a good regularizer in that setting - see point 5 below) 

2) the proposed methodology is limited to discrete-time temporal graphs, and further work would need to be done to extend the model to a continuous-time setting (this is also highlighted by the authors themselves in the paper). 

From a different perspective, a weakness of the paper is probably its exposition, which in some parts is a little confusing and could be improved. In this direction:

1) lines 67, I believe there is a typo and the authors should replace “features the dynamics” with “the features dynamics”
2) line 74 to 77 are quite unclear and should probably be refactored
3) line 105, “recurrence model” should probably be replaced with “recurrent model”
4) Table 1, while it is true that classic transformers have a quadratic space and time complexity, there are also efficient transformers (e.g. Performs) that show linear space and time complexity. I believe these are not mentioned in the paper and probably they should, in order to provide to the reader a complete picture of where we stand with transformer architectures.
5) I believe the smoothness in node state imposed in equation 3 with the dirichlet energy regularizer makes sense only on homophilic graphs, and might not be a meaningful prior to use when dealing with heterophily. This is not pointed out in the paper and probably it should be highlighted as a possible limitation.
6) Lines 180-181, the authors mention that what minimizes equation 3 is a piecewise vector. I believe what they mean is that the solution is piecewise over time, I think this should probably be clarified in the paper as it might be a source of confusion
7) line 184-185, how the nodes' memories are parameterized as polynomials is very unclear here. How are those polynomials defined? What is the input to such polynomials? It might be that a reader familiar with the HIPPO framework would be aware of such details, but the lack of an explanation in the text makes it hard to draw a complete picture here.
8) lines 335-336, GraphSSM with S5 doesn’t always achieve poor performance (it is indeed the best performing model on DBLP-3 in table 4). This should be clarified in the text.
9) line 339, I believe there is a typo and the “and” at the beginning of the sentence should be fixed.
questions: questions: One thing that I think would be helpful to improve the strength of the paper is to run some experiments where the proposed SSM module is swapped with a RNN or transformer architecture. As we can see from the ablation study of table 5, it appears that the mixing strategy one might decide to use is relevant to outperform prior art. As a result of this, I wonder what would be the performance that a comparable solution would have with the very same architecture used by the authors but without using any SSM approach. This will clearly provide an indication on whether using SSMs is beneficial for achieving the desired performance (which is the main claim of the paper) or whether the results we see in Table 2 and 3 are only due to the mixing strategy.
limitations: limitations: Yes

Comments: title: Thank you for raising the score
Gentle Reminder
Response to Reviewer PUCW (Part 2/2)
Response to Reviewer PUCW (Part 1/2)
Reply to authors' rebuttal
comment: We are very happy to hear that you have raised the score. Your recognition of our efforts in the paper is greatly appreciated! We thank you again for the invaluable feedback and insights.
Thank you for the diligent response. I have raise my score accordingly.
Dear Reviewer PUCW,

We sincerely appreciate your insightful review and feedback comments.  As the author-reviewer discussion deadline (Aug 13) is approaching, we would like to check if you have any other remaining concerns about our paper. 

**We have faithfully responded to ALL your comments. If our responses have adequately addressed your concerns, we kindly hope that you can consider increasing the score.** We understand that you have a demanding schedule, and we appreciate the time and effort you dedicate to reviewing our paper.

Kind regards,

Authors
The comparison CTDG baselines include memory-based methods TGN and TIGER (as you mentioned in the initial review), as well as attention-based methods TGAT and DyGFormer. Results for the node classification tasks are presented below:

||**DBLP-3**||**Brain**||**Reddit**||**DBLP-10**||
|---------|:-----------:|:-----------:|:-----------:|:-----------:|:-----------:|:-----------:|:-----------:|:-----------:|
||Micro-F1|Macro-F1|Micro-F1|Macro-F1|Micro-F1|Macro-F1|Micro-F1|Macro-F1|
|TGN|80.43|80.08|89.50|89.59|40.66|40.50|70.01|68.83|
|TIGER|81.53|81.74|90.49|90.33|41.17|40.98|72.48|71.14|
|TGAT|82.59|82.14|90.98|90.25|42.79|42.40|73.15|72.02|
|DyGFormer|84.77|84.51|92.01|91.32|45.24|44.51|72.40|71.86|
|GraphSSM|**86.29±1.0**|**85.78±0.9**|**93.80±0.3**|**94.47±0.6**|**49.21±0.5**|**49.05±0.7**|**76.80±0.3**|**76.00±0.4**|

Due to time and space limitation, we are currently presenting the experimental results of the node classification task on the four datasets. Further experiments are being conducted and can be presented in a few days if it is necessary. As observed, CTDG methods demonstrate relatively poor performance on the four datasets, which is inline with our above claims. Indeed, many CTDG methods require fine-grained time information to capture the evolution patterns between edges, which is impractical for DTDGs with coarsened snapshot information. 

**Response to Q2:** Firstly, we acknowledge that the learning setting in ROLAND differs from that of other works [1, 2]. The main differences arise from diverse contexts of dynamic graphs, as detailed below:

+ ROLAND focus on the **discrete-time dynamic graphs (DTDG)**, which is also the main focus of our work. In this setting, our goal is to predict edges in snapshot 𝑡 + 1 based on the snapshots accumulated up to time 𝑡. This is the general setting in discrete-time temporal graphs and is widely used in literature including several representative works, e.g., EvolveGCN[4]. 
+ For TGB[1] and DGB[2], both benchmarks focus on the **continuous-time dynamic graphs (CTDG)**, whose goal is to predict the occurenc of a link between two given nodes at a specific time point. 
+ While we agree that "Continuous vs. discrete are often modeling choices most of the time", there are differences in terms of the associated time information in the graphs for modeling. In the context of CTDG, link prediction is **fine-grained**, with each link being associated with an exact timestamp. This allows for predicting evolving edge streams at specific time points. However, in DTDG, link prediction is **coarse-grained**. Each edge in a snapshot has relative time information, and we can only predict whether a link will occur in the next time span, such as in the next week or month. In this regard, we were unable to perform link prediction under the setting of [1,2].

Secondly, while ROLAND studied and proposed a framework that is related to continual/incremental learning, the experimental setup for the link prediction task exactly follows the classic DTDG settings in previous works[4,5,6,7]. This is detailed in the ROLAND paper ("Task" section in Section 4.1). 

Finally, to avoid any potential confusion, we detail our experimental setting for the link prediction task:

+ **Train-val-test splits:** For each dataset, we use $G_1, G_2, · · · , G_{T −1}$ as the training data, $G_{T}$ as test graph. Also, we use $G_{T-1}$ as the validation graph during training for hyperparameter tuning. Each model leverages the T-1 visible graphs to prediction the occurence of edges in $G_T$.
+ **Positive and negative edges:** For the link prediction problem, the positive test set consists of the edges that appear in $G_T$ and do not present in $G\_{T-1}$, while the negative test set consists of the edges that do not appear in $G\_{T-1}$ and $G_T$.
+ **Evaluation metrics:** AUC and AP, two standard metrics in link prediction tasks.
---

If any concern still remains that might prohibit a positive recommendation of this work, we would appreciate if you could let us know.


### Reference

[3] Kazemi, Seyed Mehran, et al. "Representation learning for dynamic graphs: A survey." Journal of Machine Learning Research 21.70 (2020): 1-73.

[4] Pareja, Aldo, et al. "Evolvegcn: Evolving graph convolutional networks for dynamic graphs." Proceedings of the AAAI conference on artificial intelligence. Vol. 34. No. 04. 2020.

[5] You, Jiaxuan, Tianyu Du, and Jure Leskovec. "ROLAND: graph learning framework for dynamic graphs." *Proceedings of the 28th ACM SIGKDD conference on knowledge discovery and data mining*. 2022.

[6] Zhang, Kaike, et al. "Dyted: Disentangled representation learning for discrete-time dynamic graph." *Proceedings of the 29th ACM SIGKDD Conference on Knowledge Discovery and Data Mining*. 2023.

[7] Sankar, Aravind, et al. "Dysat: Deep neural representation learning on dynamic graphs via self-attention networks." Proceedings of the 13th international conference on web search and data mining. 2020.
Thanks for acknowledging our efforts in the rebuttal. We greatly appreciate your prompt response during the rebuttal period. We will now address your additional concerns and questions below.

**Response to Q1:** Thank you for the question regarding TGNNs. In our initial rebuttal, due to limited characters we only tried to answer your question in a somewhat conventional way by pointing out modeling on discrete-time dynamic graphs (DTDG) and continuous-time dynamic graphs (CTDG) are usually considered as different tasks[3]. We agree with you that there are more profound connections between the modeling of DTDG, CTDG and the GHiPPO abstraction we proposed. Our perspectives are three-fold: 

+ (i) We will show that DTDG and CTDG are two different lossy observation schemes of the underlying graph process in the GHiPPO abstraction. 
+ (ii) We discuss why the current SSM framework based on discretizations are not yet readily suitable for handling CTDGs. 
+ (iii) We return to CTDG models as you have mentioned, and present some empirical results suggesting that they might indeed be not an ideal solution to DTDG problems.

### (i) DTDG, CTDG and GHiPPO

Recall that in our formulation of the underlying graph process, the node features evolve continuously and the topological relations among nodes allow finite (countable) mutations. Let's rewrite the process in equation (2) of the paper:
$G(0) \overset{\mathcal{E}_1}{\longrightarrow} G(t_1) \overset{\mathcal{E}_2}{\longrightarrow} G(t_2) \longrightarrow \cdots \longrightarrow  G(t\_{M-1}) \overset{\mathcal{E}_M}{\longrightarrow}  G(t_M) = G(T).$

Now we take a closer look at DTDG and CTDG representations of the above process:

- **DTDG**: In DTDG representations, we do not directly observe the events, but we observe the entire graph at certain time spots resulting in a serious of snapshots. In this spirit, **DTDG has complete latitudinal information, but is lossy regarding longitudinal information**.
- **CTDG**: In CTDG representations, we have complete observations of events, but upon each event information, we do not observe the features of the rest of the nodes (that do not participate in those specific events). Therefore, **CTDG has complete longitudinal information, but is lossy regarding latitudinal information**.

Next we proceed to why we choose to only study the former one (DTDG).

### (ii) Handing CTDGs using SSM discretizations is challenging

In section 3.2 of our paper (especially theorem 2), we established the discretization scheme upon an ideal, discrete observation. (We observe the graph snapshot at each mutation events) We believe that this one reasonably hints tha gap between possible empirical approximations in either DTDG or CTDG scenarios. In DTDGs, we believe approximations using available snapshots are possible since from hindsight, the ideal representation is a convex combination of the snapshot representations at the mutation times. The approximation bias mostly comes from fewer snapshots, and we use mixing strategies to mitigate the biases. 

However, in CTDG scenarios, we miss the majority of information in each snapshot. Up till now we have not yet figured out any practical and sound solutions to this issue. Besides, consturcting snapshots from CTDGs is itself a very impractical method. It is possible that the GHiPPO abstraction might not be a very good fit for a principled study of CTDGs (might be too stringent, in our opinion). Hence, we regard the modeling of CTDG to be beyond the scope of GraphSSM.

### (iii) CTDG models and their performance on DTDG graphs

As you have pointed out, many CTDG models are based on a **Message Passing combined with Recurrent State Update (MPRSU)** scheme, which is closely related to GraphSSM (and indeed, related to many other DTDG methods like ROLAND as well). Recall that, SSMs are also a special case of linear RNNs from the computational viewpoint. However, **SSMs utilize a finer-grained definition of memory (via deriving them from an online approximation problem) that yields more stable and performant RNN variants as compared to GRU or LSTM.** This is indeed what we tried to establish in our paper---We equip the node-wise memory in a principled way that is formalized as the optimal approximation coefficients of the Laplacian-regularized online approximation problem. We have shown that this formulation leads to better DTDG modeling frameworks, and in our ablation study we also show the efficacy of using SSM schemes rather than ordinary RNNs or even softmax transformers. So far as we have noticed, most TGNN methods have no mechanistic interpretations of their memory mechanisms. Besides, as we have discussed earlier, it is unclear if SSM schemes are readily applicable to CTDG frameworks. 

Finally, we have gained some additional empirical insights from additional experiments. We conducted evaluations of CTDG methods on the DBLP-3, Brain, Reddit, and DBLP-10 datasets.
Thank you for your response. We greatly appreciate your thoughtful feedback and the time you dedicated to reviewing our paper!
I would like to thank the authors for the detailed response. I don't have further question.
I would like first to thank the authors for their diligent response. However, I still found my main concerns not fully addressed after reading the response.

1. I don't think "they focus on modelling continuous graph" is a sufficient reason for not comparing them with the proposed method. Continuous VS discrete are often modelling choices most of the time, and memory-based TGNN have been applied to the dataset used in the paper [1,2]. I also think the modelling choice of this paper is also assuming the underlying process is continuous?

The most appealing reason why I think a detailed comparison is necessary between memory-based TGNN and the proposed method is their memory mechanism. Memory-based TGNN captures dynamic/temporal information of the graph by incrementally processing information and distilling it into a so-called memory vector, which in a sense can be viewed as the state of the graph/node. Therefore, I think the mechanisms between memory-based TGNN and the method proposed are very similar. I am very keen to understand the fundamental difference and why we need this potentially new family of models for dynamic graphs.

2. I am very confused with the choice of experimental setting for the link prediction task. Based on my understanding, ROLAND actually focused on the efficiency of TGNN and studied a learning setting that is different from the typical ones (e.g., [1,2]). They studied and proposed a framework that is related to continual/incremental learning. 


[1] Huang, Shenyang, et al. "Temporal graph benchmark for machine learning on temporal graphs." Advances in Neural Information Processing Systems 36 (2024).

[2] Poursafaei, Farimah, et al. "Towards better evaluation for dynamic link prediction." Advances in Neural Information Processing Systems 35 (2022): 3
Thank you for your response and for raising the score. We are pleased to address the concerns you have raised.
I would like to thank the authors for their response. I don't have further question and I'm happy to raise my score to a 6