Adversarial Robustness in Graph Neural Networks: A Hamiltonian Approach
Graph neural networks (GNNs) are vulnerable to adversarial perturbations,
including those that affect both node features and graph topology. This paper
investigates GNNs derived from diverse neural flows, concentrating on their
connection to various stability notions such as BIBO stability, Lyapunov
stability, structural stability, and conservative stability. We argue that
Lyapunov stability, despite its common use, does not necessarily ensure
adversarial robustness. Inspired by physics principles, we advocate for the use
of conservative Hamiltonian neural flows to construct GNNs that are robust to
adversarial attacks. The adversarial robustness of different neural flow GNNs
is empirically compared on several benchmark datasets under a variety of
adversarial attacks. Extensive numerical experiments demonstrate that GNNs
leveraging conservative Hamiltonian flows with Lyapunov stability substantially
improve robustness against adversarial perturbations. The implementation code
of experiments is available at
https://github.com/zknus/NeurIPS-2023-HANG-Robustness.
Kai Zhao, Qiyu Kang, Yang Song, Rui She, Sijie Wang, Wee Peng Tay
# Adversarial Robustness in Graph Neural Networks:

###### Abstract

Graph neural networks (GNNs) are vulnerable to adversarial perturbations, including those that affect both node features and graph topology. This paper investigates GNNs derived from diverse neural flows, concentrating on their connection to various stability notions such as BIBO stability, Lyapunov stability, structural stability, and conservative stability. We argue that Lyapunov stability, despite its common use, does not necessarily ensure adversarial robustness. Inspired by physics principles, we advocate for the use of conservative Hamiltonian neural flows to construct GNNs that are robust to adversarial attacks. The adversarial robustness of different neural flow GNNs is empirically compared on several benchmark datasets under a variety of adversarial attacks. Extensive numerical experiments demonstrate that GNNs leveraging conservative Hamiltonian flows with Lyapunov stability substantially improve robustness against adversarial perturbations. The implementation code of experiments is available at [https://github.com/zknus/NeurIPS-2023-HANG-Robustness](https://github.com/zknus/NeurIPS-2023-HANG-Robustness).

## 1 Introduction

Graph neural networks (GNNs) [1; 2; 3; 4; 5; 6; 7; 8] have achieved great success in inference tasks involving graph-structured data, including applications from social media networks, molecular chemistry, and mobility networks. However, GNNs are known to be vulnerable to adversarial attacks [9]. To fool a trained GNN, adversaries can either add new nodes to the graph during the inference phase or remove/add edges from/to the graph. The former is called an injection attack [10; 11; 12; 13], and the latter is called a modification attack [14; 15; 16]. In some works [9; 17], node feature perturbations are also considered to enable stronger modification attacks.

Neural ordinary differential equation (ODE) networks [18] have recently gained popularity due to their inherent robustness [19; 20; 21; 22; 23; 24]. Neural ODEs can be considered as a continuous analog of ResNet [25]. Many neural ODE networks have since been proposed, including but not limited to [19; 20; 26; 27; 28]. Using neural ODEs, we can constrain the input and output of a neural network to follow certain physics laws. Injecting physics constraints to black-box neural networks improves neural networks' explainability. More recently, neural ODEs have also been successfully applied to GNNs by modeling
the way nodes exchange information given the adjacency structure of the underlying graph. We call them _graph neural flows_ which enables the interpretation of GNNs as evolutionary dynamical systems. These system equations can be learned by instantiating them using neural ODEs [29, 30, 31, 32, 33, 34]. For instance, [29, 30] models the message-passing process, i.e., feature exchanges between nodes, as the heat diffusion, while [31, 32] model the message passing process as the Beltrami diffusion. The reference [35] models the graph nodes as coupled oscillators with a coupled oscillating ODE guiding the message-passing process.

Although adversarial robustness of GNNs has been investigated in various works, including [9, 13, 17], robustness study of graph neural flows is still in its fancy. To the best of our knowledge, only the recent paper [32] has started to formulate theoretical insights into why graph neural diffusion is generally more robust against topology perturbation than conventional GNNs. The concept of Lyapunov stability was used in [20, 21]. However, there are many different notions of stability in the dynamic system literature [36, 37]. In this paper, we focus on the study of different notions of stability for graph neural flows and investigate which notion is most strongly connected to adversarial robustness. We impose an energy conservation constraint on graph neural flows that lead to a Hamiltonian graph neural flow. We find that energy-conservative graph Hamiltonian flows endowed with Lyapunov stability improve the robustness the most as compared to other existing stable graph neural flows.

**Main contributions.** This research is centered on examining various stability notions within the realm of graph neural flows, especially as they relate to adversarial robustness. Our main contributions are summarized as follows:

1. We revisit the definitions of stability from the perspective of dynamical systems as applicable to graph neural flows. We argue that vanilla Lyapunov stability does not necessarily confer adversarial robustness and provide a rationale for this observation.
2. We propose Hamiltonian-inspired graph neural ODEs, noted for their energy-conservative nature. We perform comprehensive numerical experiments to verify their performance on standard benchmark datasets. Crucially, our results demonstrate that Hamiltonian flow GNNs present enhanced robustness against various adversarial perturbations. Moreover, it is found that the effectiveness of Lyapunov stability becomes pronounced when layered on top of Hamiltonian flow GNNs, thereby fortifying their adversarial robustness.

The rest of this paper is organized as follows. We introduce various stability notions from a dynamical system perspective in Section 2. A review of existing graph neural flows is presented in Section 3 with links to the stability notions defined in Section 2. We present a new type of graph neural flow inspired by the Hamiltonian system with energy conservation in Section 4. Two different variants of this model are proposed in Section 5. Section 6 details our extensive experimental outcomes. The supplementary section provides an overview of related studies, an exhaustive outline of the algorithm, further insights into model robustness, supplementary experimental data, and the proofs for the theoretical propositions made throughout the paper.

## 2 Stability in Dynamical Systems

It is well known that a small perturbation at the input of an unstable dynamical system will result in a large distortion in the system's output. In this section, we first introduce various types of stability in dynamical physical systems and then relate them to graph neural flows. We consider the evolution of a dynamical system that is described as the following autonomous nonlinear differential equation:

\[\frac{\mathrm{d}\mathbf{z}(t)}{\mathrm{d}t}=f_{\boldsymbol{\theta}}(\mathbf{z }(t)), \tag{1}\]

where \(f_{\boldsymbol{\theta}}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}\) denotes the system dynamics, which may be non-linear in general, \(\boldsymbol{\theta}\) denotes the system parameters, and \(\mathbf{z}:[0,\infty)\rightarrow\mathbb{R}^{n}\) represents the \(n\)-dimensional system state.

