Graph Neural Networks Gone Hogwild

## Abstract

Abstract Message passing graph neural networks (GNNs) would appear to be powerful tools to learn distributed algorithms via gradient descent, but generate catastrophically incorrect predictions when nodes update asynchronously during inference.
This failure under asynchrony effectively excludes these architectures from many potential applications, such as learning local communication policies between resource-constrained agents in, e.g., robotic swarms or sensor networks.
In this work we explore why this failure occurs in common GNN architectures, and identify “implicitly-defined” GNNs as a class of architectures which is provably robust to partially asynchronous “hogwild” inference, adapting convergence guarantees from work in asynchronous and distributed optimization, e.g., Bertsekas ( 1982 ); Niu et al. ( 2011 ) .
We then propose a novel implicitly-defined GNN architecture, which we call an energy GNN . We show that this architecture outperforms other GNNs from this class on a variety of synthetic tasks inspired by multi-agent systems, and achieves competitive performance on real-world datasets.

## 1 Introduction

Message-passing graph neural networks (GNNs) compute embeddings for nodes in graph-structured data by using repeated message passing within nodal neighborhoods (Hamilton, ; Gilmer et al., ; Kipf and Welling, ; Veličković et al., ).
This effectively produces a deep neural network architecture whose computation graph reflects the structure of the input graph.
The inherent locality in the generation of node embeddings appears at first glance to lend itself to distributed and decentralized execution, which is desirable in a variety of contexts.

Consider a group of agents (e.g., robots or “motes” in a sensor network) which need to collectively perform a task, but that might be unreliable, have limited range of communication, possess scarce local computational resources, and lack central control.
GNNs are appealing as a way to learn local communication policies that solve the distributed problem, where each agent corresponds to a node in the graph and each edge to a local communication channel.
One could imagine learning an algorithm that enables a swarm of robots to localize themselves, or a collection of resource-constrained edge devices to collectively estimate environmental conditions.
Another application where distributed computation is attractive is in GNN inference over large graphs, where nodes or node collections are managed on distinct machines (respecting the graph connectivity).
This is especially relevant for GNN deployment on resource-constrained devices.
Distributed inference may also facilitate privacy in settings where nodes correspond to entities such as individuals in a social network, by enabling local inference and precluding excessive data transmission.

There is a significant obstacle in applying GNNs in these kinds of applications: distributed and decentralized computation is often asynchronous, but most GNN architectures implicitly assume synchrony in the execution of message passing at inference time.
When nodes update asynchronously, messages received from neighbors may be out of date, and the effective architecture of the GNN diverges catastrophically from the training architecture; this means the output of the GNN can exhibit arbitrarily large errors.
Figure 1 illustrates the issue on a toy problem.
If it were possible to perform such “hogwild” inference (nodding to the asynchronous optimization work of Niu et al. ()), then GNNs would unlock solutions for a significantly wider set of problem domains.

Figure: Figure 1: Modifications to the computation graph of an $L$-layer message-passing GNN resulting from asynchronous, distributed, per-node inference. (a) An undirected linear graph with binary features written in the node body. (b) Synchronous per-layer execution of a 3-layer GNN on the graph from (a); all nodes update at the same time, using neighbor information from the previous layer. Arrows represent network weights. For demonstration, weights are all equal to 1, and the value of a node embedding at the next layer is the sum of incoming values. (c) Asynchronous inference of the GNN from (b). Nodes update at random times and can use neighbor information corresponding to the incorrect layer, which introduces modifications to the computation graph. Gray arrows correspond to connections that were removed from the original computation graph, and red arrows are unintended connections resulting from asynchrony. (d) To demonstrate the effect of asynchrony, we show the output of the GNN varies significantly over asynchronous runs with different node update orderings.
Refer to caption: extracted/5699996/figures/fig1.png

In this work, we characterize a class of “implicitly-defined” GNNs which do not exhibit this catastrophic behavior under asynchronous inference.
In these GNNs, node embeddings are implicitly defined as the solution to an optimization problem (Gu et al., ; Liu et al., ; Yang et al., ; Zhu et al., ; Ma et al., ; Scarselli et al., ).
This is in contrast to GNNs in which node embeddings are computed as the output of a specific forward computation, such as graph attention networks (GAT, Veličković et al. ()) and graph convolutional networks (GCN, Kipf and Welling ()); we refer to these architectures as “explicitly-defined”.
Implicitly-defined GNNs are provably robust to asynchrony given lenient assumptions on staleness and per-node update frequency, making them ideally suited to distributed, asynchronous inference. However, the diversity of existing implicitly-defined GNN architectures is limited.
We contribute a novel architecture, which we call energy GNN, which exposes a rich parameterization of optimization problems by using input-convex neural networks defined on nodal neighborhoods.
We show that energy GNNs outperform other implicitly-defined GNNs on a variety of synthetic tasks motivated by problems which are relevant to multi-agent systems, where synchronous inference may be undesirable.
We also achieve competitive performance on tasks with benchmark graph datasets, certifying the merit of our approach even as a stand-alone GNN architecture.

Our contributions are as follows:

- •
We present the first unified framework for GNNs which can execute inference asynchronously. We refer to this class of GNNs as being *implicitly defined*, and articulate sufficient conditions that such architectures must satisfy. We then contextualize all major GNN architectures and identify which are implicitly defined, and which are not.
- •
We propose a novel implicitly defined architecture called the energy GNN. Unlike existing implicitly-defined GNNs, energy GNNs 1) incorporate edge features, a crucial property in many applications, 2) enable the generation of neighbor-specific messages, 3) enable the utilization of neighborhood attention, and 4) use neural networks for generation of messages and embeddings.
- •
We provide conditions under which there are guarantees for the convergence of our architecture under partially asynchronous inference, adapting classical results from the distributed optimization literature.
- •
We present several empirical studies which span synthetic datasets crafted to tease out important distinctions between GNNs which are/aren’t amenable to asynchronous inference, and several widely used GNN benchmarks to assess performance in a cross-referenced manner with existing results published in the literature.

## 2 Preliminaries

### 2.1 Partially Asychronous Algorithms

Computational models for asynchronous algorithms vary depending on the constraints imposed on the sequencing or frequency of computation or communication.
We consider partial asynchronism as defined by Bertsekas and Tsitsiklis (), which, informally, places bounds on two key characteristics: the time between updates for each node, and the amount by which data retained at any node can be out of date from the perspective of some other node(s).

Concretely, consider a collection of $n$ nodes, which could be a set of processors carrying out a distributed computation.
Each node has hidden state denoted ${{\bm{h}}_{i}\in\mathbb{R}^{k}}$, which corresponds to one “block” (row) of the aggregate state ${{\bm{H}}:=({\bm{h}}_{1},\ldots,{\bm{h}}_{n})^{T}\in\mathbb{R}^{n\times k}}$.
The distributed computation we aim to execute is an asynchronous iterative algorithm, where each node $i=1,2,\dots,n$ iteratively updates its state according to ${\bm{h}}_{i}:=f^{i}({\bm{H}})$ using some node-specific update function $f^{i}:\mathbb{R}^{n\times k}\to\mathbb{R}^{k}$.

Asynchrony in the execution of the algorithm is accounted for as follows. We are given a set $T^{i}\subseteq\{0,1,2,\dots\}$ of times at which each node $i$ is updated.
Additionally, for each ${t\in T^{i}}$ we are given variables $0\leq\tau_{j}^{i}(t)\leq t$ which represent the time of node $j$’s last update from the perspective of node $i$, at time $t$.
Note that $\tau^{i}_{j}(t)\in T^{j}$.
The quantities ${t-\tau_{j}^{i}(t)\in[0,t]}$ can then be interpreted as the amount (in time) by which information associated with node $j$ is outdated or “stale” when used in the update of ${\bm{h}}_{i}$ at time $t$.
For instance, the state ${T^{3}=\{0,4,7\}}$, ${T^{1}=\{1,2,3\}}$, ${\tau_{1}^{3}(4)=2}$ means that when node 3 executes its second update (at ${t=4}$), it uses the stale value of node 1 that was computed at ${t=2}$, rather than the most recent value at ${t=3}$.

We can now describe the asynchronous algorithm using the following node update equations:

$$ $\displaystyle{\bm{h}}_{i}(t\!+\!1)$ $\displaystyle={\bm{h}}_{i}(t)\;\;\text{if}\;t\notin T^{i},$ (1) $\displaystyle{\bm{h}}_{i}(t\!+\!1)$ $\displaystyle=f^{i}({\bm{h}}_{1}(\tau_{1}^{i}(t)),\dots,{\bm{h}}_{n}(\tau_{n}^ {i}(t)))\;\;\text{if}\;t\in T^{i},$ (2) $$

for $t\geq 0$.
Partial asynchronism corresponds to the following assumptions:

###### Assumption 2.1 (Partial Asynchronism).

(Bertsekas and Tsitsiklis () Assumption 5.3)
There exists an integer ${B>0}$ such that:

- (a)
(Bounded time to next update) For every node $i$ and for every $t\geq 0$, at least one of the elements of the set $\{t,t+1,\dots,t+B-1\}$ belongs to $T^{i}$.
- (b)
(Bounded staleness) There holds $t-B<\tau_{j}^{i}(t)\leq t$,
for all $i,j$, and all $t\geq 0$ belonging to $T^{i}$.
- (c)
There holds $\tau_{i}^{i}(t)=t$ for all $i=1,2,\dots,n$ and $t\in T^{i}$.

###### Assumption 2.1 (Partial Asynchronism).

### 2.2 Graph Neural Networks

In this section, we provide a brief overview of GNNs.
Consider a directed graph with a collection of $n$ vertices ${V=\{1,...,n\}}$, and edges ${E\subseteq V\times V}$.
The connectivity of the graph is contained in its adjacency matrix ${\bm{A}}\in\{0,1\}^{n\times n}$, where ${\bm{A}}_{i,j}=1$ if there is an edge from node $i$ to node $j$, and 0 0 otherwise.
The graph may also have associated node and edge features ${{\bm{X}}\in\mathbb{R}^{n\times p}}$ and ${{\bm{E}}\in\mathbb{R}^{|E|\times r}}$, so we use ${\mathcal{G}}=({\bm{A}},{\bm{X}},{\bm{E}})$ to denote the graph (and its features).
We use $\mathcal{N}(i)$ to denote the set of neighbors of node $i$, and $|{\mathcal{G}}|$ to denote the number of nodes in a graph ${\mathcal{G}}$.

We are given a dataset ${\mathcal{D}=\{({\mathcal{G}},{\bm{Y}})^{q}\}_{q=1}^{|\mathcal{D}|}}$, where each ${\mathcal{G}}^{q}$ is associated with a node-level prediction target ${{\bm{Y}}^{q}\in\mathbb{R}^{n_{q}\times J}}$, with $n_{q}=|{\mathcal{G}}^{q}|$.
We focus on node-level predictions, since they often underlie graph-level predictions.

In a GNN, each node $i=1,\dots,n$ is associated with a $k$-vector embedding ${\bm{h}}_{i}$ which is updated through iterations (or layers) of message passing; we use ${\bm{H}}:=({\bm{h}}_{1},\dots,{\bm{h}}_{n})^{T}\in\mathbb{R}^{n\times k}$ to denote the aggregate embedding of nodes in the graph.
The dimensionality $k$ can vary across layers but, without loss of generality, we assume in this work that it is fixed.
Each layer ${\ell\in\{1,...,L\}}$ in a GNN defines a parameterized embedding function $f_{\theta}^{\ell}:{\mathcal{G}}_{i}\times\mathbb{R}^{|{\mathcal{G}}_{i}|\times
k
}\rightarrow\mathbb{R}^{k}$. $f_{\theta}^{\ell}$ is applied independently for each node $i=1,\dots,n$ and takes as argument a subgraph ${\mathcal{G}}_{i}\subseteq{\mathcal{G}}$ corresponding to the neighborhood of node $i$ and the subgraph’s node embeddings from the previous layer ${\bm{H}}^{\ell-1}\in\mathbb{R}^{|{\mathcal{G}}_{i}|\times k}$,
and returns an embedding ${\bm{h}}^{\ell}_{i}$ for node $i$.
The new embeddings ${\bm{H}}^{\ell}$ are used in the next layer.
Typically, ${\bm{H}}^{0}$ is initialized as the node features ${\bm{X}}$.
For convenience, we omit the layer index $\ell$ from $f_{\theta}^{\ell}$ in contexts where there is only one parameterized layer, or when referring generally to a message passing layer.