We first introduce the stability notions from a dynamical systems perspective that is related to the GNN robustness against _node feature perturbation_.

**Definition 1** (Bibo stability).: _The system is called BIBO (bounded input bounded output) stable if for any bounded input, there exists a constant \(M\) s.t. the output \(\|\mathbf{z}(t)\|<M,\forall\,t\geq 0\)._

Suppose \(f\) has an equilibrium at \(\mathbf{z}_{e}\) so that \(f_{\boldsymbol{\theta}}\left(\mathbf{z}_{e}\right)=0\). We can define the stability notion for \(\mathbf{z}_{e}\).
**Definition 2** (Lyapunov stability and asymptotically stable [38]).: _The equilibrium \(\mathbf{z}_{e}\) is Lyapunov stable if for every \(\epsilon>0\), there exists a \(\delta>0\) such that, if \(\|\mathbf{z}(0)-\mathbf{z}_{e}\|<\delta\), then for every \(t\geq 0\) we have \(\|\mathbf{z}(t)-\mathbf{z}_{e}\|<\epsilon\). Furthermore, the equilibrium point \(\mathbf{z}_{e}\) is said to be asymptotically stable if it is Lyapunov stable and there exists a \(\delta^{\prime}>0\) such that if \(\|\mathbf{z}(0)-\mathbf{z}_{e}\|<\delta^{\prime}\), then \(\lim_{t\to\infty}\|\mathbf{z}(t)-\mathbf{z}_{e}\|=0\)._

**Remark 1**.: _Lyapunov stability indicates that the solutions whose initial points are near an equilibrium point \(\mathbf{z}_{e}\) stay near \(\mathbf{z}_{e}\) forever. For the special linear time-invariant system \(\,\mathrm{d}\mathbf{z}(t)/\,\mathrm{d}t=\mathbf{A}\mathbf{z}(t)\) with a constant matrix \(\mathbf{A}\), it is Lyapunov stable if and only if all eigenvalues of \(\mathbf{A}\) have non-positive real parts and those with zero real parts are the simple roots of the minimal polynomial of \(\mathbf{A}\)[39, 40]. Asymptotically stable means that not only do trajectories stay near \(\mathbf{z}_{e}\) for all time (Lyapunov stability), but trajectories also converge to \(\mathbf{z}_{e}\) as time goes to infinity (asymptotic stability)._

We next introduce the concept of structural stability from dynamical systems theory, which is related to the robustness of GNNs against _graph topological perturbation_. It describes the sensitivity of the qualitative features of a solution to changes in parameters \(\boldsymbol{\theta}\). The definition of structural stability requires the introduction of a topology on the space of \(\mathbf{z}\) in (1), which we do not however present here rigorously due to space constraints and not to distract the reader with too much mathematical details. Instead, we provide a qualitative description of structural stability to elucidate how it can indicate the robustness of a graph neural flow against topology perturbations.

**Definition 3** (Structural stability).: _Unlike Lyapunov stability, which considers perturbations of initial conditions for a fixed \(f_{\boldsymbol{\theta}}\), structural stability deals with perturbations of the dynamic function \(f_{\boldsymbol{\theta}}\) by perturbing the parameter \(\boldsymbol{\theta}\). The qualitative behavior of the solution is unaffected by small perturbations of \(f_{\boldsymbol{\theta}}\) in the sense that there is a homeomorphism that globally maps the original solution to the solution under perturbation._

**Remark 2**.: _In the graph neural flows to be detailed in Section 3 and Section 4, the parameter \(\boldsymbol{\theta}\) includes the graph topology (i.e., the adjacency matrix) and learnable neural network weights. Unlike adversarial attacks on other deep learning neural networks where the attacker targets only the input \(\mathbf{z}\), it is worth noting that adversaries for GNNs can also attack the graph topology, which forms part of \(\boldsymbol{\theta}\). If there are different Lyapunov stable equilibrium points, one for each class of nodes, one intuitive example of breaking structural stability in graph neural flows is by perturbing the graph topology in such a way that there are strictly fewer equilibrium points than the number of classes._

In this study, we will propose GNNs drawing inspiration from Hamiltonian mechanics. In a Hamiltonian system, \(\mathbf{z}=(q,p)\in\mathbb{R}^{2n}\) refers to the generalized coordinates, with \(q\) and \(p\) corresponding to the generalized position and momentum, respectively. The dynamical system is characterized by the following nonlinear differential equation:

\[\frac{\mathrm{d}\mathbf{z}(t)}{\mathrm{d}t}=J\nabla H(\mathbf{z}(t)), \tag{2}\]

where \(\nabla H(\mathbf{z})\) is the gradient of a scalar function \(H\) at \(\mathbf{z}\) and \(J=\left(\begin{array}{cc}0&I\\ -I&0\end{array}\right)\) is the \(2n\times 2n\) skew-symmetric matrix with \(I\) being the \(n\times n\) identity matrix.

We now turn our attention to the notion of conservative stability in dynamical systems. It is worth noting that a general dynamical system, as characterized in (1), might not consistently resonate with traditional perspectives on energy and conservation, especially when compared to physics-inspired neural networks, like (2).

**Definition 4** (Conservative stability).: _In a dynamical system inspired by physical principles, such as (2), a conserved quantity might be present. This quantity, which frequently embodies the notion of the system's energy, remains invariant along the system's evolution trajectory \(\mathbf{z}(t)\)._

Our focus in this work is on graph neural flows that can be described by either (1) or (2). Chamberlain et al. [29] postulate that many GNN architectures such as GAT can be construed as discrete versions of (1) via different choices of the function \(f_{\boldsymbol{\theta}}\) and discretization schemes. Therefore, the stability definitions provided above can offer additional insights into many popular GNNs. Most existing graph neural flows only scrutinize the BIBO/Lyapunov stability of their system. For instance, GRAND [29] proposes BIBO/Lyapunov stability against node feature perturbation over \(\mathbf{z}\). However, the more fundamental structural stability in graph neural flows, which is related to the robustness against graph topological changes, remains largely unexplored. Some models, such as GraphCON [29, 35], exhibit conservative stability under certain conditions. We direct the reader to Section 3 and Table 1 for a comprehensive discussion of the stability properties of each model.
## 3 Existing Graph Neural Flows and Stability

Consider an undirected, weighted graph \(\mathcal{G}=(\mathcal{V},\mathcal{E})\) where \(\mathcal{V}\) is a finite set of vertices and \(\mathcal{E}\subset\mathcal{V}\times\mathcal{V}\) denotes the set of edges. The adjacency matrix of the graph is denoted as \((\mathbf{W}[u,v])=\mathbf{W}([v,u])\) for all \([u,v]\in\mathcal{E}\). Let \(\mathbf{X}(t)\in\mathbb{R}^{|\mathcal{V}|\times r}\) represent the features associated with the vertices at time \(t\). In this section, we introduce several graph neural flows on \(\mathcal{G}\), categorizing them according to the stability concepts outlined in Section 2.

**GRAND:** Inspired by the heat diffusion equation, GRAND [29] employs the following dynamical system:

\[\frac{\mathrm{d}\mathbf{X}(t)}{\mathrm{d}t}=\overline{\mathbf{A}}_{G}( \mathbf{X}(t))\mathbf{X}(t)\coloneqq(\mathbf{A}_{G}(\mathbf{X}(t))-\alpha \mathbf{I})\mathbf{X}(t), \tag{3}\]

with the initial condition \(\mathbf{X}(0)\). Within this model, \(\mathbf{A}_{G}(\mathbf{X}(t))\) is either a time-invariant static matrix, represented as GRAND-l, or a trainable time-variant attention matrix \((a_{G}\left(\mathbf{x}_{i}(t),\mathbf{x}_{j}(t)\right))\), labeled as GRAND-nl, reflecting the graph's evolutionary features. The function \(a_{G}(\cdot)\) calculates similarity for pairs of vertices, and \(\mathbf{I}\) is an identity matrix with dimensions that fit the context. In [29], \(\alpha\) is set to be 1. Let \(\mathbf{D}\) be the diagonal node degree matrix where \(\mathbf{D}[u,u]=\sum_{v\in\mathcal{V}}\mathbf{W}[u,v]\).

**Theorem 1**.: _We can prove the following stability:_

1. _For GRAND-nl, if the attention matrix_ \(\mathbf{A}_{G}(\mathbf{X}(t))\) _is set as a doubly stochastic attention_ _[_41_]__, we have BIBO stability and Lyapunov stability for any_ \(\alpha\geq 1\)_. When_ \(\alpha>1\)_, it reaches global asymptotic stability under any perturbation._
2. _Within the GRAND-l setting, if_ \(\mathbf{A}_{G}\) _is set as a constant column- or row-stochastic matrix, such as the normalized adjacency matrices_ \(\mathbf{WD}^{-1}\) _or_ \(\mathbf{D}^{-1}\mathbf{W}\)_, global asymptotic stability is achieved for_ \(\alpha>1\) _under any perturbation. If the graph is additionally assumed to be strongly connected_ _[_42_]__[_Sec.6.3_]__, BIBO and Lyapunov stability are realized for_ \(\alpha=1\)_._
3. _Furthermore, when_ \(\mathbf{A}_{G}\) _is specifically a constant column-stochastic matrix like_ \(\mathbf{WD}^{-1}\) _and_ \(\alpha=1\)_, GRAND conserves a quantity that can be interpreted as energy. Furthermore, in this setting, asymptotic stability is attained when the graph is aperiodic and strongly connected and the perturbations on_ \(\mathbf{X}(0)\) _ensure unaltered column summations._

**BLEND:** _In comparison to GRAND, BLEND [31] introduces the use of positional encodings. Following a similar line of reasoning to that used for GRAND, BLEND also exhibits BIBO/Lyapunov stability as stated in Theorem 1. Moreover, it is noteworthy that if positional features_ \(\mathbf{U}(t)\) _are eliminated, for instance by setting them as a constant, BLEND simplifies to the GRAND model._

**GraphCON:** Inspired by oscillator dynamical systems, GraphCON is a graph neural flow proposed in [35] and defined as