A layer can be described as performing message passing, where each node $i$ constructs a “message” ${\bm{m}}_{i}$ by aggregating information (e.g., embeddings from the previous layer and node/edge features) from nodes $j\in\mathcal{N}(i)$ in its local neighborhood,
then uses this aggregate message to update its embedding.
In the final layer, a readout function ${o_{\phi}:\mathbb{R}^{k}\rightarrow\mathbb{R}^{J}}$ is applied to each embedding, which results in node predictions ${\hat{{\bm{Y}}}=(o_{\phi}({\bm{h}}^{L}_{1}),...,o_{\phi}({\bm{h}}^{L}_{n}))^{T
}\in\mathbb{R}^{n\times J}}$.

Different GNN architectures are distinguished by the parameterization of the update function $f_{\theta}$.
In Appendix A, we describe GAT and GCN, two widely-used architectures.

## 3 Explicitly-defined vs. Implicitly-defined GNNs

In a conventional feed-forward message passing architecture that consists of $L$ layers, it is assumed that all update functions for a given layer are applied in a synchronized manner.
That is, each node embedding is derived from the previous layer embeddings.
In distributed systems synchrony is either impossible or comes with substantial overhead.
If nodes update asynchronously according to Equations 1 and 2, GNN architectures that assume synchronicity will fail catastrophically and nondeterministically because the architecture at inference time is different than it was during training.
This effect is illustrated in Figure 1.

This failure motivates us to distinguish between two types of GNN architectures: *explicitly-defined* GNNs which specify a specific feed-forward layer-wise computation, and *implicitly-defined* GNNs in which the layer-wise message passing updates correspond to iterations toward a fixed point.
We further subdivide implicitly-defined GNNs into two types: *fixed-point* GNNs and the important special case of *optimization-based* GNNs.
Explicitly-defined GNNs are widely-used but are susceptible to failure under asynchrony, while implicitly-defined GNNs are robust (as we discuss in Section 4).

### 3.1 Fixed-point GNNs

Fixed-point GNNs obtain node embeddings as the fixed point of a contractive message passing function.
Using ${\bm{h}}\in\mathbb{R}^{nk}$ to denote the unrolled embeddings ${\bm{H}}$ and taking ${F_{\theta}:{\mathcal{G}}\times\mathbb{R}^{nk}\to\mathbb{R}^{nk}}$ to be the aggregate update of the hidden state from applying $f_{\theta}$ to the neighborhood subgraph of each node $i=1,...n$ in a graph ${\mathcal{G}}$, a fixed-point GNN iterates $F_{\theta}$ until numerical convergence of node embeddings, i.e., ${F_{\theta}({\bm{h}})\approx{\bm{h}}}$.

Convergence to a unique fixed point is guaranteed provided the update function $F_{\theta}$ is a contraction map with respect to the embeddings, i.e. $||F_{\theta}({\bm{h}})-F_{\theta}({\bm{h}}^{\prime})||\leq\mu||{\bm{h}}-{\bm{h
}}^{\prime}||$ holds for $0<\mu<1$ and $\forall{\bm{h}},{\bm{h}}^{\prime}\in\mathbb{R}^{nk}$.
Since it is difficult to design non-trivial parameterizations of $F_{\theta}$ which are contractive by construction,
existing fixed-point GNN architectures (e.g. IGNN (Gu et al., ), EIGNN (Liu et al., ), and APPNP (Gasteiger et al., )) are limited in diversity. These architectures use linear transformations of the node embeddings to form messages ${\bm{m}}_{i}$, where the parameters can easily be constrained such that $F_{\theta}$ is contractive.
Gu et al. () apply a component-wise non-expansive nonlinearity to compute updated node embeddings from these messages, while Liu et al. (); Gasteiger et al. () take ${\bm{m}}_{i}$ as the node embeddings directly.

The original ‘nonlinear GNN’ proposed by Scarselli et al. (), in which $f_{\theta}$ contains general multi-layer neural networks operating on the input data, represents another approach to parameterizing a fixed-point GNN.
Their method encourages rather than guarantees contraction via
a multi-objective problem that includes the norm of the Jacobian of the update function as a quantity to be minimized.
This heuristic can work in practice, but the sequence of iterates does not definitively converge, particularly if node embeddings are initialized far from the fixed point (as the norm of the Jacobian is only penalized at the fixed point).
We thus do not consider this non-linear GNN to be a true fixed-point GNN.
We provide more details on existing fixed-point GNN architectures in appendix A.

### 3.2 Optimization-based GNNs

Optimization-based GNNs are an important special case of fixed-point GNNs, obtaining node embeddings by minimizing a convex scalar-valued graph function ${E_{\theta}:{\mathcal{G}}\times\mathbb{R}^{n\times k}\to\mathbb{R}}$ with respect to the aggregate graph node embeddings ${\bm{H}}$, where $\theta$ are parameters:

$$ $\displaystyle{\bm{H}}^{*}=\textstyle\operatorname*{arg\,min}_{{\bm{H}}}E_{ \theta}({\mathcal{G}},{\bm{H}}).$ (3) $$

In this work, we assume $E_{\theta}$ is separable per node, with each node computing its contribution to $E_{\theta}$ using local information. That is, $E_{\theta}$ can be expressed as:

$$ $\displaystyle E_{\theta}({\mathcal{G}},{\bm{H}})=\textstyle\sum_{i=1}^{n}e_{ \theta}\left({\bm{h}}_{i},{\bm{x}}_{i},\{{\bm{h}}_{j},{\bm{x}}_{j},{\bm{e}}_{ ij}\mid j\in\mathcal{N}(i)\}\right)=\textstyle\sum_{i=1}^{n}e^{i}_{\theta}.$ (4) $$

Here $e^{i}_{\theta}$ depends only on information from node $i$’s neighborhood, and we further assume it is twice differentiable with a Hessian of bounded norm.
Crucially, due to the dependence of each $e^{i}_{\theta}$ on local information, gradient-based minimization of $E_{\theta}$ can be expressed per node via message passing:

$$ $\displaystyle{\bm{h}}_{i}(t+1)$ $\displaystyle={\bm{h}}_{i}(t)-\alpha\textstyle\sum_{j\in\mathcal{N}(i)\cup\{i \}}{\bm{g}}_{ji}(t)$ (5) $\displaystyle{\bm{g}}_{ji}(t)$ $\displaystyle:=\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}_{j}(t),\{{\bm{h}}_ {j^{\prime}}(t)|j^{\prime}\in\mathcal{N}(j)\}),$ (6) $$

where $\alpha\in\mathbb{R}_{>0}$, and node and edge features are omitted for clarity.
We assume that at time $t$, node $i$ obtains the values ${\bm{g}}_{ji}$ and ${\bm{h}}_{j}$ (needed to compute ${\bm{g}}_{ii}$) from neighbors $j\in\mathcal{N}(i)$. The number of iterations is dictated by the convergence of the embeddings to a fixed point.

Existing optimization-based GNNs (Yang et al., ; Zhu et al., ; Ma et al., ) use an objective $E_{\theta}$ where node embeddings arise as:

$$ $\displaystyle{\bm{H}}^{*}=\textstyle\operatorname*{arg\,min}_{{\bm{H}}}\gamma| |{\bm{H}}-g_{\theta}({\bm{X}})||_{F}^{2}+\beta\text{tr}({\bm{H}}^{T}{\bm{L}}{ \bm{H}}),$ (7) $$

where ${g_{\theta}:\mathbb{R}^{n\times p}\rightarrow\mathbb{R}^{n\times k}}$, ${\bm{L}}$ is a function of ${\bm{A}}$ and can be viewed as a generalized incidence matrix (assumed to be symmetric and positive semi-definite), $\gamma$ and $\beta$ are constants, and $\text{tr}()$ is the trace.
The first term drives node embeddings to approximate some function of the node features ${\bm{X}}$, and the second term is a regularizer that rewards smoothness of neighboring node embeddings.
We refer to optimization-based GNNs which use this form of objective as GSDGNNs, since minimizing the objective can be interpreted as performing graph signal denoising (GSD).

Previous work has shown that the node embeddings obtained by the fixed-point message passing scheme defined by APPNP (Gasteiger et al., ) and by EIGNN (Liu et al., ) correspond directly to minimization of this form of objective (Ma et al., ; Yang et al., ).
This correspondence between fixed-point GNNs and optimization-based GNNs does not always exist; while the node embeddings in optimization-based GNNs can always be expressed as satisfying a fixed-point equation (i.e., the gradient of the convex graph function is equal to zero at the solution), not all fixed-point equations correspond with convex optimization problems.
For this reason, we distinguish between these two types of implicit GNNs.

## 4 Hogwild GNN Inference

Executing the layer updates during GNN inference is typically performed synchronously, as depicted in Figure 1; each node embedding is derived from the previous layer embeddings.
In this section, we discuss GNN inference under the asynchronous, “hogwild” (Niu et al., ) execution model presented in Section 2.1.

In the general asynchronous update in Equation 2 it is assumed that node $i$ has access to all of the values ${\bm{h}}_{1},...{\bm{h}}_{n}$ required for performing its update.
In general these values could be obtained, for example, by providing nodes access to a shared memory structure which contains these values.
However, in this work, we are interested in decentralized, distributed, and asynchronous GNN embedding updates, where neither centralized memory nor a centralized controller are used.
We assume that individual nodes store their own values ${\bm{h}}_{i}$ (and ${\bm{x}}_{j}$, ${\bm{e}}_{ij}$), and if these values are needed by other nodes to perform their update, they must obtain them by communicating with neighbors.
Our aim in this section is to demonstrate that node updates can be performed only using information that can be obtained through local communication.

### 4.1 Explicitly-defined GNN and Fixed-point GNN Inference under Partial Asynchronism

Without loss of generality, assume that the embedding dimension $k$ is fixed for all layers of parameterized node update functions.
We do not write $f_{\theta}$ indexed by layer, but this is straightforwardly generalized to the case of layer-specific parameters and functions described in Section 2.2.
Without loss of generality, assume the embedding update function $f_{\theta}$ is continuously differentiable, so that the following restriction can be stated:

$$ $\displaystyle j\notin\mathcal{N}(i)\;\implies\frac{\partial f_{\theta}}{ \partial{\bm{h}}_{j}}({\bm{z}})=0\quad\forall{\bm{z}}\in\mathbb{R}^{k}.$ (8) $$

With this restriction, the general node update from Equation 2 can be adapted to describe partially asynchronous message passing, in which node updates are performed using only information from a node’s local neighborhood. In particular, for node $i$ and for $t\geq 0$ and $t\in T^{i}$, the update equation is:

$$ $\displaystyle{\bm{h}}_{i}(t+1)$ $\displaystyle=f_{\theta}({\bm{h}}_{i}(t),\{{\bm{h}}_{j}(\tau_{j}^{i}(t))\mid j \in\mathcal{N}(i)\})\,$ (9) $$

where we omit node and edge features for clarity.
Note the crucial difference introduced by asynchrony: the neighbor data ${\bm{h}}_{j}(\tau_{j}^{i}(t))$ may be out of date, corresponding to a different layer in the network.
For explicitly-defined GNNs such as GCN or GAT, the number of iterations $|T_{i}|$ executed by each node is fixed and equal to the number of layers $L$ in the GNN.
For fixed point GNNs, the number of iterations is not pre-specified.

Since explicitly-defined GNNs implement a specific feed-forward neural network architecture, partially asynchronous inference using Equation 9 corrupts the computation performed by the network.
We illustrate this in Figure 1, where partial asynchrony results in a (different) computation graph with some connections removed, and new connections that are not present in the original synchronous computation graph.
This means that there are no convergence guarantees under partial asynchrony, and in general the final node embeddings may vary significantly with respect to the particular node update sequence.

In contrast, fixed-point GNN architectures in which $F_{\theta}$ is contractive with respect to node embeddings ${\bm{h}}$ are provably robust to partially asynchronous inference.
In particular, they satisfy the assumptions of the following proposition (Bertsekas, ).

###### Proposition 4.1 .