\[\left\{\begin{array}{c}\frac{\mathrm{d}\mathbf{Y}(t)}{\mathrm{d}t}=\sigma( \mathbf{F}_{\theta}(\mathbf{X}(t),t))-\gamma\mathbf{X}(t)-\alpha\mathbf{Y}(t), \\ \frac{\mathrm{d}\mathbf{X}(t)}{\mathrm{d}t}=\mathbf{Y}(t),\end{array}\right. \tag{4}\]

where \(\mathbf{F}_{\theta}(\cdot)\) is a learnable \(1\)-neighborhood coupling function, \(\sigma\) denotes an activation function, and \(\gamma\) and \(\alpha\) are tunable parameters.

_As described in [35, Proposition 3.1], under specific settings where \(\sigma\) is the identity function and \(\mathbf{F}_{\theta}(\mathbf{X}(t),t)=\mathbf{AX}(t)\) with \(\mathbf{A}\) being a constant matrix, GraphCON conserves Dirichlet energy_ (11)_, thereby demonstrating conservative stability._

**GraphBel:** Generalizing the Beltrami flow, mean curvature flow and heat flow, a stable graph neural flow [32] is designed as

\[\frac{\mathrm{d}\mathbf{X}(t)}{\mathrm{d}t}=(\mathbf{A}_{\mathbf{S}}(\mathbf{ X}(t))\odot\mathbf{B}_{\mathbf{S}}(\mathbf{X}(t))-\Psi(\mathbf{X}(t)))\mathbf{X}(t), \tag{5}\]

where \(\odot\) is the element-wise multiplication. \(\mathbf{A}_{\mathbf{S}}(\cdot)\) and \(\mathbf{B}_{\mathbf{S}}(\cdot)\) are learnable attention function and normalized vector map, respectively. \(\mathbf{\Psi}(\mathbf{X}(t))\) is a diagonal matrix in which \(\Psi(\mathbf{x}_{i},\mathbf{x}_{i})=\sum_{\mathbf{x}_{j}}(\mathbf{A}\odot \mathbf{B})(\mathbf{x}_{i},\mathbf{x}_{j})\).

_Analogous to BLEND, under certain conditions with \(\Psi(\mathbf{X}(t))=\mathbf{B}_{\mathbf{S}}(\mathbf{X}(t))=\mathbf{I}\), GraphBel simplifies to the GRAND model. Consequently, it exhibits BIBO/Lyapunov stability in certain scenarios._

The incorporation of ODEs via graph neural flows may enhance the stability of graph feature representations. A summarized relationship between model stability and these graph neural flows can be found in Table 1.
### Lyapunov Stability vs. Node Classification Robustness:

At first glance, Lyapunov stability has a strong correlation with node classification robustness against feature perturbations. However, before diving into experimental evidence, we point out an important conclusion: **Lyapunov stability _by itself_ does not necessarily imply adversarial robustness.** Consider a scenario where a graph neural flow has only one equilibrium point \(\mathbf{z}_{e}\), while the node features are derived from more than one class. In a Lyapunov asymptotically stable graph neural flow, such as GRAND (as shown in Theorem 1), all node features across different classes would inevitably converge to a single point \(\mathbf{z}_{e}\) due to global contraction. We note that this is why the model in [20] requires a diversity-promoting layer to ensure that different classes converge to different Lyapunov-stable equilibrium points.

**Example 1**.: _We provide an example to demonstrate our claim. Consider the following Lyapunov stable ODE_

\[\dot{\mathbf{x}}(t)=\begin{pmatrix}-1&0\\ 0&-5\end{pmatrix}\mathbf{x}(t) \tag{6}\]

_with initial condition \(\mathbf{x}(0)=\left[x_{1}(0),x_{2}(0)\right]^{\intercal}\). The solution to this ODE is given by \(\mathbf{x}(t)=x_{1}(0)e^{-t}[1,0]^{\intercal}+x_{2}(0)e^{-5t}[0,1]^{\intercal}\). For all initial points in \(\mathbb{R}^{2}\), we have \(\mathbf{x}(t)\rightarrow\mathbf{0}\) as \(t\rightarrow\infty\). Furthermore, as \(t\rightarrow\infty\), the trajectory \(\mathbf{x}(t)\) for any initial point is approximately parallel to the \(x\)-axis. We draw the phase plane in Fig. 0(a)._

_Assume that the points on the upper half y-axis belongs to class 1 while we have a linear classifier that seperates class 1 and class 2 as shown in Fig. 0(a). We observe that for the initial point \(A\) belonging to class 1, the solution from a small perturbed initial point \(A+\epsilon\) is misclassified as class 2 for a large enough \(t\) for any linear classifier. We see from this example that Lyapunov stability iteself does not imply adversarial robustness in graph neural flow models._

_This example indicates that Lyapunov stability does not guarantee node classification robustness. Additionally, for a system exhibiting global contraction to a single equilibrium point, structural stability may also be ensured. For instance, in the case of GRAND, even if the edges are perturbed, the system maintains the same number of equilibrium points with global contraction._ We conclude that even an amalgamation of both Lyapunov stability and structural stability may not help the graph's adversarial robustness for node classification.

In the example shown in Fig. 1, we observe that in the case of GRAND when \(\alpha>1\), the node features from different classes tend to become closer to each other as time progresses. This phenomenon can potentially create more vulnerability to adversarial attacks.

## 4 Hamiltonian-Inspired Graph Neural Flow

Drawing inspiration from the principles of Hamiltonian classical mechanics, we introduce a novel graph neural flow paradigm, namely HamiltoniAN Graph diffusion (HANG). In Hamiltonian mechanics, the notation \((q,p)\) is traditionally used to represent the properties of nodes (position and

Figure 1: (a): We plot the vector field and the system solution trajectories of Example 1. (b) and (c): In GRAND, the node features’ energy tends to converge towards each other. In HANG, we observe the node features’ energy remains relatively stable over time. Two nodes are from different classes.
momentum). Hence, in this section, we adopt this notation instead of \(\mathbf{X}\) (as used in Section 3) to denote node features.

### Physics Philosophy Behind HANG

In Hamiltonian mechanics, the state evolution of a multi-object physical system, such as an electromagnetic field, double pendulum, or spring network [43, 44], adheres to well-established physical laws. For instance, in a system of charged particles, each particle generates an electromagnetic field that influences other particles. The Hamiltonian for such a system includes terms representing the kinetic and potential energy of each particle, along with terms representing the interactions between particles via their electromagnetic fields. In this paper, we propose a novel concept of information propagation between graph nodes, where interactions follow a similar Hamiltonian style.

In a Hamiltonian system, the position \(q\) and momentum \(p\) together constitute the phase space \((q,p)\), which comprehensively characterizes the system's evolution. In our HANG model, we process the raw node features of the graph using a linear input layer, yielding \(2r\)-dimensional vectors. Following the methodologies introduced in the GNN work [31, 35], we split the \(2r\) dimensions into two equal halves, with the first half serving as the feature (position) vector and the second half as "momentum" vector that guides the system evolution. Concretely, each node \(k\) is associated with an \(r\)-dimensional feature vector \(\mathbf{q}_{k}(0)=(q_{k}^{1},\ldots,q_{k}^{r})\) and an \(r\)-dimensional momentum vector \(\mathbf{p}_{k}(0)=(p_{k}^{1},\ldots,p_{k}^{r})\). Subsequently, \(\mathbf{q}_{k}(t)\) and \(\mathbf{p}_{k}(t)\) will evolve along with the propagation of information between graph nodes, with \(\mathbf{q}_{k}(0)\) and \(\mathbf{p}_{k}(0)\) serving as the initial conditions.

Following the modeling conventions in physical systems, we concatenate the feature positions of all \(|\mathcal{V}|\) vertices into a single vector, treating it as the system's generalized coordinate within a \(r|\mathcal{V}|\)-dimensional manifold, a process that involves index relabeling3.

Footnote 3: In multilinear algebra and tensor computation, vector components employ upper indices, while covector components use lower indices. We adhere to this convention.

\[q(t)=\left(q^{1}(t),\ldots q^{r|\mathcal{V}|}(t)\right)=\left(\mathbf{q}_{1}(t ),\ldots,\mathbf{q}_{|\mathcal{V}|}(t)\right). \tag{7}\]

This \(r|\mathcal{V}|\)-dimensional coordinate representation at each time instance provides a snapshot of the state of the graph system. Similarly, we concatenate all the "momentum" vectors at time \(t\) to construct an \(r|\mathcal{V}|\)-dimensional vector:

\[p(t)=\left(p_{1}(t),\ldots p_{r|\mathcal{V}|}(t)\right)=\left(\mathbf{p}_{1}(t ),\ldots,\mathbf{p}_{|\mathcal{V}|}(t)\right), \tag{8}\]

which can be interpreted as a generalized momentum vector for the entire graph system.

In physics, the system evolves in accordance with fundamental physical laws, and a conserved quantity function \(H(q,p)\) remains constant along the system's evolution trajectory. This conserved quantity is typically interpreted as the "system energy". In our HANG model, instead of defining an explicit energy function \(H(p,q)\) from a fixed physical law, we utilize a learnable energy function \(H_{\mathrm{net}}:\mathcal{G}\rightarrow\mathbb{R}^{+}\) parameterized by a neural network, referred to as the _Hamiltonian energy function_:

\[H_{\mathrm{net}}:\mathcal{G}\rightarrow\mathbb{R}^{+} \tag{9}\]

We allow the graph features to evolve according to a learnable Hamiltonian law analogous to basic physical laws. More specifically, we model the feature evolution trajectory as the following canonical Hamilton's equations, which is a restatement of (2):

\[\dot{q}(t)=\frac{\partial H_{\mathrm{net}}}{\partial p},\quad\dot{p}(t)=- \frac{\partial H_{\mathrm{net}}}{\partial q}, \tag{10}\]

with the initial features \((q(0),p(0))\in\mathbb{R}^{2r|\mathcal{V}|}\) at time \(t=0\) being the vectors after the raw node features transformation.

The neural ODE given by (10) can be trained and solved through integration to obtain the trajectory \((q(t),p(t))\). At the terminal time point \(t=T\), the system's solution is represented as \((q(T),p(T))\). We then apply the canonical projection map \(\pi\) to extract the nodes' concatenated feature vector \(q(T)\) as follows: \(\pi((q(T),p(T)))=q(T)\). This concatenated feature vector \(q(T)\) is subsequently decompressed into individual node features for utilization in downstream tasks. For this study, we employ backpropagation to minimize the cross-entropy in node classification tasks. The complete model architecture is depicted in Fig. 2, while a comprehensive summary of the full algorithm can be found in the Appendix G.
### Hamiltonian Energy Conservation

Referring to [45], it is known that the total Hamiltonian energy \(H_{\mathrm{net}}\) remains constant along the trajectory of its induced Hamiltonian flow. This principle is recognized as the law of energy conservation in a Hamiltonian system.

**Theorem 2**.: _If the graph system evolves in accordance with (10), the total energy \(H_{\mathrm{net}}(q(t),p(t))\) of the system remains constant. BIBO stability is achieved if \(H_{\mathrm{net}}\) remains bounded for all bounded inputs and, as \((q,p)\to\infty\), \(H_{\mathrm{net}}(q,p)\to\infty\)._

In light of Theorem 2, if our system evolves following (10), it adheres to the law of energy conservation. As a result, our model guarantees conservative stability.

**Definition 5** (Dirichlet energy [35]).: _The Dirichlet energy is defined on node features \(q(t)\) at time \(t\) of an undirected graph \(\mathcal{G}\) as_

\[\mathcal{E}(q(t))=\frac{1}{|\mathcal{V}|}\sum_{i}\sum_{j\in\mathcal{N}(i)}\| \mathbf{q}_{i}(t)-\mathbf{q}_{j}(t)\|^{2}. \tag{11}\]

Compared to GraphCON, which conserves Dirichlet energy over time \(t\)_under specific conditions_ as detailed in Section 3, the notion of Hamiltonian energy conservation is broader in scope. This is due to the fact that \(H_{\mathrm{net}}\) can always be defined as \(\mathcal{E}(q)\). Therefore, under the settings delineated in Section 3, GraphCON can be considered a particular instance of HANG.

## 5 Different Hamiltonian Energy Functions

In physical systems, the system is often depicted as a graph where two neighboring vertices with mass are connected by a spring of given stiffness and length [44]. The system's energy is thus related to the graph's topology. Similarly, in our graph system, the energy function \(H_{\mathrm{net}}\) involves interactions between neighboring nodes, signifying the importance of the graph's topology. There exist multiple ways to learn the energy function, and we present two examples below.

### Vanilla HANG

We define \(H_{\mathrm{net}}\) as a composition of two graph convolutional layers:

\[H_{\mathrm{net}}=\left\|\left(g_{\mathrm{gcn_{2}}}\circ\tanh\circ g_{\mathrm{ gcn_{1}}}\right)\left(q,p\right)\right\|_{2}, \tag{12}\]

where \(g_{\mathrm{gcn_{1}}}:\mathbb{R}^{2r\times|\mathcal{V}|}\to\mathbb{R}^{d\times| \mathcal{V}|}\) and \(g_{\mathrm{gcn_{2}}}:\mathbb{R}^{d\times|\mathcal{V}|}\to\mathbb{R}^{|\mathcal{ V}|}\) are two GCN [3] layers with different hidden dimensions. A \(\tanh\) activation function is applied between the two GCN layers, and \(\|\cdot\|_{2}\) denotes the \(\ell_{2}\) norm. We concatenate \(\mathbf{q}_{k}\) and \(\mathbf{p}_{k}\) for each node \(k\) at the input of the above composite function, resulting in an input dimension of \(2r\) per node. From Theorem 2, it follows that HANG exhibits BIBO stability. If \((q(t),p(t))\) were unbounded, the value of \(H_{\mathrm{net}}\) would also become unbounded, contradicting the energy conservation principle. In subsequent discussions, this invariant is referred to as HANG.

### Quadratic HANG (HANG-quad)

The general vanilla HANG with conservative stability does not possess Lyapunov stability. Other additional conditions may be required for Lyapunov stability. The following Lagrange-Dirichlet Theorem provides a sufficient condition.

**Theorem 3** (Lagrange-Dirichlet Theorem [46]).: _Let \(\mathbf{z}_{e}\) be a locally quadratic equilibrium of the natural Hamiltonian (2) with_

\[H=T(q,p)+U(q), \tag{13}\]

_where \(T\) is a positive definite, quadratic function of \(p\). Then \(\mathbf{z}_{e}\) is Lyapunov stable if the position of it is a strict local minimum of \(U(q)\)._

Theorem 3 implies that we can design an energy function \(H_{\mathrm{net}}\) such that the induced graph neural flow is both Lyapunov stable and energy conservative. For instance, we can define \(T\) as

\[T(q,p)=\sum_{k}\mathbf{p}_{k}^{\intercal}\left(\mathbf{A}_{G}(\mathbf{q}_{k}, \mathbf{q}_{k})\mathbf{A}_{G}^{\intercal}(\mathbf{q}_{k},\mathbf{q}_{k})+ \sigma\mathbf{I}\right)\mathbf{p}_{k}. \tag{14}\]
\(\sigma\) is a small positive number ensuring positive definiteness. The matrix \(\mathbf{A}_{G}\) can be the adjacency matrix or one learnable attention matrix. If \(U(q)\) is chosen to possess only a single, global minimum--such as \(U(q)=\|q\|\) with an \(\ell_{2}\) norm, as indicated in Section 3.1--this kind of stability does not necessarily guarantee adversarial robustness. One potential approach to this issue could be to set \(U(q)=\|\sin(q)\|\), thereby promoting Lyapunov stability at many local equilibrium points. However, this choice considerably restricts the form that \(U\) can take and may consequently limit the model's capacity. In the implementation, we set the function \(U(q)\) in \(H\) to be a single GAT [6] layer \(g_{\mathrm{gat}}:\mathbb{R}^{r\times|\mathcal{V}|}\rightarrow\mathbb{R}^{r \times\mathcal{V}}\) with a sin activation function followed by an \(\ell_{2}\) norm, i.e.,

The relationship between the stability of HANG and HANG-quad is summarized in Table 1.

\begin{table}
\begin{tabular}{c c c c c c} \hline \hline Graph Neural Flows & BIBO stability & Lyapunov stability & Structural stability & Conservative stability \\ \hline GRAND & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) \\ BLEND & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) \\ GraphCON & \(\mathcal{\prec}\) & \(\times\) & \(\times\) & \(\mathcal{\prec}\) \\ GraphBLE & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) & \(\mathcal{\prec}\) \\ HANG & \(\mathcal{\times}\) & \(\mathcal{\times}\) & \(\mathcal{\times}\) & \(\mathcal{\check{\prec}}\) \\ HANG-quad & \(\mathcal{\times}\) & \(\mathcal{\check{\vee}}\) & \(\mathcal{\times}\) & \(\mathcal{\check{\vee}}\) \\ \hline \hline \end{tabular}
\end{table}
Table 1: The stability summary for different graph neural flows where the \(\mathcal{\check{\prec}}\)denotes that stability is affirmed under additional conditions.