If $F_{\theta}:{\mathcal{G}}\times\mathbb{R}^{nk}\rightarrow\mathbb{R}^{nk}$ is contractive with respect to node embeddings ${\bm{h}}$, then under the bounded staleness conditions in 2.1, the fixed-point iteration of Equation 9 converges.

###### Proposition 4.1 .

### 4.2 Optimization-based GNN Inference under Partial Asynchrony

In order to examine inference of optimization-based GNNs under partial asynchrony, we assume gradient-based optimization is used in computing node embeddings.
Recall that optimization-based GNNs are performing a minimization of a separable objective $E_{\theta}$, as defined in Equation 4.
We thus state the following restriction, analogous to the restriction in Equation 8:

$$ $\displaystyle j\notin\mathcal{N}(i)\implies\frac{\partial e^{i}_{\theta}}{ \partial{\bm{h}}_{j}}({\bm{z}})=0\;\text{for all}\;\;{\bm{z}}\in\mathbb{R}^{k}.$ (10) $$

We could then naively adapt the general node update equations from Equation 1 to describe partially asynchronous, gradient-based minimization of $E_{\theta}$ as follows:

$$ $\displaystyle{\bm{h}}_{i}(t+1)$ $\displaystyle={\bm{h}}_{i}(t)-\alpha\textstyle\sum_{j\in\mathcal{N}(i)\cup\{i \}}{\bm{g}}_{ji}(\tau^{i}_{j}(t))$ (11) $\displaystyle{\bm{g}}_{ji}(\tau^{i}_{j}(t))$ $\displaystyle:=\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}_{j}(\tau^{i}_{j}(t )),\{{\bm{h}}_{j^{\prime}}(\tau^{i}_{j^{\prime}}(t))|j^{\prime}\in\mathcal{N}( j)\}),$ (12) $$

where $\alpha\in\mathbb{R}_{>0}$ is the step size, and node and edge features are omitted for clarity.
This formulation would allow us to directly cite a proof of convergence from Bertsekas and Tsitsiklis (), as we did in Section 4.1 for partially asynchronous fixed-point GNN inference.

However, recall that for a node to perform an update using only local communication, we previously assumed that ${\bm{g}}_{ji}$ was obtained by node $i$ from its neighbor $j$.
With the formulation in Equations 11 and 12, node $j$ cannot provide node $i$ with ${\bm{g}}_{ji}$ since the value depends on node $i$’s view of the embeddings (and features) of the neighbors of node $j$, rather than deferring to node $j$’s view of its neighbors.
That is, node $i$ needs access to information about its 2-hop neighbors in addition to its 1-hop neighbors. Since we assume communication is only possible with 1-hop neighbors, 2-hop neighbor information would need to be forwarded by direct neighbors of node $i$.
This inflates the cost of communication, requiring a number of bits per transmission which is proportional to $|\mathcal{N}(j)|$
Furthermore, as the number of neighbors that are shared between node $i$ and node $j$ increases, the contents of these messages become increasingly redundant.

Instead, we preserve fully local communication and fixed-size transmissions where neighbors $j$ send fixed-size packets containing only (${\bm{h}}_{j}$, ${\bm{g}}_{ji}$) to node $i$ by defining ${\bm{g}}_{ji}(\tau^{i}_{j}(t))$ as follows:

$$ $\displaystyle{\bm{g}}_{ji}(\tau^{i}_{j}(t))$ $\displaystyle:=\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}_{j}(\tau^{i}_{j}(t )),\{{\bm{h}}_{j^{\prime}}(\tau^{j}_{j^{\prime}}(\tau^{i}_{j}(t))):j^{\prime} \in\mathcal{N}(j)\}).\,$ (13) $$

The crucial difference between equations Equation 12 and Equation 13 is that instead of ${\bm{g}}_{ji}(\tau^{i}_{j}(t))$ depending on node $i$’s view of 2-hop neighbors $j^{\prime}$ at time $t$, it now depends on neighbor $j$’s view of its neighbors at time $\tau^{i}_{j}(t)$, the time corresponding to node $i$’s view of $j$ at time $t$.

###### Proposition 4.2 .

If $E_{\theta}$ is strongly convex and separable per node, then for a sufficiently small step size and the bounded staleness conditions in 2.1, the optimization procedure of Equations 11 and 13 will converge when executed in a distributed manner underp partial asynchrony.

See Appendix D for a proof adapting results from Bertsekas and Tsitsiklis ().

###### Proposition 4.2 .

## 5 Energy GNNs

Under the assumptions of Section 2.1, implicit GNNs in which node embeddings are updated iteratively using local information are well suited for partially asynchronous, decentralized, and distributed inference.
However, as we discussed in Section 3, the diversity of existing implicit GNN architectures (both fixed-point and optimization-based) is limited.
Fixed-point GNNs employ simple message-passing schemes that can be enforced to be contractive, where $f_{\theta}$ consists of a linear function to compute messages ${\bm{m}}_{i}$, possibly followed by the application of a component-wise non-decreasing non-linearity to compute updated embeddings.
Optimization-based GNNs use a restricted class of graph objectives $E_{\theta}$ as shown in Equation 7.
One limitation with this objective is that it assumes informative node features exist, which may not be the case.
Edge features are not easily incorporated, which may be present and relevant for prediction.
Finally, smoothness across embeddings for neighboring nodes may be inappropriate for some datasets or tasks where neighboring nodes are expected to generate distinct embeddings.

We propose a novel implicitly-defined, optimization-based GNN architecture which we call the energy GNN.
Energy GNNs compute node embeddings that minimize a parameterized, convex graph function $E_{\theta}$, which we refer to as the ‘energy’ function.
In contrast to previous work on optimization-based GNNs, our energy function makes use of partially input-convex neural networks (PICNNs, Amos et al. ()).
PICNNs are scalar-valued neural networks that constrain the parameters in such a way that the network is convex with respect to a specified subset of the inputs.
This exposes a rich and flexible class of convex energy functions $E_{\theta}$ of the form:

$$ $\displaystyle E_{\theta}({\mathcal{G}},{\bm{H}})$ $\displaystyle=\textstyle\sum_{i=1}^{n}e^{i}_{\theta}\qquad\text{where}$ (14) $\displaystyle e^{i}_{\theta}$ $\displaystyle=u\left({\bm{m}}_{i},{\bm{h}}_{i},{\bm{x}}_{i};\theta_{u}\right)+ (\nicefrac{{\beta}}{{2}})||{\bm{h}}_{i}||_{2}^{2}\,.$ (15) $\displaystyle{\bm{m}}_{i}$ $\displaystyle=\textstyle\sum_{j\in\mathcal{N}(i)}m({\bm{h}}_{i},{\bm{h}}_{j},{ \bm{x}}_{i},{\bm{x}}_{j},{\bm{e}}_{ij};\theta_{m})$ (16) $$

Here ${\bm{m}}_{i}$ is the aggregation of state from neighboring nodes, where $m$ is a function that is both convex and nondecreasing (in each dimension) with respect to ${\bm{h}}_{j}$ and ${\bm{h}}_{i}$.
The function $u$ is convex with respect to ${\bm{m}}_{i}$ and ${\bm{h}}_{i}$.
These functions are both implemented as PICNNs and their composition is convex.
Summing these functions along with the squared norm penalty results in $E_{\theta}$ being strongly convex with respect to the node embeddings ${\bm{H}}$.
This architecture for the energy can be described as a (single layer) partially input-convex GNN; we provide more details in Appendix B.

This formulation for $E_{\theta}$ offers significantly more flexibility than the architectures of other implicitly-defined GNNs.
Any combination of inputs to $m$ or $u$ are valid, provided that $E_{\theta}$ remains a convex function of ${\bm{H}}$.
This means that edge features are easily incorporated into the model, and neighbor-specific or neighbor agnostic messages can be used (e.g., $m$ can take in information from just a node’s neighbor, or information pertaining to both a node and its neighbor); in our experiments in Section 6.3 we show that this translates empirically to improved performance on various tasks.
Since the aggregation in Equation 16 is only constrained to be a non-negative sum over neighbors, it can be replaced with, for example, a mean or a sum weighted by the entries of the symmetric renormalized adjacency matrix (as in GCNs, see Appendix A).
Alternatively, Equation 16 can easily incorporate a neighbor attention mechanism (as in GATs), where neighbor contributions to the sum are scaled by neighbor-specific attention weights (see Appendix B).
The attention weights can depend on any of the non-convex inputs (i.e., the features).

## 6 Experiments

### 6.1 Synthetic Multi-Agent Tasks

We perform experiments on a number of synthetic datasets, motivated by tasks which are of interest for multi-agent systems and where distributed, asynchronous inference is desirable.
We describe each task and associated dataset below.

#### Chains

The ability to communicate information across long distances in a group of agents is important when agent predictions depend on global information.
This communication is made more difficult in the absence of a central controller (as is the case for distributed, asynchronous inference).

The chains dataset, used in Gu et al. (); Liu et al. (), is meant to evaluate the ability to capture long-range dependencies between nodes.
The dataset consists of $p$ undirected linear graphs with $l$ nodes, with each graph having a label $k\in\{1,...,p\}$.
The task is node classification of the graph label, where class information is contained only in the feature of the first node in the chain; the node feature matrix ${\bm{X}}\in\mathbb{R}^{n\times p}$ for a graph with class $k$ has ${\bm{X}}_{1,k}=1$ and zeros at all other indices.
Perfect classification accuracy indicates that information is successfully propagated from the first node to the final node in the chain.
For our dataset, we use chains with length $l=100$.

#### Counting

Counting the number of agents in a group may be important in various multi-agent tasks.
This value can be used, for example, to calculate means, or in agent decisions which rely on group size.
We construct a dataset meant to evaluate the ability of GNNs to count in undirected chain graphs.
Our dataset consists of $50$ graphs with $1$-$50$ nodes.
Since no informative node features are present for this task, we set node features as one-hot embeddings of node degrees.
The prediction target for each node in a given graph is the total number of nodes in that graph.

#### Sums

For this task, we consider summation, a basic functional building block relevant for many multi-agent tasks.
For instance, in reinforcement learning tasks, agents might aim to perform actions that optimize their collective rather than individual rewards, requiring each agent to sum the rewards associated with all other agents.
Many distributed and asynchronous algorithms exist for summation (Kempe et al., ).

We construct a dataset to evaluate the ability of GNNs to perform binary sums in undirected chain graphs.
Our data are $2000$ graphs with $50$ nodes each, with different instantiations of binary node features ${\bm{x}}_{i}\in\{0,1\}$.
The prediction target for each node in a given graph is $\hat{{\bm{y}}}_{i}:=\sum_{i}{\bm{x}}_{i}$.

#### Coordinates

A common task for multi-agent collectives such as robot swarms is localization.
This problem has previously been tackled in various ways that all employ a bespoke algorithm tailored for the task (Todescato et al., ; Huang and Tian, , ).

We test the ability of GNNs to solve this problem on static graphs.
We construct a dataset where each node has a position in $\mathbb{R}^{2}$ and neighbors within some radius are connected by an edge.
We do not assume a global coordinate system; instead, we focus on relative localization, where pairwise distances between nodes are maintained.
Each node predicts a position in $\mathbb{R}^{2}$, and the objective is the mean squared error between true pairwise node distances, and distances between their predicted positions.
In order to break symmetries, each node has a unique ID which is one-hot encoded and used as the node feature. Distances to connected neighbors are provided as edge features.
We generate $1500$ random graphs where all graphs consist of $20$ nodes.
We sample uniformly in the unit square to get node positions and connect nodes by an edge if they are within a distance of 0.5.

#### MNIST “Terrain"

The final synthetic task we consider is “terrain” classification.
Suppose a number of agents are placed in an environment where each agent performs some local measurement, and the agents must collectively make predictions about some global state of the environment using only local communication.
For this experiment, we use MNIST images (those with 0/1 labels only) to represent the environment (LeCun et al., ).
Agents are placed at random locations in the image, and use the coordinates and pixel value at their location as their node features.
The prediction target for each agent is the image label.
We resize the images to $10\times 10$ pixels, and sample 10 random pixels for agent locations.
Nodes share an edge if they are within 5 pixels of each other.

Unlike the previous synthetic tasks, in which existing distributed algorithms can be applied, no bespoke algorithm exists for the MNIST terrain task.
This is precisely the type of problem which motivates the development of GNNs which are robust to asynchronous and distributed inference.

### 6.2 Experimental Setup

For each synthetic task, we compare performance of energy GNNs to other implicitly-defined GNNs we identified in Section 3.
We employ three energy GNN architecture variants; node-wise, where messages are constructed using information from individual nodes, edge-wise, where messages are constructed using information pertaining to both nodes on an edge, including edge features, and edge-wise energy GNN with neighborhood attention.
In terms of other implicitly-defined GNNs, we focus on the fixed-point GNN architecture IGNN defined by Gu et al. () which is described in Appendix A, and against GSDGNN, an optimization-based GNN which uses the objective from Equation 7.
Two fixed point GNN architectures identified in Section 3 are excluded: EIGNN (Liu et al., ) and the GNN proposed by Scarselli et al. ().
The former is excluded because the fixed point is solved for directly in the forward pass using global information rather than iteratively using local information; the latter is excluded because fixed point convergence may not be achieved.
In addition to implicitly-defined GNNs, we also compare against two common explicitly-defined GNN architectures; GCN (Kipf and Welling, ), and GAT (Veličković et al., ).

The cost of the forward and backward pass (in terms of computation and/or memory) for implicitly-defined GNNs is variable, depending on the number of iterations required for convergence.
We mitigate this cost during training in two ways.
Since convergence for both fixed-point GNNs and optimization-based GNNs is guaranteed (since for fixed-point GNNs, $F_{\theta}$ is contractive and for optimization-based GNNs, $E_{\theta}$ is convex) implicit differentiation can be used to obtain gradients of the task-specific loss function $\mathcal{L}$ with respect to parameters $\theta$.
This avoids unrolling the fixed-point iterations in the backward pass, and requires a fixed amount of computation and memory.
We derive the gradient in Appendix C.
Furthermore, since the solution of the forward pass is unique (i.e. not dependent on the initialization of ${\bm{H}}$), the number of iterations in the forward pass can be reduced by initializing ${\bm{H}}$ to be the solution from the previous epoch of training.
In our work we employ both of these strategies during training.

We additionally perform experiments on benchmark datasets MUTAG (Srinivasan et al., ), PROTEINS (Borgwardt et al., ), and PPI (Hamilton et al., ) for node and graph classification to evaluate energy GNNs as a synchronous GNN architecture.
Although our objective is not performance in the synchronous setting, we show that they are nevertheless competitive on each dataset.
Details are provided in Appendix G.

Additional training details for the synthetic experiments are in Appendix F.

### 6.3 Results

In Table 1, we report synchronous performance of each GNN architecture on the synthetic tasks.
For regression tasks (counting, sums, coordinates) task performance is calculated as the root mean squared error over the test dataset normalized by the root mean value of the test dataset prediction targets.
For classification tasks (chains, MNIST) task performance is calculated as the mean test dataset classification error.
Table 1 reports performance for each task, with mean and standard deviation taken across 10 dataset folds and 5 random parameter seeds.

The results of the synthetic experiments empirically demonstrate the superiority of our energy GNN architecture compared to other implicitly-defined GNNs.
The node-wise energy GNN architecture improves performance over IGNN and GSDGNN, which we attribute to the use of PICNNs.
When edge-wise rather than node-wise information is used in constructing messages to neighbors, further improvements in performance are observed.
The strong performance on tasks requiring long-distance communication between nodes for correct predictions (chains, counting, and sums) shows that our architecture is capable of capturing long-range dependencies between node predictions.

**Table 1: Task performance on test data, reported as percentage error (relative root mean squared error for COUNT, SUM, COORDINATES). Mean and standard deviation are across 10 random seeds and 5 train/test splits. Although inference on these experiments is done synchronously, the poor performance of the explicitly-defined GCN and GAT can be attributed to their depth-limited ability to propagate information.**
| MODEL | CHAINS | COUNT | SUM | COORDINATES | MNIST |
| --- | --- | --- | --- | --- | --- |
| IGNN | $26.9\pm 13.3$ | $40.2\pm 5.5$ | $12.8\pm 0.8$ | $52.0\pm 5.5$ | $30.4\pm 0.7$ |
| GSDGNN | $35.9\pm 6.3$ | $40.3\pm 0.6$ | $13.3\pm 0.2$ | $44.0\pm 1.0$ | $29.3\pm 0.6$ |
| Energy GNN node-wise | $15.8\pm 17.9$ | $19.5\pm 1.7$ | $12.2\pm 0.7$ | $41.3\pm 2.4$ | $13.8\pm 0.8$ |
| Energy GNN edge-wise | $1.2\pm 2.2$ | $4.0\pm 3.6$ | $\mathbf{4.9\pm 3.3}$ | $33.5\pm 3.2$ | $\mathbf{13.0\pm 0.6}$ |
| Energy GNN + attention | $\mathbf{0.25\pm 0.5}$ | $\mathbf{3.6\pm 3.8}$ | $6.0\pm 4.0$ | $\mathbf{30.9\pm 1.8}$ | $13.8\pm 0.8$ |
| GCN | $47.0\pm 0.0$ | $40.7\pm 1.0$ | $13.1\pm 0.3$ | $53.2\pm 0.9$ | $29.7\pm 0.5$ |
| GAT | $47.0\pm 0.0$ | $41.5\pm 0.8$ | $12.9\pm 0.5$ | $39.3\pm 1.0$ | $15.3\pm 3.2$ |

**Table 2: Decrease in task performance (decrease in accuracy for CHAINS, MNIST, and increase in relative RMSE for COUNT, SUM, COORDINATES) observed from switching from synchronous to asynchronous inference on sub-sample of test data (10 samples) using one trained model instance. Mean and standard deviation are across 5 asynchronous runs. The poor performance of GCN and GAT are consistent with the expected unreliability of explicitly-defined GNNs with asynchronous inference.**
| MODEL | CHAINS | COUNT | SUM | COORDINATES | MNIST |
| --- | --- | --- | --- | --- | --- |
| IGNN | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ |
| GSDGNN | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ |
| Energy GNN node-wise | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ |
| Energy GNN edge-wise | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ |
| Energy GNN + attention | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ | $0.0\pm 0.0$ |
| GCN | $38.8\pm 3.3$ | $584.6\pm 42.4$ | $2.6\pm 0.2$ | $63.4\pm 0.0$ | $37.4\pm 10.3$ |
| GAT | $6.6\pm 1.8$ | $250.3\pm 58.6$ | $45.1\pm 2.1$ | $97.0\pm 34.4$ | $50.4\pm 3.5$ |

In Table 2, we demonstrate empirically that the explicitly-defined architectures GCN and GAT perform poorly and unreliably under asynchrony, with task performance decreasing for all experiments.
In our experiments, we simulate asynchronous inference; our algorithm is in Appendix E.
In most cases the variance of the performance decrease is large as a result of inconsistent predictions under different random node update schedules.
In cases where the variance is low, we observe that predictions for different random schedules collapse to the same/similar values due to non-linearities in the architecture.
We reiterate that explicitly-defined GNNs are not suitable for the tasks we consider, in which we assume asynchronous, distributed inference and aim for reliable predictions under this regime.
For implicitly-defined GNNs, the decrease in performance under asynchronous inference is less than 0.1%, empirically confirming the convergence guarantees given by Proposition 4.1 and Proposition 4.2.

## 7 Related work

#### Asynchronous algorithms

Asynchronous algorithmic models (sometimes called chaotic relaxation models) date back at least to the 1960s (Chazan and Miranker, ), and were explored extensively into the 1970s and 1980s (Donnelly, ; Miellou, ; Robert et al., ; Baudet, ; Bertsekas, , ; Bojańczyk, ; Mitra, ; Üresin and Dubois, ).
In Section 4 we study GNNs under partial asynchrony, a model which imposes constraints on the sequencing of computations and frequency of communication between distributed elements.
Our analysis of partially asynchronous GNN inference draws directly from prior work which analyzes the sufficient conditions for convergence (Tsitsiklis, ; Tsitsiklis et al., ; Bertsekas and Tsitsiklis, ).

#### Asynchronous GNN training

There has been a plethora of work in distributed training of GNNs, where data are computed on by separate workers (Besta and Hoefler, ; Shao et al., ).
In this regime, workers may not be operating on independent data, e.g., when data come from a single connected graph.
Since ignoring this dependence reduces performance, workers exchange embedding information associated with “boundary vertices” which are logically connected but delegated to different workers.
Some distributed training assumes workers operate asynchronously, using stale
embedding (Md et al., ; Peng et al., ; Wan et al., ) or gradient information (Thorpe et al., ) corresponding to previous training epochs.
However, across all these frameworks, the forward pass is executed synchronously, and to our knowledge there is no work examining asynchronous GNN inference.

#### GNNs applied to multi-agent tasks

Several works have identified GNNs as a useful tool in parameterizing decentralized policies for multi-agent systems, with applications including flock formation, target tracking, path planning, and channel allocation in wireless networks.
These works employ GNNs to generate neighborhood-aggregated state representations for agent control policies, which are optimized through imitation learning (Grattarola et al., ; Zhou et al., ; Gama et al., ) or multi-agent reinforcement learning (Nakashima et al., ; Khan et al., ; Jiang and Lu, ; Jiang et al., ).
A pitfall of these works is that they do not consider asynchronous execution of the GNNs which generate state representations.
Consistent with the observations in our work, Blumenkamp et al. () show that performance of multi-agent policies that use GNNs suffer in real-world deployment where synchronous execution is not possible.

## 8 Conclusion

GNNs have the potential to provide learning frameworks for distributed systems, with applications to privacy, robotics, remote sensing, and other domains.
However, most conventional GNN architectures are not compatible with asynchronous inference and this hinders their deployment on these applications.

In this work, we characterize the class of implicitly-defined GNNs as being provably robust to partially asynchronous inference.
Motivated by the lack of diversity in this class of architectures, we contribute a novel addition in the form of energy GNNs.
We show that energy GNNs achieve better performance than other implicitly-defined GNNs on a number of synthetic tasks motivated by decentralized and constrained multi-agent systems, where distributed, asynchronous inference is desirable.
In addition to its robustness to asynchrony, our method is comparable in generalization performance (on benchmarks) with other modern GNN architectures that do not offer these guarantees.

The positive results of our synthetic experiments motivates additional work in applying GNN architectures to multi-agent related tasks.
A specific line of work which we expect to be interesting is real-time inference of dynamic graphs, both because of the relevance to problems in, e.g., robotics, but also due to the “anytime” nature of implicitly-defined GNN architectures.
Distributed and asynchronous inference over large graphs is another application for implicitly-defined GNNs which should be explored.
We are optimistic that the energy GNN framework itself, and implicitly-defined GNNs in general, will provide a path forward for “learning to learn” distributed, asynchronous algorithms.

There are several limitations of this work.
The forward pass in implicitly-defined GNNs is more involved computationally, and may require an amount of time that cannot be predetermined.
The condition number of the PICNN has a significant effect on convergence rate and is difficult to control.
The backward pass requires implicit gradients rather than simple backpropagation.
Our model for partial asynchrony, although standard, makes assumptions that may not match real-world systems.
Finally, operationalizing the theoretical results may require constants that are not readily available.

### Acknowledgments

This work was partially supported by NSF grants IIS-2007278 and OAC-2118201.

## Appendix A GNN architectures

#### Graph Convolutional Networks

GCNs
(Kipf and Welling, ) replace the adjacency matrix ${\bm{A}}$ with the symmetric normalized adjacency matrix with added self-loops, $\tilde{{\bm{A}}}=({\bm{D}}+{\bm{I}})^{-\frac{1}{2}}({\bm{A}}+{\bm{I}})({\bm{D}
}+{\bm{I}})^{-\frac{1}{2}}$. With node embeddings initialized to be equal to the node features, $f_{\theta}^{\ell}$ is defined as:

$$ $\displaystyle{\bm{m}}_{i}^{\ell}$ $\displaystyle:=\textstyle\sum_{j\in\mathcal{N}(i)}\tilde{{\bm{A}}}_{i,j}\theta _{m}^{\ell}{\bm{h}}_{j}^{\ell}$ $\displaystyle{\bm{h}}_{i}^{\ell+1}$ $\displaystyle:=\text{ReLU}({\bm{m}}_{i}^{\ell}),$ (17) $$