\begin{table}
\begin{tabular}{c c c c c c c c c} \hline \hline Dataset & Attack & HANG & HANG-quad & GraphCON & GraphBLE & GRAND & GAT & GraphSAGE & GCN \\ \hline \multirow{4}{*}{Cens} & _clean_ & 87.13\(\pm\)0.86 & 79.68\(\pm\)0.62 & 86.27\(\pm\)0.51 & 86.13\(\pm\)0.51 & 87.53\(\pm\)0.59 & 87.58\(\pm\)0.64 & 86.65\(\pm\)1.51 & 88.31\(\pm\)0.48 \\  & PGD & 78.37\(\pm\)1.84 & 79.05\(\pm\)0.42 & 42.81\(\pm\)0.30 & 40.04\(\pm\)0.68 & 39.65\(\pm\)1.32 & 38.27\(\pm\)2.73 & 35.43\(\pm\)5.05 & 35.83\(\pm\)0.71 \\  & TDIGA & 79.76\(\pm\)0.99 & 78.94\(\pm\)0.65 & 43.02\(\pm\)0.34 & 39.10\(\pm\)0.80 & 41.77\(\pm\)2.70 & 39.39\(\pm\)5.77 & 34.38\(\pm\)2.25 & 33.05\(\pm\)1.09 \\  & MetaCIA & 77.48\(\pm\)1.02 & 78.28\(\pm\)0.56 & 42.30\(\pm\)0.33 & 39.93\(\pm\)0.59 & 39.36\(\pm\)1.26 & 39.49\(\pm\)1.97 & 38.14\(\pm\)2.23 & 35.84\(\pm\)0.73 \\ \hline \multirow{4}{*}{Citeseer} & _clean_ & 74.11\(\pm\)0.62 & 71.85\(\pm\)0.48 & 74.84\(\pm\)0.49 & 69.62\(\pm\)0.56 & 74.98\(\pm\)0.45 & 67.87\(\pm\)4.97 & 63.22\(\pm\)9.14 & 72.63\(\pm\)1.14 \\  & PGD & 72.31\(\pm\)1.16 & 71.07\(\pm\)0.41 & 40.56\(\pm\)0.36 & 55.67\(\pm\)3.53 & 36.68\(\pm\)1.05 & 32.65\(\pm\)3.80 & 32.70\(\pm\)5.11 & 30.69\(\pm\)2.33 \\  & TDIGA & 71.21\(\pm\)0.52 & 71.09\(\pm\)0.40 & 36.67\(\pm\)1.25 & 34.74\(\pm\)4.68 & 36.67\(\pm\)1.25 & 30.53\(\pm\)3.57 & 36.11\(\pm\)2.94 & 21.10\(\pm\)2.35 \\  & MetaCIA & 72.92\(\pm\)0.66 & 71.60\(\pm\)0.48 & 43.86\(\pm\)2.12 & 45.60\(\pm\)3.41 & 46.23\(\pm\)2.01 & 37.68\(\pm\)4.04 & 35.75\(\pm\)5.50 & 35.86\(\pm\)0.68 \\ \hline \multirow{4}{*}{ConstroCS} & _clean_ & 96.16\(\pm\)0.09 & 59.27\(\pm\)1.12 & 95.10\(\pm\)0.12 & 93.93\(\pm\)0.48 & 95.08\(\pm\)0.12 & 92.84\(\pm\)0.41 & 93.09\(\pm\)0.39 & 93.33\(\pm\)0.73 \\  & PGD & 94.80\(\pm\)0.33 & 95.08\(\pm\)0.23 & 42.68\(\pm\)1.13 & 73.15\(\pm\)2.90 & 74.96\(\pm\)1.24 & 43.22\(\pm\)32.58 & 7.92\(\pm\)2.58 & 11.02\(\pm\)5.04 \\  & TDIGA & 95.40\(\pm\)0.13 & 95.09\(\pm\)0.09 & 7.92\(\pm\)4.11 & 43.13\(\pm\)3.58 & 5.05\(\pm\)1.43 & 16.08\(\pm\)15.74 & 6.47\(\pm\)4.25 & 3.61\(\pm\)1.77 \\  & MetaCIA & 94.85\(\pm\)0.31 & 94.83\(\pm\)0.28 & 67.79\(\pm\)4.04 & 73.98\(\pm\)8.26 & 84.31\(\pm\)4.26 & 52.01\(\pm\)24.21 & 7.82\(\pm\)2.69 & 15.70\(\pm\)3.95 \\ \hline \multirow{4}{*}{Pubuned} & _clean_ & 89.91\(\pm\)0.27 & 88.10\(\pm\)0.33 & 88.78\(\pm\)0.46 & 86.97\(\pm\)0.37 & 88.44\(\pm\)0.34 & 87.41\(\pm\)1.73 & 88.71\(\pm\)0.37 & 88.46\(\pm\)0.20 \\  & PGD & 81.81\(\pm\)1.94 & 89.67\(\pm\)0.57 & 45.06\(\pm\)0.51 & 46.06\(\pm\)0.97 & 44.61\(\pm\)2.78 & 49.84\(\pm\)12.99 & 44.62\(\pm\)6.49 & 99.03\(\pm\)0.10 \\  & TDIGA & 86.62\(\pm\)1.05 & 87.55\(\pm\)0.60 & 46.30\(\pm\)1.56 & 52.24\(\pm\)0.68 & 44.99\(\pm\)1.10 & 47.56\(\pm\)3.11 & 47.61\(\pm\)0.91 & 42.64\(\pm\)1.41 \\  & MetaCIA & 87.58\(\pm\)0.75 & 87.04\(\pm\)0.62 & 45.53\(\pm\)1.18 & 50.04\(\pm\)0.64 & 44.36\(\pm\)1.20 & 44.75\(\pm\)2.53 & 24.29\(\pm\)0.53 & 40.42\(\pm\)0.17 \\ \hline \hline \end{tabular}
\end{table}
Table 2: Node classification accuracy (%) on graph **injection, evasion, non-targeted** attack in **inductive** learning. The best and the second-best result for each criterion are highlighted in red and blue respectively.