where ${\theta_{m}^{\ell}\in\mathbb{R}^{k^{(\ell)}\times k^{(\ell)}}}$.
This update can be succinctly described at the graph level as ${{\bm{H}}^{\ell+1}=\text{ReLU}(\tilde{{\bm{A}}}{\bm{H}}^{\ell}{\bm{W}}^{\ell})}$.
Note that for explicitly-defined message passing GNNs which have $L$ layers, such as GCN, it is impossible to propagate information farther than $L$ hops.

#### Graph Attention Networks

GATs
(Veličković et al., ) apply an attention mechanism to determine the weighting of information from different neighbors. With node embeddings initialized to be equal to the node features, $f_{\theta}^{\ell}$ is defined as:

$$ $\displaystyle\alpha_{i,j}$ $\displaystyle=\frac{\exp\left(\text{LeakyReLU}((\theta_{a}^{\ell})^{T}[\theta_ {m}^{\ell}{\bm{h}}_{i}||\theta_{m}^{\ell}{\bm{h}}_{j}])\right)}{\textstyle\sum _{k\in\mathcal{N}(i)}\exp\left(\text{LeakyReLU}((\theta_{a}^{\ell})^{T}[\theta _{m}^{\ell}{\bm{h}}_{i}||\theta_{m}^{\ell}{\bm{h}}_{j}])\right)}$ (18) $\displaystyle{\bm{m}}_{i}^{\ell}$ $\displaystyle:=\textstyle\sum_{j\in\mathcal{N}(i)}\alpha_{i,j}\theta_{m}^{\ell }{\bm{h}}_{j}^{\ell}$ (19) $\displaystyle{\bm{h}}_{i}^{\ell+1}$ $\displaystyle:=\text{ReLU}({\bm{m}}_{i}^{\ell}),$ (20) $$

where ${\theta_{m}^{\ell}\in\mathbb{R}^{q\times k^{(\ell)}}}$, and ${\theta_{a}^{\ell}\in\mathbb{R}^{2q}}$ are additional parameters used in aggregation.
Multiple attention heads can be used, in which there are $k$ aggregation functions with their own parameters. The messages generated by each of the attention heads are either concatenated or averaged to generate a single message ${\bm{m}}_{i}^{\ell}$.

#### Implicit Graph Neural Networks (IGNNs)

IGNNs Gu et al. () are a fixed point GNN architecture, in which each parameterized node update layer is contractive with respect to the node embeddings.
For simplicity, we assume a single parameterized layer $f_{\theta}$ and exclude the layer superscript from our notation.
The parameterized node embedding update is repeated for steps $t=0,...,T-1$, where the stopping point $T$ is determined by when node embeddings converge (within some numerical tolerance).
IGNNs use a similar embedding update function as GCN, but add node features as an additional input to the update function.
A layer $f_{\theta}$ is defined as:

$$ $\displaystyle{\bm{m}}_{i}(t)$ $\displaystyle:=\textstyle\sum_{j\in\mathcal{N}(i)}\tilde{{\bm{A}}}_{i,j}\theta _{m}{\bm{h}}_{j}(t)$ $\displaystyle{\bm{h}}_{i}(t+1)$ $\displaystyle:=\phi({\bm{m}}_{i}(t)+g({\bm{x}}_{i};\theta_{u})),$ (21) $$

where ${\theta_{m}\in\mathbb{R}^{k\times k}}$,  $g_{\theta}:\mathbb{R}^{n\times p}\rightarrow\mathbb{R}^{n\times k}$ and $\phi$ is a component-wise non-expansive function such as ReLU.
Convergence is guaranteed by constraining $||\theta_{m}||_{\infty}<\lambda_{max}(\tilde{{\bm{A}}})^{-1}$, where $\lambda_{max}$ is the maximum eigenvalue of $\tilde{{\bm{A}}}$. This ensures that the update is contractive, a sufficient condition for convergence. Since the fixed point is unique, ${\bm{h}}_{i}(0)$ can be initialized arbitrarily (although convergence time will vary).

#### Efficient Implicit Graph Neural Networks (EIGNNs)

EIGNNs Liu et al. () are another fixed point GNN architecture which are very similar to IGNNs.
The message passing iteration is constructed such that a closed-form solution can be obtained for the fixed point of a layer, which is more efficient than iterating message passing to convergence.
A layer $f_{\theta}$ is defined as:

$$ $\displaystyle{\bm{m}}_{i}(t)$ $\displaystyle:=\textstyle\sum_{j\in\mathcal{N}(i)}\gamma\alpha\tilde{{\bm{A}}} _{i,j}(\theta_{m})^{T}(\theta_{m}){\bm{h}}_{j}(t)$ $\displaystyle{\bm{h}}_{i}(t+1)$ $\displaystyle:={\bm{m}}_{i}(t)+{\bm{x}}_{i},$ (22) $$

where ${\theta_{m}\in\mathbb{R}^{k\times k}}$, $\alpha>0$ is a scaling factor equal to $\frac{1}{||(\theta_{m})^{T}(\theta_{m})||_{F}+\epsilon}$ with arbitrarily small $\epsilon$, and $\gamma\in(0,1]$ is an additional scaling factor.
The overall scaling factor $\gamma\alpha$ is chosen to ensure that the update is contractive, from which it follows that the sequence of iterates converges.

#### Nonlinear Fixed Point GNN

Scarselli et al. () introduce a general nonlinear fixed point GNN whose update can be written as

$$ $\displaystyle{\bm{m}}_{i}(t)$ $\displaystyle=\sum_{j\in\mathcal{N}(i)}m({\bm{h}}_{i}(t),{\bm{h}}_{j}(t),{\bm{ x}}_{i},{\bm{x}}_{j},{\bm{e}}_{ij};\theta_{m})$ (23) $\displaystyle{\bm{h}}_{i}(t+1)$ $\displaystyle=u({\bm{m}}_{i}(t),{\bm{h}}_{i}(t),{\bm{x}}_{i};\theta_{u}).$ (24) $$

where $m$ and $u$ are multi-layer neural networks.
As we discuss in Section 3, this flexible parameterization of message passing comes at a cost: it is difficult to enforce that the overall update is definitely contractive.

## Appendix B Input-convex GNN architecture details

As in Amos et al. (), we construct a parametric family of neural networks $f_{\theta}({\bm{x}},{\bm{y}})$ with inputs ${\bm{x}}\in\mathbb{R}^{n},{\bm{y}}\in\mathbb{R}^{m}$ which are convex with respect to ${\bm{y}}$ (i.e. a subset of the inputs).
We define a $k$-layer partially convex neural network by the recurrences:

$$ $\displaystyle{\bm{u}}_{i+1}$ $\displaystyle=\tilde{g}_{i}(\tilde{\bm{W}}_{i}{\bm{u}}_{i}+\tilde{\bm{b}}_{i})$ $\displaystyle{\bm{z}}_{i+1}$ $\displaystyle=g_{i}($ $\displaystyle{\bm{W}}_{i}^{(z)}({\bm{z}}_{i}\circ[{\bm{W}}_{i}^{(zu)}{\bm{u}}_ {i}+{\bm{b}}_{i}^{(z)}]_{+})+$ $\displaystyle{\bm{W}}_{i}^{(y)}({\bm{y}}\circ({\bm{W}}_{i}^{(yu)}{\bm{u}}_{i}+ {\bm{b}}_{i}^{(y)}))+$ $\displaystyle{\bm{W}}_{i}^{(u)}{\bm{u}}_{i}+{\bm{b}}_{i})$ $$

$$ $\displaystyle f_{\theta}({\bm{x}},{\bm{y}})={\bm{z}}_{k},\quad{\bm{u}}_{0}={ \bm{x}},\quad{\bm{z}}_{0}=\mathbf{0}$ $$

Provided the $W^{(z)}$ are elementwise nonnegative, and the activation functions $g_{i}$ are non-decreasing in each argument, it follows that $f_{\theta}$ is convex in ${\bm{y}}$.

In the context of an energy GNN, $E_{\theta}$ is comprised of two PICNNs, specialized to operate on graph structures; we refer to the resulting architecture as a partially input-convex GNN (PICGNN).
The message function $m$ in Equation 16 corresponds to a PICNN which has node embeddings ${\bm{h}}_{i},i=1,\dots,n$ as its convex inputs (all other features are non-convex inputs).
A second PICNN corresponds with the update function $u$ in Equation 15, which is convex in the node embeddings and the messages computed by the message function, and nonconvex in the node features.
The graph energy $E_{\theta}$ can then be written as the sum of the outputs of the second updating PICNN.

### B.1 Neighborhood Attention

We can incorporate a neighborhood attention mechanism in the PICGNN by modifying the aggregation of messages ${\bm{m}}_{i}$ as follows:

$$ $\displaystyle{\bm{m}}_{i}$ $\displaystyle=\textstyle\sum_{j\in\mathcal{N}(i)}\alpha_{i,j}m({\bm{h}}_{i},{ \bm{h}}_{j},{\bm{x}}_{i},{\bm{x}}_{j},{\bm{e}}_{ij};\theta_{m})\qquad\text{where}$ (25) $\displaystyle\alpha_{i,j}$ $\displaystyle=\frac{\exp\left(\text{LeakyReLU}(\theta_{a})^{T}[\theta_{m_{x}}{ \bm{x}}_{i}||\theta_{m_{x}}{\bm{x}}_{j}||\theta_{m_{e}}{\bm{e}}_{ij}])\right)} {\textstyle\sum_{k\in\mathcal{N}(i)}\exp\left(\text{LeakyReLU}((\theta_{a})^{T }[\theta_{m_{x}}{\bm{x}}_{i}||\theta_{m_{x}}{\bm{x}}_{j}||\theta_{m_{e}}{\bm{e }}_{ij}])\right)}$ (26) $$

where ${\theta_{m_{x}}\in\mathbb{R}^{q\times p}}$,
${\theta_{m_{e}}\in\mathbb{R}^{s\times r}}$ and ${\theta_{a}\in\mathbb{R}^{2q+s}}$ are additional parameters used in aggregation.
As in GATs, multiple attention heads can be used, in which there are $k$ aggregation functions with their own parameters. The messages generated by each of the attention heads are either concatenated or averaged to generate a single message ${\bm{m}}_{i}$.

## Appendix C Implicit Differentiation

Since we use an optimization procedure to compute the node embeddings within the forward pass in implicitly-defined GNNs, we need to obtain derivatives of the node embeddings with respect to the parameters of the model.
We compute derivatives by implicitly differentiating the optimality conditions.
For fixed-point GNNs, we use the fact that at the fixed point, we have parameters $\theta^{*}\in\mathbb{R}^{p}$ and node embeddings ${\bm{H}}^{*}\in\mathbb{R}^{n\times k}$ such that:

$$ $g({\bm{H}}^{*},\theta)=f_{\theta}({\bm{H}}^{*})-{\bm{H}}^{*}=\mathbf{0}.$ (27) $$

For optimization-based GNNs, we use the fact that at the solution of the minimization problem, we have:

$$ $g({\bm{H}}^{*},\theta)=\frac{\partial E_{\theta^{*}}}{\partial{\bm{H}}}({\bm{H }}^{*})=\mathbf{0}.$ (28) $$

Let $h^{*}(\theta^{*})={\bm{H}}^{*}$, so that we can write the optimality conditions in terms of the parameters only:

$$ $g(h^{*}(\theta^{*}),\theta^{*})=\mathbf{0}.$ (29) $$

Given an objective $\mathcal{L}:\mathbb{R}^{p}\mapsto\mathbb{R}$, the desired quantity is the total derivative of $\mathcal{L}$ with respect to the parameters. By the chain rule,

$$ $\displaystyle\frac{d\mathcal{L}}{d\theta}=\frac{\partial\mathcal{L}}{\partial h ^{*}}\frac{dh^{*}}{d\theta}+\frac{\partial\mathcal{L}}{\partial\theta}.$ (30) $$

We compute $\frac{\partial\mathcal{L}}{\partial h^{*}}$ and $\frac{\partial\mathcal{L}}{\partial\theta}$ using normal automatic differentiation, and the solution Jacobian $\frac{dh^{*}}{d\theta}$ using implicit differentiation.
Notice that, at the fixed point (where the optimality constraint is satisfied), we have:

$$ $\displaystyle\frac{d}{d\theta}g(h^{*}(\theta),\theta)$ $\displaystyle=\mathbf{0}$ (31) $\displaystyle\frac{\partial g}{\partial h^{*}}\frac{dh^{*}}{d\theta}+\frac{ \partial g}{\partial\theta}=\mathbf{0}$ (32) $\displaystyle\frac{\partial g}{\partial h^{*}}\frac{dh^{*}}{d\theta}$ $\displaystyle=-\frac{\partial g}{\partial\theta}.$ (33) $$

This is the primal (or tangent) system associated with the constraint function $g$.
In our setup we utilize reverse mode automatic differentiation, since $p\gg 1$ parameters are mapped to a single scalar objective.
Provided $\frac{\partial g}{\partial h^{*}}$ is invertible, we can rewrite the solution Jacobian as:

$$ $\frac{dh^{*}}{d\theta}=-(\frac{\partial g}{\partial h^{*}})^{-1}\frac{\partial g }{\partial\theta},$ (34) $$

and substitute this expression into Equation 30 as follows:

$$ $\frac{d\mathcal{L}}{d\theta}=-\frac{\partial\mathcal{L}}{\partial h^{*}}(\frac {\partial g}{\partial h^{*}})^{-1}\frac{\partial g}{\partial\theta}+\frac{ \partial\mathcal{L}}{\partial\theta}.$ (35) $$

For reverse mode, we compute the dual (or adjoint) of this equation,

$$ $\frac{d\mathcal{L}}{d\theta}^{T}=-\frac{\partial g}{\partial\theta}^{T}(\frac{ \partial g}{\partial h^{*}})^{-T}\frac{\partial\mathcal{L}}{\partial h^{*}}^{T }+\frac{\partial\mathcal{L}}{\partial\theta}^{T},$ (36) $$

And solve the dual system:

$$ $\frac{\partial g}{\partial h^{*}}^{T}\lambda=-\frac{\partial\mathcal{L}}{ \partial h^{*}}^{T}$ (37) $$

for the dual variable $\lambda$.

## Appendix D Proof of Convergence Under Partial Asynchrony

Our guarantee of convergence uses the classical result from Bertsekas and Tsitsiklis (, Chapter 7.5), which depends on several formal assumptions.
We reproduce these assumptions here, in the notation used in the present paper, and address how they are satisfied in our problem setting.
In the following, $s_{i}(t)$ is the “search direction” taken by node $i$, i.e., the vector used by node $i$ to take optimization step.
Ideally, this would be the negative gradient $-\nabla_{{\bm{h}}_{i}}E_{\theta}$, but partial asynchrony means it may be a vector constructed from stale information.
Our goal is to show that convergence of the optimization is nevertheless guaranteed.
We use ${\bm{h}}\in\mathbb{R}^{nk}$ to denote the unrolled embeddings ${\bm{H}}$.

###### Assumption D.1 ( Bertsekas and Tsitsiklis ( 1989 ) Assumption 5.1).

- (a)
There holds $E_{\theta}({\bm{h}})\geq 0$ for every ${\bm{h}}\in\mathbb{R}^{nk}$.
- (b)
(Lipschitz Continuity of $\nabla E_{\theta}$) The function $E_{\theta}$ is continuously differentiable and there exists a constant $K_{1}$ such that
$||\nabla E_{\theta}({\bm{h}})-\nabla E_{\theta}({\bm{h}}^{\prime})||\leq K_{1}
||{\bm{h}}-{\bm{h}}^{\prime}||,\quad\forall{\bm{h}},{\bm{h}}^{\prime}\in
\mathbb{R}^{nk}$

For part *(a)*, since both convexity and absolute continuity are preserved under nonnegative summation, the graph energy $E$ given in Equation 4 is a smooth, strictly convex function.
Without loss of generality, we can assume the node energies $e^{i}_{\theta}$ described in Equation 4 satisfy $e^{i}_{\theta}({\bm{h}}_{i})\geq 0$ for all ${\bm{h}}_{i}\in\mathbb{R}^{d_{i}}$.
This follows because the $e^{i}_{\theta}$ are strictly convex, therefore the optimal value $p_{i}^{*}=\inf\{e^{i}_{\theta}({\bm{h}}_{i})\}$ is achieved and thus $\tilde{e}_{i}=e^{i}_{\theta}+p_{i}^{*}$ is nonnegative.
Thus, we have that the graph energy is the sum of nonnegative terms: $E_{\theta}({\bm{h}}_{1},{\bm{h}}_{2},\dots,{\bm{h}}_{n})\geq 0$ for all $({\bm{h}}_{1},{\bm{h}}_{2},\dots,{\bm{h}}_{n})\in\mathbb{R}^{nk}$.

For part *(b)*, the assumption of the $e^{i}_{\theta}$ having a bounded Hessian implies that their sum also has a bounded Hessian, which further implies Lipschitz continuity of the gradient of $E_{\theta}$.

###### Assumption D.2 ( Bertsekas and Tsitsiklis ( 1989 ) Assumption 5.5).

- (a)
(Block-Descent) There holds $s_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))\leq-||s_{i}(t)||^{
2}/K_{3}$ for all $i$ and all $t\in T^{i}$.
- (b)
There holds $||s_{i}(t)||\geq K_{2}||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||$ for all $i$ and all $t\in T^{i}$.

For part *(a)*, the situation is slightly more complex than the conventional optimization setup, as the gradients for the search direction $s_{i}(t)$ are being computed from potentially-outdated neighbor embeddings, rather than it being the gradients themselves that are outdated.
Writing the negative search direction for node $i$ at time $t$ in terms of the update times $\tau(t)$, we have

$$ $\displaystyle\bar{s}_{i}(t)$ $\displaystyle:=-s_{i}(t)$ (38) $\displaystyle=\nabla_{{\bm{h}}_{i}}e^{i}_{\theta}({\bm{h}}_{1}(\tau^{i}_{1}(t) ),\ldots,{\bm{h}}_{n}(\tau^{i}_{n}(t)))+\!\!\!\sum_{j\in\mathcal{N}(i)}\!\!\! \nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}_{1}(\tau^{j}_{1}(\tau^{i}_{j}(t)) ),\ldots,{\bm{h}}_{n}(\tau^{j}_{n}(\tau^{i}_{j}(t)))),$ (39) $$

where the the gradient communicated from node $j$ may have used stale versions of the embedding both for node $i$ itself and for other nodes connected to $j$.
Contrast this with the “true” gradient computed at $i$ which would be computed from its current estimate of the complete state of the graph:

$$ $\displaystyle\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))$ $\displaystyle=\nabla_{{\bm{h}}_{i}}e^{i}_{\theta}({\bm{h}}_{1}(\tau^{i}_{1}(t) ),\ldots,{\bm{h}}_{n}(\tau^{i}_{n}(t)))+\!\!\!\sum_{j\in\mathcal{N}(i)}\nabla_ {{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}_{1}(\tau^{i}_{1}(t)),\ldots,{\bm{h}}_{n} (\tau^{i}_{n}(t)))\,.$ (40) $$

We introduce the following notation to simplify the exposition:

$$ $\displaystyle g^{k/k^{\prime}}_{j/i}(t)$ $\displaystyle:=\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}_{1}(\tau^{k}_{1}( \tau^{k^{\prime}}_{j}(t))),\ldots{\bm{h}}_{n}(\tau^{k}_{n}(\tau^{k^{\prime}}_{ j}(t))))\,,$ (41) $$

which can be read as “gradient of $e^{j}_{\theta}$ with respect to ${\bm{h}}_{i}$ from the perspective of node $k$ at the time corresponding to some node $k^{\prime}$’s view of node $j$ at time $t$.”.
With this notation, we define

$$ $\displaystyle\bar{s}_{i}(t)$ $\displaystyle=g^{i/i}_{i/i}(t)+\sum_{j\in\mathcal{N}(i)}g^{j/i}_{j/i}(t)$ $\displaystyle\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))$ $\displaystyle=g^{i/i}_{i/i}(t)+\sum_{j\in\mathcal{N}(i)}g^{i/j}_{j/i}(t)\,.$ (42) $$

We wish to show that the inner product between $\bar{s}_{i}(t)$ and $\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))$ is greater than $||s_{i}(t)||^{2}/K_{3}$, for some $K_{3}>0$ and all $i$ and $t$.

###### Lemma D.3 .

There exists an $\alpha>0$ such that $\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))\geq||s_{i}(t
)||^{2}/K_{3}$.

###### Proof.

Starting with the squared error in the negative search direction:

$$ $\displaystyle||\bar{s}_{i}(t)-\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{ 2}_{2}$ $\displaystyle=(\bar{s}_{i}(t)-\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t)))^{ \top}(\bar{s}_{i}(t)-\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t)))$ (43) $\displaystyle=||\bar{s}_{i}(t)||^{2}_{2}+||\nabla_{{\bm{h}}_{i}}E_{\theta}({ \bm{h}}(t))||^{2}_{2}-2\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({ \bm{h}}(t))$ (44) $$

we find an expression for the inner product:

$$ $\displaystyle\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))$ $\displaystyle=\frac{1}{2}\left(||\bar{s}_{i}(t)||^{2}_{2}+||\nabla_{{\bm{h}}_{ i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}-||\bar{s}_{i}(t)-\nabla_{{\bm{h}}_{i}}E_{ \theta}({\bm{h}}(t))||^{2}_{2}\right)\,.$ (45) $$

and so we require the following to be greater than or equal to zero:

$$ $\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))-||\bar{s}_{i }(t)||^{2}/K_{3}\\ =\frac{1}{2}\left((1\!-\!\frac{2}{K_{3}})||\bar{s}_{i}(t)||^{2}_{2}+||\nabla_{ {\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}-||\bar{s}_{i}(t)-\nabla_{{\bm{h }}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}\right)$ (46) $$

We can use the triangle inequality to find an upper bound on the term being subtracted:

$$ $\displaystyle||\bar{s}_{i}(t)-\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||_{2}$ $\displaystyle=\left\|\left(g^{i/i}_{i/i}(t)+\sum_{j\in\mathcal{N}(i)}g^{j/i}_{ j/i}(t)\right)-\left(g^{i/i}_{i/i}(t)+\sum_{j\in\mathcal{N}(i)}g^{i/j}_{j/i}(t )\right)\right\|_{2}$ (47) $\displaystyle=\left\|\sum_{j\in\mathcal{N}(i)}g^{j/i}_{j/i}(t)-g^{i/j}_{j/i}(t )\right\|_{2}$ (48) $\displaystyle\leq\sum_{j\in\mathcal{N}(i)}||g^{j/i}_{j/i}(t)-g^{i/j}_{j/i}(t)| |_{2}$ (49) $\displaystyle=\sum_{j\in\mathcal{N}(i)}||\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({ \bm{h}}_{1}(\tau^{j}_{1}(\tau^{i}_{j}(t))),\ldots,{\bm{h}}_{n}(\tau^{j}_{n}( \tau^{i}_{j}(t))))-\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}(\tau^{i}(t)))| |_{2}\,.$ (50) $$

Any difference between the states ${\bm{h}}(\tau^{j}(\tau^{i}_{j}(t)))$ and ${\bm{h}}(\tau^{i}(t))$ would arise because node $i$ and node $j$ observe different staleness states of one or more of their shared neighbors $j^{\prime}$; these different staleness states correspond to differences in the number of gradient steps taken by shared neighbors $j^{\prime}$, as observed by node $i$ and $j$.
Note that we can assume that node $i$ and $j$ agree on the values of neighbors $j^{\prime}$ of either node $i$ or $j$ which are not shared between them.
Assuming the norm of the gradient is bounded, i.e., $||\nabla_{{\bm{h}}_{j^{\prime}}}e^{k}_{\theta}||_{2}\leq K_{0}$, and the number of neighbors a node has is bounded, i.e., $|\mathcal{N}(j^{\prime})|\leq n_{\text{max}}$, then the staleness bound $B$ and the step size $\alpha$ imply

$$ $\displaystyle||{\bm{h}}(\tau^{j}(\tau^{i}_{j}(t)))-{\bm{h}}(\tau^{i}(t))||_{2}$ $\displaystyle\leq\alpha B_{0}K_{0}\qquad\text{where}\quad B_{0}:=2n_{\text{max }}|\mathcal{N}(i)\cap\mathcal{N}(j)|B\,,$ (51) $$

where we assume $\mathcal{N}(k)$ includes node $k$ itself.
The constant $2$ in the inequality above arises because the staleness of the gradient received by node $i$ from node $j$ is outdated by at most $B$ time units, and node $j$’s view of its neighbor embeddings is also outdated by at most $B$ time units; this means the staleness of the embeddings in the gradient received by node $i$ is stale by at most $2B$ time units.
The Lipschitz continuity condition then gives

$$ $\displaystyle||\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}(\tau^{j}(\tau^{i}_ {j}(t))))-\nabla_{{\bm{h}}_{i}}e^{j}_{\theta}({\bm{h}}(\tau^{i}(t)))||_{2}$ $\displaystyle\leq\alpha B_{0}K_{0}K_{1}$ (52) $$

and therefore

$$ $\displaystyle\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t)) -||\bar{s}_{i}(t)||^{2}/K_{3}$ (53) $\displaystyle\qquad=\frac{1}{2}\left((1-\frac{2}{K_{3}})||\bar{s}_{i}(t)||^{2} _{2}+||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}-||\bar{s}_{i}(t) -\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}\right)$ (54) $\displaystyle\qquad\geq\frac{1}{2}\left((1-\frac{2}{K_{3}})||\bar{s}_{i}(t)||^ {2}_{2}+||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}-(\alpha| \mathcal{N}(i)|B_{0}K_{0}K_{1})^{2}\right)\,.$ (55) $$

We can therefore satisfy the assumption by choosing $\alpha>0$ such that

$$ $\displaystyle\frac{(1-\frac{2}{K_{3}})||s_{i}(t)||^{2}_{2}+||\nabla_{{\bm{h}}_ {i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}}{(|\mathcal{N}(i)|B_{0}K_{0}K_{1})^{2}} \geq\alpha^{2}\,.$ (56) $$

For part *(b)* we require there to exist a lower bound on the magnitude of $s_{i}(t)$ relative to the magnitude of the true gradient; this prevents the step from being too small.

###### Lemma D.4 .

There exists an $\alpha>0$ such that $\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))\geq K_{2}||
\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}$.

###### Proof.

We can use a nearly identical argument to that done above in part *(a)*, but instead of Equation 53 we write

$$ $\displaystyle\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t)) -K_{2}||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}$ (57) $\displaystyle\qquad=\frac{1}{2}\left(||\bar{s}_{i}(t)||^{2}_{2}+(1-2K_{2})|| \nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}-||\bar{s}_{i}(t)-\nabla _{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}\right)$ (58) $\displaystyle\qquad\geq\frac{1}{2}\left(||\bar{s}_{i}(t)||^{2}_{2}+(1-2K_{2})| |\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}_{2}-(\alpha|\mathcal{N}(i) |B_{0}K_{0}K_{1})^{2}\right)\,.$ (59) $$

We can then choose $\alpha$ to be

$$ $\displaystyle\frac{||s_{i}(t)||^{2}_{2}+(1-2K_{2})||\nabla_{{\bm{h}}_{i}}E_{ \theta}({\bm{h}}(t))||^{2}_{2}}{(|\mathcal{N}(i)|B_{0}K_{0}K_{1})^{2}}\geq \alpha^{2}\,.$ (60) $$

###### Lemma D.5 .

If $\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))\geq K_{2}||
\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}$ then $||s_{i}(t)||\geq K_{2}||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||$.

###### Proof.

Noting that $||s_{i}(t)||=||\bar{s}_{i}(t)||$, we have

$$ $\displaystyle||s_{i}(t)||\cdot||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))|| \geq\bar{s}_{i}(t)^{\top}\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))\geq K_{2 }||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||^{2}\,.$ (61) $$

Dividing both the left and right sides by $||\nabla_{{\bm{h}}_{i}}E_{\theta}({\bm{h}}(t))||$ gives the desired result.

Having satisfied the assumptions, we can now apply the result that guarantees convergence.

###### Proposition D.1 (Bertsekas and Tsitsiklis (1989), Proposition 5.1).

Under Assumptions 2.1, D.1, and D.2, there exists some $\alpha_{0}>0$ (depending on $n$, $B$, $K_{1}$, and $K_{3}$) such that if $0<\alpha<\alpha_{0}$ then $\lim_{t\to\infty}\nabla E_{\theta}({\bm{h}}(t))=0$.

###### Assumption D.1 ( Bertsekas and Tsitsiklis ( 1989 ) Assumption 5.1).

###### Assumption D.2 ( Bertsekas and Tsitsiklis ( 1989 ) Assumption 5.5).

###### Lemma D.3 .

###### Proof.

###### Lemma D.4 .

###### Proof.

###### Lemma D.5 .

###### Proof.

###### Proposition D.1 (Bertsekas and Tsitsiklis (1989), Proposition 5.1).

## Appendix E Asynchronous GNN implementation

In our asynchronous inference experiments, we simulate partially asynchronous execution (see Algorithm 1). We fix the maximum staleness bound to $B=7$.
Staleness is achieved by staggering node updates across time and introducing delay in messages sent.
In particular, when a node updates, its next update time is selected randomly between time $t+1$ and $t+S$, where $S$ is the “stagger” time, and one of the last $D$ values of its neighbors is chosen for performing the node update.
This satisfies assumptions 1 and 2 of partial asynchronism; nodes update at least every $S$ time units, and a node’s view of its neighbors is stale by at most $B=S+D$ time steps.

Figure: Algorithm 1 Simulated asynchronous GNN inference

## Appendix F Experiment Details (Synthetic Experiments)

### F.1 Architecture details

For all implicitly-defined GNN architectures, we use the same node embedding size with ${\bm{h}}_{i}\in\mathbb{R}^{2}$.
The architectures are chosen such that the number of parameters is approximately equal between models (with the constraint of using the same node embedding size).
All architectures employ an output function $o_{\phi}$ which is parameterized as an MLP with layers $(4,4,1)$ ($(4,4,2)$ for the coordinates experiment, as node predictions are positions in $\mathbb{R}^{2}$).

#### Energy GNN

For energy GNN, we use a PICNN with layer sizes $(4,4,2)$ for the message function $m$ in Equation 16 and a PICNN with layer sizes $(4,4,1)$ for the update function $u$ in Equation 15.
The aggregation in Equation 16 uses the entries of the unnormalized adjacency matrix ${\bm{A}}$ with no self-loops added.
We add an independently parameterized self-loop to the message passing function $m$ (instead of using the same parameters as for neighbors).
We set $\beta=0.04$.

#### Energy GNN + Attention

For energy GNN with the attention mechanism, the attention weights $\alpha_{ij}$ are computed as in Equation 18; however, to maintain convexity with respect to the embeddings, the node and edge features (if present) are used in rather than using the embeddings ${\bm{h}}$.
We use 2 attention heads, where the attention head outputs are concatenated to form the message ${\bm{m}}_{i}$ in Equation 16.

#### GSDGNN

For GSDGNN (where Equation 7 describes the optimization objective for obtaining embeddings), we parameterize $g_{\theta}$ as an MLP with layer sizes $(16,16,16,2)$.
We use the symmetric renormalized Laplacian matrix $\tilde{{\bm{L}}}={\bm{I}}-\tilde{{\bm{A}}}={\bm{I}}-({\bm{D}}+{\bm{I}})^{-
\frac{1}{2}}({\bm{A}}+{\bm{I}})({\bm{D}}+{\bm{I}})^{-\frac{1}{2}}$ for the Laplacian regularization term, and set $\gamma=1.0,\beta=5.0$.
(Note that gradient-based optimization of this objective function has a direct correspondence to the embedding update function of APPNP Gasteiger et al. ().)

#### IGNN

For IGNN (described in Equation 21), we parameterize $g_{\theta}:\mathbb{R}^{n\times p}\rightarrow\mathbb{R}^{n\times k}$ as an MLP with layer sizes $(16,16,16,2)$.

#### GCN and GAT

For GCN, we use $5$ layers of message passing with layer sizes $(10,10,10,10,10)$.
For GAT, we use $5$ layers of message passing with layer sizes $(3,3,3,3,3)$, and concatenate the output of $3$ attention heads at each layer.

### F.2 Training details

For binary classification experiments, we use binary cross entropy loss for training.
For regression experiments, we use mean squared error.
We use the Adam optimizer with weight decay, where we set the optimizer parameters as $\alpha=0.001,\beta_{1}=0.9,\beta_{2}=0.999$.
We set the learning rate to $0.002$, and use exponential decay with rate $0.98$ ever $200$ epochs.
We train for a maximum of 5000 epochs.
These experiments were performed on a single NVIDIA RTX 2080 Ti.

## Appendix G Experiment Details (Benchmark Datasets)

### G.1 Dataset Details

The benchmark datasets we report performance for are MUTAG and PROTEINS, where the prediction task is graph classification, and PPI, where the prediction task is node classification.

#### MUTAG

MUTAG is a dataset consisting of 188 graphs, each of which corresponds to a nitroaromatic compound (Srinivasan et al., ). The goal is to predict the mutagenicity of each compound on Salmonella typhimurium. Nodes in the graphs correspond to atoms (and are associated with a one-hot encoded feature in $\mathbb{R}^{7}$ corresponding to the atom type), and edges correspond to bonds. The average number of nodes in a graph is 17.93, and the average number of edges is 19.79.

#### PROTEINS

The PROTEINS dataset Borgwardt et al. () consists of 1113 graphs, each of which corresponds to a protein. The task is predicting whether or not the protein is an enzyme. Nodes in the graph correspond to amino acids in the protein (and are associated with node features in $\mathbb{R}^{3}$ representing amino acid properties). Edges connect amino acids that are less than some threshold distance from one another in the protein. The average number of nodes in a graph is 39.06, and the average number of edges is 72.82.

#### PPI

The PPI dataset (Hamilton, ) consists of 24 graphs, each of which corresponds to a protein-protein interaction network found in different areas of the body. Each node in the graph corresponds to a protein, with edges connecting proteins that interact with one another. Nodes are associated with features in $\mathbb{R}^{50}$, representing some properties of the protein. Each protein has $121$ binary prediction targets, each of which corresponds to some ontological property that the protein may or may not have. We use a 20/2/2 train/valid/test split consistent with Hamilton et al. ().

### G.2 Results

For all experiments with benchmark datasets, we use the same training procedure and architectures as described in Appendix F.
For node classification tasks, the final layer of output function $o_{\phi}$ is modified to use layer sizes $(4,4,\text{(num_classes)})$ where num_classes is the number of distinct class labels.
For graph classification tasks, we obtain graph-level predictions by passing the mean of the node-level predictions through a graph readout function parameterized by an MLP with layers $(4,4,\text{(num_classes)})$.

For PROTEINS and MUTAG, We perform 10-fold cross validation and report average classification accuracy and standard deviations in Table 3. We include performance reported by Gu et al. () (marked by an asterisk), as well as performance of our own implementation of IGNN with a single layer and GCN with two layers.
For PPI, we use a 20/2/2 train/valid/test split consistent with Hamilton et al. ().

Table 4 shows average micro-f1 scores for energy GNNs compared to other GNN architectures. We include performance reported by Gu et al. () (marked by an asterisk).
Where layer specifications are not included (for asterisked values), we were unable to determine them from Gu et al. ().

**Table 3: Graph classification accuracy (%). Results are averaged (and standard deviations are computed) using 10 fold cross validation with 5 random parameter seeds. Asterisked values are from Gu et al. ().**
|  | DATASET |  |
| --- | --- | --- |
| MODEL | MUTAG | PROTEINS |
| GCN* (5 layer) | 85.6 $\pm$ 5.8 | 76.0 $\pm$ 3.2 |
| IGNN* (3 layer) | 89.3 $\pm$ 6.7 | 77.7 $\pm$ 3.4 |
| Energy GNN | 87.6 $\pm$ 3.4 | 72.5 $\pm$ 0.3 |
| Energy GNN + attention | 79.5 $\pm$ 1.8 | 72.5 $\pm$ 0.5 |
| IGNN (1 layer) | 73.0 $\pm$ 0.6 | 71.8 $\pm$ 1.3 |
| GSD GNN | 78.4 $\pm$ 2.2 | 72.8 $\pm$ 0.8 |
| GCN (5 layer) | 78.0 $\pm$ 1.4 | 73.7 $\pm$ 0.5 |
| GAT (5 layer) | 76.1 $\pm$ 1.5 | 71.7 $\pm$ 3.2 |