Figure 2: The model architecture: each node is assigned a learnable “momentum” vector at time \(t=0\) which initializes the evolution of the system together with node features. The graph features evolve on following a _learnable_ law (10) derived from the \(H_{\mathrm{net}}\). At the time \(t=T\), we use \(q(T)\) as the final node feature. \(H_{\mathrm{net}}(q(t),p(t))\) is a learnable graph energy function.
## 6 Experiments

In this section, we conduct a comprehensive evaluation of our theoretical findings and assess the robustness of two conservative stable models: HANG and HANG-Quad. We compare their performance against various benchmark GNN models, including GAT [47], GraphSAGE [48], GCN [3] and other prevalent graph neural flows. We incorporate different types of graph adversarial attacks as described in Section 6.1 and Section 6.3. These attacks are conducted in a black-box setting, where a surrogate model is trained to generate perturbed graphs or features. For more experiments, we direct readers to Appendix C.

### Graph Injection Attacks (GIA)

We implement various GIA algorithms following the methodology in [49]. This framework consists of node injection and feature update procedures. Node injection involves generating new edges for injected nodes using gradient information or heuristics. We use the PGD-GIA method [49] to randomly inject nodes and determine their features with the PGD algorithm [50]. TDGIA [12] identifies topological vulnerabilities to guide edge generation and optimize a smooth loss function for feature generation. MetaGIA [49] performes iterative updates of the adjacency matrix and node features using gradient information. Our datasets include citation networks (Cora, Citeseer, Pubmed) [51], the Coauthor academic network [52], an Amazon co-purchase network (Computers) [52], and the Ogbn-Arxiv dataset [53]. For inductive learning, we follow the data splitting method in the GRB framework [54], with 60% for training, 10% for validation, and 20% for testing. Details on data statistics and attack budgets can be found in Appendix C.1. Targeted attacks are applied to the Ogbn-Arxiv and Computers datasets [49]. Additional results for various attack strengths and white-box attacks can be found in Appendix C.4 and Appendix C.3, respectively.