**Table 4: Mean micro-F1 score on PPI dataset (%). Asterisked values are from Gu et al. ().**
| MODEL | micro f1 |
| --- | --- |
| MLP* | 46.2 |
| GCN* | 59.2 |
| GraphSAGE* | 78.6 |
| GAT* (3 layer) | 97.3 |
| IGNN* (5 layer) | 97.6 |
| Energy GNN | 76.2 |
| Energy GNN + attention | 76.0 |
| IGNN (1 layer) | 75.5 |
| GSD GNN | 76.0 |
| GAT (5 layer) | 74.3 |
| GCN (5 layer) | 76.2 |

## Peer Reviews

confidence: 2
contribution: 3
initial_rating: 6
presentation: 2
questions: - Could you provide a more intuitive explanation of implicitly defined GNNs, specifically highlighting the mechanisms contributing to their robustness under asynchronous updates?
- Do you have any insight into why the performance on real-world datasets is not as great as on synthetic data?
- Can you offer an explanation or hypothesis as to why explicitly defined GNNs did not perform as poorly as expected on the PPI dataset under asynchronous conditions?

Line 369. There appears to be a typo in the index notation within the description of targets for the Sums dataset.
soundness: 3
strengths: - Although I am not deeply familiar with the literature on asynchronous inference in GNNs, the proposed Energy GNN model introduces a novel and meaningful contribution to this field.
- The paper is of high quality. The authors provide comprehensive mathematical proofs for the convergence of Energy GNNs. They also provide a detailed description of the experimental setup and results, which are well-organized and easy to follow.
- The paper is generally well-structured, with a clear definition of the problem space and a concise summary of related GNN methods. The paper clearly defines the problem of asynchronous inference in explicitly-defined GNNs.
- The proposed Energy GNN model demonstrates significant improvements on synthetic datasets, and competitive performance on real-world datasets.
summary: This paper investigates the performance of Graph Neural Networks (GNNs) in partially asynchronous inference settings, where nodes update in a staggered or asynchronous manner. The authors categorize GNNs into two types: "explicitly-defined" and "implicitly-defined." They demonstrate that explicitly-defined GNNs are highly vulnerable to asynchronous updates. In contrast, implicitly-defined GNNs are shown to be robust. Additionally, the authors propose a novel explicitly-defined model, termed Energy GNN, which achieves notable improvements over existing methods on synthetic datasets.
weaknesses: The paper lacks a clear and intuitive explanation of implicitly-defined GNNs, which is essential for understanding their robustness to asynchronous updates. While the authors offer detailed explanations for explicitly-defined GNNs, which are more straightforward, they do not provide the same depth of insight into implicitly-defined GNNs. This makes it difficult for readers unfamiliar with the topic to understand how implicitly-defined GNNs work and why they are resilient to asynchronous inference.

Additionally, although the proposed Energy GNN shows strong results on synthetic datasets, its performance on real-world datasets is rather competitive. On the PPI dataset, in particular, its performance is comparable to that of explicitly-defined GNNs, which were expected to fail under asynchrony. This discrepancy between synthetic and real dataset performance is not explained. A broader evaluation across various real-world datasets would increase the credibility of Energy GNNs as a robust solution.

---

confidence: 4
contribution: 2
initial_rating: 6
presentation: 4
questions: Stemming from the weaknesses:

**Q1 (W1)**. Is there a different way to evaluate the differences among implicit GNNs other than on sync tasks?  
**Q2 (W1)**. How strong is EnergyGNN compared to more long-range optimized GNNs like Half-Hop, DRew, and graph transformers?  
**Q3 (W2)**. Does the synthetic EnergyGNN performance correlate with the tasks from more real-world benchmarks like LRGB and TGB?
soundness: 3
strengths: **S1**. A theoretical study on the async inference with GNNs is timely and important - many real-world tasks are of that nature, so having a principled, robust approach for such problems (instead of tinkering standard sync GNN architectures) might be of interest to the graph learning community.

**S2**. The paper is well-written - complex concepts are properly introduced and explained (which is often a challenge in the literature on implicit GNNs), the story and motivation are easy to follow.
summary: The paper studies GNNs in the asynchronous “hogwild” setup – when node states are not updated simultaneously at the same time. Such setup is common in distributed environments, agentic systems, and temporal learning where enforcing synchronous node updates might be impossible or too expensive to maintain. Standard, explicitly-defined GNNs trained in the synchronous mode fail in the async mode, so the authors turn their attention to implicitly-defined GNNs (further categorized into fixed-point and optimization-based GNNs). 

It is shown that implicit GNNs are robust to partial asynchronicity. Motivated by the theoretical findings, the authors propose EnergyGNN - an optimization-based implicit GNN where node states minimize a convex energy function. The space of possible energy functions is rather wide and allows for node-wise, edge-wise, and attention-based parameterizations. Experimenting on a range of synthetic tasks, the proposed EnergyGNN outperforms other implicit GNNs in the synchronous regime, is robust to the delayed node update setup, and is on par with synchronous baselines on MUTAG and Proteins datasets.
weaknesses: The main problem of the work is in the experiments - it is hard to judge the claimed effectiveness of the proposed EnergyGNN using only synthetic experiments and basic GCN / GAT as baselines.

**W1**. In the proposed suite of tasks, all implicit GNNs are robust in the async setup (the main goal of the work), and the main difference lies in the performance in the sync setup. Is there a different way to evaluate the differences among implicit GNNs other than on sync tasks? EnergyGNN is better than other implicit models but comparing against vanilla GCN and GAT on benchmarks defined to be of a long-range nature (where vanilla GNNs are bound to fail) is questionable. There is a variety of explicitly-defined GNNs that might be stronger baselines in such setups like Half-Hop [1], DRew [2] and other graph rewiring methods, as well as various graph transformers from the dedicated Long-Range Graph Benchmark [3]. 

**W2**. The proposed synthetic benchmarks are rather small and might not correlate with the performance on real datasets where async inference is important (or MUTAG and Proteins with sync inference); it seems to be a stretch to attribute chains, node counting, and node sums to “agentic” tasks. Instead, experiments on more real benchmarks might be more informative and evidential: 
* Since some tasks focus on the long-range dependencies, LRGB [3] is a suitable choice;
* As async inference implies nodes appearing and disappearing at some moments of time, temporal graph benchmarks with many snapshots [4] might be a great choice.

Generally, I am willing to increase the score if the authors add more modern baselines and/or real-world datasets.

[1] Azabou et al. Half-Hop: A graph upsampling approach for slowing down message passing. ICML 2023.   
[2] Gutteridge et al. DRew: Dynamically Rewired Message Passing with Delay. ICML 2023.  
[3] Dwivedi et al. Long Range Graph Benchmark. NeurIPS 2022.   
[4] Huang et al. Temporal Graph Benchmark for Machine Learning on Temporal Graphs. NeurIPS 2023.

---

confidence: 4
contribution: 2
initial_rating: 5
presentation: 2
questions: Please consider changing the term energy GNN. My first reaction was that this solution was energy-efficient which is a very big concern in AI at the moment. However this work is not about energy efficient AI. 

The experiments are on toy examples that do not need machine learning. Please expand the experimental range to other problems and to larger problem scale. 

Can you be clearer about what you summarize from existing work and what are the new contributions in sections 2, 3, 4 

The description of energy GNN is very sparse, please explain the computation of the neuron and node, architecture, and the method for training in the main manuscript. Please address summarize the properties of this model in Section 5. 

Can you say anything about the performance and energy consumption of these models?
soundness: 3
strengths: By focusing on implicitly-defined GNNs, the work addresses a major limitation of conventional GNNs in handling asynchronous and unreliable communication, making it highly relevant for real-world multi-agent systems.

Experimental results show that energy GNNs outperform other implicitly-defined GNNs on synthetic tasks, providing empirical validation for the architecture's effectiveness in multi-agent tasks (although most experiments are toy-ish.)
summary: This paper introduces Graph Neural Networks (GNNs) for distributes multi-agent systems with asynchronous communication. Traditional GNNs struggle with asynchronous execution and unreliable communication. This results in unreliable predictions. To address this, the authors focus on implicitly-defined GNNs, a class of models that can handle partial asynchrony. They introduce a model called energy GNNs, which outperform existing implicitly-defined GNNs on synthetic multi-agent tasks.

The experimental results highlight potential applications for GNNs in control tasks and real-time inference on dynamic graphs. This is an important property for multi-robot systems. The paper also notes the training limitations for implicitly-defined GNNs, particularly in achieving convergence, as it requires complex and often unpredictable computations. Strategies like warm-start initialization and implicit differentiation could help with these challenges, but challenges remain in achieving stability.
weaknesses: The experiments are conducted on toy examples. 

For the experiments other than the "terrain" examples there are great solutions that do not require machine learning. The results on the benchmark datasets in the supplementary show small improvements as compared to the more toy examples in the main manuscript.

---

confidence: 3
contribution: 3
initial_rating: 6
presentation: 3
questions: 1) How does the choice of the staleness bound B and stagger time S affect the performance and convergence of energy GNNs in practice? Is there a principled way to select these parameters?
2) The authors mention the potential for numerical instability due to ill-conditioned Hessians or Jacobians. Have they explored techniques like preconditioning to address this issue?
3) How does the performance of energy GNNs scale with the size of the graph and the number of nodes? Would these be applicable to very large graphs? (I see a comment on this in conclusion section but would be interesting to know the authors intuition about this)
soundness: 4
strengths: 1) Adoption of GNNs in real-world distributed systems.
2) Theoretical guarantees for the robustness of implicitly-defined GNNs to asynchrony, drawing on concepts from distributed optimization.
3) The proposed energy GNN offers a flexible and expressive way to define convex energy functions, potentially leading to more powerful GNN models.
4) Experiments on both synthetic (chains, counting, sums, coordinates) and benchmark datasets (MUTAG, PROTEINS, PPI) demonstrate the effectiveness of the proposed approach.
summary: In this paper, the authors tackle the issue of asynchrony in graph neural networks (GNNs). Traditional GNNs, particularly multi-layer variants, assume synchronized node updates, which is often unrealistic in real-world distributed systems like robotic swarms or sensor networks. 

The authors identify a class of GNNs called "implicitly-defined" GNNs that are robust to asynchronous or "hogwild" inference. They propose a novel implicitly-defined GNN architecture called "energy GNN", which leverages input-convex neural networks to parameterize a convex energy function. The results from experiments on synthetic multi-agent tasks and benchmark graph datasets demonstrate the superior performance of energy GNNs under asynchrony and their competitiveness even in synchronous settings.
weaknesses: 1) Even though the paper proposes mitigation strategies for trraining implicitly-defined GNNs. They are computationally expensive due to the iterative nature of the forward pass. This would be a challenge in practical scenarios.
2) The experiments primarily focus on single-layer energy GNNs. The performance of multi-layer variants and their scalability to larger graphs are unclear.
3) The convergence of implicitly-defined GNNs can be sensitive to the choice of hyperparameters like step size and convergence tolerance. Have the authors investigated the robustness perspective?
4) While the paper motivates the problem with real-world scenarios(in abstract) the experiments are primarily on synthetic datasets. Unsure how well these would be effective in real world scenarios due to few of the weaknesses pointed above

Meta-review: The paper addresses the issue of graph networks (GNNs) that this networks generate incorrect predictions if nodes update asynchronously. This is an issue in distributed, decentralized multi-agent systems, where it is not feasible to update the nodes at the same time. The authors propose implicitly-defined GNNs that are robust to asynchronous updates and can be used in distributed systems. The paper provides comprehensive mathematical proofs for the convergence of Energy GNNs and demonstrates significant, or competitive, improvements on synthetic datasets and real-world datasets. As a weakness, most experiments were provided on small synthetic datasets, though the target graphs in multi-robot systems are usually small as well. I recommend for Accept (poster).

Decision: Accept (Poster)