### Performance Results Under GIAs

Upon examining the results in Table 2 and Table 3 pertaining to experiments conducted under GIA conditions, the robustness of our proposed HANG and HANG-quad, is notably prominent across different GIA scenarios. Interestingly, GRAND, despite its Lyapunov stability as analyzed in Theorem 1, does not significantly outperform GAT under certain attacks. In contrast, HANG consistently displays robustness against attacks. Notably, HANG-quad exhibits superior performance to HANG on the Pubmed dataset under GIA perturbations, underscoring the effectiveness of integrating both Lyapunov stability and Hamiltonian mechanics to boost robustness. Although other graph neural

\begin{table}
\begin{tabular}{c c c c c c c c c c c c} \hline \hline Dataset & Phb Rate(\%) & HANG & HANG-quad & GraphCON & GraphBled & GRAND & GAT & GCN & RGCN & GCN-SVD & Pro-GNN \\ \hline \multirow{8}{*}{Publogs} & 0 & 94.77\(\pm\)1.07 & 94.63\(\pm\)1.06 & 93.14\(\pm\)0.84 & 85.13\(\pm\)2.22 & 95.57\(\pm\)0.44 & 95.35\(\pm\)0.20 & 95.69\(\pm\)0.38 & 95.22\(\pm\)0.14 & 95.31\(\pm\)0.18 & 93.20\(\pm\)0.64 \\  & 5 & 80.19\(\pm\)2.92 & 94.38\(\pm\)0.82 & 70.24\(\pm\)2.71 & 51.84\(\pm\)3.38 & 77.58\(\pm\)3.46 & 83.69\(\pm\)1.45 & 70.77\(\pm\)0.80 & 74.34\(\pm\)0.19 & 89.09\(\pm\)0.22 & 93.92\(\pm\)0.18 \\  & 10 & 73.49\(\pm\)4.32 & 92.94\(\pm\)1.56 & 71.87\(\pm\)1.71 & 56.54\(\pm\)2.30 & 77.99\(\pm\)1.35 & 76.32\(\pm\)0.85 & 70.22\(\pm\)1.13 & 71.04\(\pm\)0.34 & 81.24\(\pm\)0.49 & 89.42\(\pm\)1.09 \\  & 15 & 71.65\(\pm\)1.34 & 90.58\(\pm\)2.43 & 69.00\(\pm\)0.90 & 53.41\(\pm\)1.08 & 73.84\(\pm\)1.46 & 68.80\(\pm\)1.14 & 64.96\(\pm\)1.91 & 67.28\(\pm\)0.38 & 68.10\(\pm\)1.33 & 86.04\(\pm\)2.71 \\  & 20 & 66.27\(\pm\)3.39 & 89.19\(\pm\)3.72 & 64.06\(\pm\)2.31 & 52.18\(\pm\)0.54 & 69.14\(\pm\)1.32 & 51.50\(\pm\)1.63 & 51.27\(\pm\)1.23 & 59.89\(\pm\)0.34 & 57.33\(\pm\)1.35 & 79.56\(\pm\)5.68 \\  & 25 & 65.30\(\pm\)2.33 & 86.89\(\pm\)3.80 & 50.56\(\pm\)1.87 & 81.39\(\pm\)1.36 & 67.55\(\pm\)1.65 & 51.19\(\pm\)1.49 & 49.23\(\pm\)1.36 & 56.02\(\pm\)0.56 & 48.66\(\pm\)0.93 & 63.18\(\pm\)4.40 \\ \hline \multirow{8}{*}{Pubmed} & 0 & 85.08\(\pm\)0.20 & 85.23\(\pm\)0.14 & 86.65\(\pm\)0.17 & 84.02\(\pm\)0.26 & 85.06\(\pm\)0.26 & 83.73\(\pm\)0.40 & 87.19\(\pm\)0.09 & 86.16\(\pm\)0.18 & 83.44\(\pm\)0.21 & 87.33\(\pm\)0.18 \\  & 5 & 85.08\(\pm\)0.18 & 85.12\(\pm\)0.18 & 86.52\(\pm\)1.41 & 89.19\(\pm\)0.26 & 84.11\(\pm\)0.30 & 78.00\(\pm\)0.44 & 83.09\(\pm\)0.13 & 81.08\(\pm\)0.20 & 84.41\(\pm\)0.15 & 87.35\(\pm\)0.09 \\ \cline{1-1}  & 10 & 85.17\(\pm\)0.23 & 85.05\(\pm\)0.19 & 86.41\(\pm\)0.13 & 84.62\(\pm\)0.26 & 84.24\(\pm\)0.18 & 74.93\(\pm\)0.38 & 81.21\(\pm\)0.09 & 77.51\(\pm\)0.27 & 83.27\(\pm\)0.21 & 87.25\(\pm\)0.09 \\ \cline{1-1}  & 15 & 85.04\(\pm\)0.22 & 85.15\(\pm\)0.17 & 86.21\(\pm\)0.15 & 84.83\(\pm\)0.20 & 83.74\(\pm\)0.34 & 71.13\(\pm\)0.51 & 78.66\(\pm\)0.12 & 73.91\(\pm\)0.25 & 83.10\(\pm\)0.18 & 87.02\(\pm\)0.09 \\ \cline{1-1}  & 20 & 85.20\(\pm\)0.19 & 85.03\(\pm\)0.19 & 86.00\(\pm\)0.18 & 84.89\(\pm\)0.45 & 83.88\(\pm\)0.20 & 68.21\(\pm\)0.96 & 73.73\(\pm\)0.19 & 71.18\(\pm\)0.31 & 83.01\(\pm\)0.22 & 87.09\(\pm\)0.10 \\ \cline{1-1}  & 25 & 85.06\(\pm\)0.17 & 84.99\(\pm\)0.16 & 86.04\(\pm\)0.14 & 85.07\(\pm\)0.15 & 83.66\(\pm\)0.25 & 65.41\(\pm\)0.77 & 75.05\(\pm\)0.17 & 67.95\(\pm\)0.15 & 82.72\(\pm\)0.18 & 86.71\(\pm\)0.09 \\ \hline \hline \end{tabular}
\end{table}
Table 4: Node classification accuracy (%) under modification, poisoning, non-targeted attack (Metattack) in **transductive** learning.

\begin{table}
\begin{tabular}{c c c c c c c c c} \hline \hline Dataset & Attack & HANG & HANG-quad & GraphCON & GraphBled & GRAND & GAT & GraphSAGE & GCN \\ \hline \multirow{3}{*}{Computers} & PGD & 90.83\(\pm\)0.53 & 87.53\(\pm\)0.99 & 74.01\(\pm\)4.87 & 89.33\(\pm\)0.56 & 65.75\(\pm\)5.00 & 35.72\(\pm\)1.07 & 18.33\(\pm\)0.57 & 17.80\(\pm\)0.06 \\  & TDGIA & 90.88\(\pm\)0.50 & 87.75\(\pm\)0.53 & 80.11\(\pm\)3.21 & 82.81\(\pm\)0.63 & 77.18\(\pm\)7.60 & 66.09\(\pm\)17.07 & 50.89\(\pm\)1.30 & 66.51\(\pm\)3.90 \\  & MetaGIA & 90.73\(\pm\)0.64 & 89.32\(\pm\)0.48 & 87.95\(\pm\)2.65 & 89.27\(\pm\)0.57 & 81.85\(\pm\)4.66 & 60.28\(\pm\)1.11 & 37.51\(\pm\)4.49 & 36.19\(\pm\)1.00 \\ \hline \multirow{3}{*
flows might demonstrate a range of improved performance compared to conventional GNN models under GIA, the degree of improvement is not consistently distinct. Despite the pronounced association between conservative stability in Hamiltonian systems and adversarial robustness, a clear relationship between adversarial robustness with the other stability of the graph neural flows outlined in Table 1 is not immediately discernible. The performance differential between HANG variants and other graph neural flows further underscores the potential of our proposed Hamiltonian models in enhancing robustness against GIA attacks.

### Graph Modification Attacks

To evaluate the robustness of our proposed conservative models, we conducted graph modification adversarial attacks using the Metattack method [17]. We followed the attack setting described in the Pro-GNN [55] and utilized the perturbed graph provided by the library [56] to ensure a fair comparison. The perturbation rate, which indicates the proportion of altered edges, was incrementally varied in 5% increments from 0% to 25%. For comparison, we also considered other defense models for GNNs, namely Pro-GNN [55], RGCN [57], and GCN-SVD [58]. We report the results of baseline models from [55].

### Performance Results Under Modification/Poisoning/Transductive Attacks

In the case of the Polblogs dataset [59], as shown in Table 4, our proposed HANG-quad model demonstrates superior performance compared to other methods, including existing defense models. This result indicates that incorporating Lyapunov stability indeed enhances HANG's robustness against graph modification and poisoning attacks. For the Pubmed dataset, we note that the impact of Meta-attacks of varying strengths on _all_ graph neural flows, including our proposed ones, is negligible. Conversely, traditional GNN models such as GAT, GCN, and RGCN are marginally affected as the attack strength escalates. This observation underlines the robustness of graph neural flows, including our proposed models, against Meta-attacks on this dataset.

### Combination with other defense mechanisms

It merits noting that our models, HANG and HANG-quad, can be readily integrated with additional defense mechanisms against adversarial attacks. These include Adversarial Training (AT) [60] and other preprocessing methods such as GNNGUARD [61]. This integration can further bolster the robustness of the HANG model. To validate this enhancement, extensive experiments are conducted, with results detailed in Appendix C.8 and Appendix C.9.

## 7 Conclusion

In this paper, we conducted a comprehensive study on stability notions in the context of graph neural flows and made significant findings. While Lyapunov stability is frequently employed, it alone may not suffice in guaranteeing robustness against adversarial attacks. With a grounding in foundational physics principles, we proposed a shift towards conservative Hamiltonian neural flows for crafting GNNs resilient against adversarial attacks. Our empirical comparisons across diverse neural flow GNNs, as tested on multiple benchmark datasets subjected to a range of adversarial attacks, have further corroborated this proposition. Notably, GNNs that amalgamate conservative Hamiltonian flows with Lyapunov stability exhibited marked enhancement in their robustness metrics. We are optimistic that our work will inspire further research into marrying physics principles with machine learning paradigms for enhanced security.

## 8 Acknowledgments and Disclosure of Funding

This research is supported by the Singapore Ministry of Education Academic Research Fund Tier 2 grant MOE-T2EP20220-0002, and the National Research Foundation, Singapore and Infocomm Media Development Authority under its Future Communications Research and Development Programme. To improve the readability, parts of this paper have been grammatically revised using ChatGPT [62].