Triggered by limitations of graph-based deep learning methods in terms of computational expressivity and model flexibility, recent years have seen a surge of interest in computational models that operate on higher-order topological domains such as hypergraphs and simplicial complexes.
While the increased expressivity of these models can indeed lead to a better classification performance and a more faithful representation of the underlying system, the computational cost of these higher-order models can increase dramatically.
To this end, we here explore a simplicial complex neural network learning architecture based on random walks and fast 1D convolutions (SCRaWl), in which we can adjust the increase in computational cost by varying the length and number of random walks considered while accounting for higher-order relationships.
Importantly, due to the random walk-based design, the expressivity of the proposed architecture is provably incomparable to that of existing message-passing simplicial neural networks.
We empirically evaluate SCRaWl on real-world datasets and show that it outperforms other simplicial neural networks.

§ INTRODUCTION

In recent years, there has been a strong interest in extending graph-based learning methods, in particular graph neural networks, to higher-order domains such as hypergraphs and simplicial complexes.
These higher-order abstractions offer a more comprehensive representation of relationships among entities than traditional graph-based approaches, which are based on pairwise interactions.
Accordingly, it has been shown that we can achieve performance gains in various learning and prediction tasks for complex systems by employing such higher-order models [Benson et al., 2018, Schaub et al., 2021, Battiston & Petri, 2022].
For instance, in the context of graph classification, it is well known that standard message-passing graph neural networks with anonymous inputs are limited in their expressiveness by the one-dimensional Weisfeiler-Leman graph isomorphism test [Morris et al., 2019, Xu et al., 2019, Morris et al., 2023] and thus cannot distinguish certain non-isomorphic graphs.
A concrete example is that such graph neural networks cannot differentiate between a cycle of length $6$ vs. two non-connected triangles.
In contrast, higher-order extensions of graph-neural networks have greater expressiveness and naturally enable us to distinguish such graphs.

However, due to the combinatorial increase in the possible interactions, there is a significant cost in terms of increased memory and compute time associated with the use of higher-order representations.
To alleviate these issues, here we take inspiration from random walk-based graph neural networks on graphs [Tönshoff et al., 2023] and explore the use of a random walk-based learning architecture for simplicial complexes.
The advantages of this strategy are particularly apparent for simplicial complexes and related architectures.
First, by choosing the number of random walks sampled, we can effectively trade off the computational demands with the expressivity of our architecture.
Second, as we can compute fast 1D convolutions on the generated random walk trajectories via fast Fourier transforms, this saves us from computing far more expensive convolutional filters on simplicial complexes.
Interestingly, it was shown in [Tönshoff et al., 2023] that there exist graphs that cannot be distinguished by the classical Weisfeiler-Leman (WL) test, that can be distinguished using such a random walk-based learning strategy, and vice versa that certain graphs cannot be distinguished that are distinguishable by the 1-WL test.
This incomparability to message-passing schemes extends to the case of simplicial complexes, and as a result, our learning architecture has an expressivity that is not comparable to message-passing-based schemes such as [Bodnar et al., 2021].

0.65 !

Sketch of SCRaWl illustrating the individual steps of the model.
We sample a collection of random walks on the simplicial complex and transform them into walk feature matrices.
These are then processed by a 1D convolutional network and pooled into updated simplex states.

We contribute SCRaWl, a novel neural network architecture for simplicial complexes that utilizes random walks on higher-order simplices to incorporate higher-order relations between entities into the neural network.
This leads to an architecture different from existing simplicial neural networks.
We prove that SCRaWl's expressiveness is incomparable to that of message-passing simplicial networks and show that our model outperforms existing approaches on real-world datasets.

Related work
Compared to corresponding graph learning problems, learning from data supported on simplicial complexes has so far received far less attention.
presented a basic convolutional neural network (SNN) for simplicial complexes, which was later extended by Yang et al., 2022 to incorporate information flows between different orders of simplices (SCNN).
proposed a convolutional architecture for trajectory predictions that is rooted in tools from algebraic topology and derived three desirable properties for simplicial neural architectures: permutation equivariance, orientation equivariance, and simplicial awareness.

proposed a message-passing simplicial network (MPSN), which adapts the proven message-passing concept from graph neural networks to simplicial complexes using messages to adjacent, upper-adjacent, and lower-adjacent simplices.
They later extended this line of work to the domain of cell complexes [Bodnar et al., 2021].
Similarly, Goh et al., 2022 and Giusti et al., 2022 independently adapted the well-known attention mechanism from graph attention networks and constructed the simplicial attention networks SAT and SAN, respectively.
Apart from this line of research focused on extending graph models to simplicial complexes, Keros et al., 2022 proposed a graph convolutional model for learning functions parametrized by the $k$-homological features of simplicial complexes.

Random walks on simplicial complexes have been considered in [Schaub et al., 2020, Alev & Lau, 2020, Parzanchevski & Rosenthal, 2017, Mukherjee & Steenbergen, 2016].
and Hacker, 2020 proposed two similar random walk-based representation approaches (simplex2vec and $k$-simplex2vec) to learn embeddings from co-appearing simplices in the walks, with the idea that co-appearing simplices should have similar embeddings while other simplices should be further apart in the embedding space.
extended this work with sc2vec, which takes more connections into account and therefore yields more effective embeddings.
These unsupervised representation learning methods compute embeddings for simplices or simplicial complexes, which can then be used for downstream tasks such as classification, while SCRaWl is a supervised learning method that is trained end-to-end.

For graphs, several neural networks based on random walks have been proposed.
proposed RWNN, a model that integrates a random walk-based kernel into a GNN architecture.
and Eliasof et al., 2022 proposed two architectures that aggregate node features along random walks to learn information based on walk distances.
Our work here extends the approach by , which proposed a random walk-based graph neural network architecture that enriches the random walks with local structural information and processes the resulting feature matrices with a 1D convolutional neural network.
connect some of these concepts and extend ideas behind RWNN to simplicial complexes.
In contrast to SCRaWl, they use walk probabilities to enrich the node and edge features for following message-passing steps, while we make use of explicit random walks as the main building block of our architecture.

Another popular abstraction tool for higher-order data is hypergraphs, for which several learning architectures have been proposed.
In a similar vein to simplicial neural networks, graph neural networks have been abstracted to hypergraphs, e.g., convolutional networks [Yadati et al., 2019], attention networks [Chen et al., 2020], and message-passing networks [Heydari & Livi, 2022, Chien et al., 2022].
Very recently, Behrouz et al., 2023 presented a walk-based hypergraph neural network that uses temporal walks to extract higher-order causal patterns.
Compared to their approach, SCRaWl can capture more structural properties as part of the feature matrices resulting from the random walks.

The remainder of this manuscript is structured as follows.
<Ref> gives an overview of the mathematical concepts of simplicial complexes and random walks on them.
In <ref>, we present our proposed architecture SCRaWl in detail.
<Ref> analyzes the expressiveness of SCRaWl compared to MPSN and <ref> evaluates SCRaWl empirically on a variety of datasets.

§ BACKGROUND

Here, we present an elementary overview of concepts used to process signals defined on simplicial complexes.
For more details, see [Bredon, 1993, Hatcher, 2002] for background in algebraic topology and [Schaub et al., 2021, Roddenberry et al., 2022] for topological signal processing.

We use $[n]$ to denote the set of integers from $0$ to $n-1$.
Vectors are denoted by boldface lowercase letters $\vec{v}$, matrices and tensors by capital letters $\mat{F}$.
For ease of notation, we index rows and columns from $0$.
$\mat{F}_{i,-}$ and $\mat{F}_{-,j}$ denote the $i$-th row and $j$-th column, respectively, of a matrix $\mat{F}$.
$\mat{0}_n$ and $\mat{1}_n$ denote an $n$-dimensional all-zero and all-one vector, respectively.

An abstract simplicial complex (SC) [Bredon, 1993, Hatcher, 2002] $X$ consists of a finite set of points $\mathcal{V}$, and a set of non-empty subsets of $\mathcal{V}$ that is closed under taking non-trivial subsets.
A $k$-simplex $\mathcal{S}^k \in X$ is a subset of $\mathcal{V}$ with $k+1$ points and if $\mathcal{S}^{k} \in X$, then for all non-empty $\mathcal{S}^{k-1} \subset \mathcal{S}^k$, $\mathcal{S}^{k-1} \in X$.
We denote the set of $k$-simplices in $X$ by $X_k$ and their cardinality by $n_k = |X_k|$.
Furthermore, for $k > 0$ the faces $\mathfrak{F}(\mathcal{S}^k) = \left\{ \mathcal{S}^{k-1} \in X_{k-1} \mid \mathcal{S}^{k-1} \subset \mathcal{S}^k \right\}$ of $\mathcal{S}^k$ are the subsets of $\mathcal{S}^k$ with cardinality $k$. We set $\mathfrak{F}(\mathcal{S}^k) = \emptyset$ for $k=0$.
If $\mathcal{S}^{k-1}$ is a face of $\mathcal{S}^k$, $\mathcal{S}^k$ is called a coface of $\mathcal{S}^{k-1}$.
We denote the set of cofaces of $\mathcal{S}^k$ by $\mathfrak{C}(\mathcal{S}^k)$.

We respectively denote the set of lower and upper adjacent simplices of $\mathcal{S}^k$ by:

\begin{align}
N^\downarrow(\mathcal{S}^k) & = \bigl\{ \mathcal{S}'^k \in X_k \mid \mathfrak{F}(\mathcal{S}^k) \cap \mathfrak{F}(\mathcal{S}'^k) \neq \emptyset \bigr\}, \\
N^\uparrow(\mathcal{S}^k) & = \bigl\{ \mathcal{S}'^k \in X_k \mid \mathfrak{C}(\mathcal{S}^k) \cap \mathfrak{C}(\mathcal{S}'^k) \neq \emptyset \bigr\}.
\end{align}

Note that SCs can be understood as an extension of graphs: $X_0$ is the set of vertices, $X_1$ is the set of edges, $X_2$ is the set of filled-in triangles (not all $3$-cliques have to be filled-in triangles), and so on.

Simplices are equipped with features $f_k : X_k \to \mathbb{R}^{d_k}$.
For flexibility, we allow the feature spaces to be different for each simplex order, as is common in real-world data sets.
If no features are supported on a simplex order, we set $d_k = 0$.

§ METHOD

At its core, SCRaWl extracts information from data supported on a simplicial complex by sampling random walks on the underlying complex.
These walks are then transformed into feature matrices, which are used to update the states of the simplices in the subsequent steps.
More precisely, SCRaWl consists of three steps, as are illustrated in <ref>:
First, we sample random walks on simplicial complexes, as described in <ref>.
Second, as discussed in <ref>, we transform these walks into feature matrices that encode the random walks, the simplex states that appear on them, and the local adjacencies between simplices along the walk.
Finally, we process these feature matrices with a simple convolutional neural network to update the hidden states of the simplices (see <ref>).

§.§ Random Walks on Simplicial Complexes

A walk of length $\ell$ on the $k$-simplices of a simplicial complex $X$ is a sequence of $k$-simplices $(v_0, \dots, v_{\ell - 1}) \in X_k^{l}$ such that $v_i$ and $v_{i+1}$ share a common face or coface $e_i$ for all $i \in [\ell]$, i.e., $v_{i+1} \in N^\downarrow(v_i) \cup N^\uparrow(v_i)$.
Note that different from the case of a graph $e_i$ is not uniquely defined given $v_{i}$ and $v_{i+1}$, since $v_{i+1}$ can be reached by either a common face or a common coface.
We keep track of the connection used to reach the next neighbor and define $W  = (v_0, e_0, v_1, e_1, \dots, e_{\ell - 2}, v_{\ell - 1})$.

While there are many possible sampling schemes for random walks on SCs, for simplicity, we only consider two elementary sampling methods:
[label= (*)]

* uniform connection sampling, and

* uniform neighbor sampling.

For both sampling strategies, we obtain a random walk on the $k$-simplices of $X$ as follows:
First, sample a random starting $k$-simplex $v_0 \sim \mathcal{U}(X_k)$ and then sample the subsequent simplices in a walk as described in the following subsections.
Using this strategy we sample $m$ random walks on the simplicial complex, which can be chosen at runtime and can vary between training and prediction.
In the case $m = |X|$, we opt to start a random walk from each simplex, i.e., we do not sample $v_0$ uniformly at random but choose every simplex once as $v_0$.
Note that by choosing a smaller $m$ we can reduce the computational cost we incur for learning and processing the data.

Uniform Connection Sampling

Starting from simplex $v_0$ we iteratively sample a face or coface $e_i \sim \mathcal{U} \bigl( \mathfrak{F}(v_i) \cup \mathfrak{C}(v_i) \bigr)$ and a neighbor

\begin{align}
v_{i+1} \sim \begin{dcases}
\mathcal{U} \bigl( \mathfrak{C}(e_i) \bigr) & \text{if } e_i \in \mathfrak{F}(v_i) \\
\mathcal{U} \bigl( \mathfrak{F}(e_i) \bigr) & \text{if } e_i \in \mathfrak{C}(v_i)
\end{dcases}
\end{align}

of $v_i$ that is connected to $v_i$ via $e_i$ uniformly at random.
This way, a neighbor with more connections to the current simplex is more likely to be sampled.

Uniform Neighbor Sampling

Starting from simplex $v_0$ we iteratively sample a neighboring simplex
$v_{i+1} \sim \mathcal{U} \bigl( N^\downarrow(v_i) \cup N^\uparrow(v_i) \bigr)$ and a suitable connection

\begin{align}
e_i \sim \mathcal{U} \Bigl( \bigl( \mathfrak{F}(v_i) \cap \mathfrak{F}(v_{i+1}) \bigr) \cup \bigl( \mathfrak{C}(v_i) \cap \mathfrak{C}(v_{i+1}) \bigr) \Bigr)
\end{align}

uniformly at random.
This way, higher-order simplices do not have a disproportionate influence on the walk due to their many low-order connections.

§.§ Walk Feature Matrices

Information flow in SCRaWl modules for simplex orders $0$, $1$, and $2$.
Using the walks shown in <ref>, we compute the feature matrices $F_0, F_1$, and $F_2$ based on simplex features $f$ and a window size of $s = 4$.
The feature matrices are convolved with a 1D CNN while keeping track of the central simplex within the window (values to the left of the convolution matrices).
Convolved matrices are pooled into simplices grouped by the central simplices.

Based on the sampled collection of random walks ${\{ W_j \}}_{j \in [m]}$ from the previous step and a local window size $s$, we transform each walk $W_j$ into feature matrix $\mat{F}_{W_j}$.
The feature matrix consists of $6$ sub-matrices:
[label= (*)]

* one matrix for the features of the current simplex,

* one matrix for the features of the face used to traverse to the current simplex,

* one matrix for the features of the coface used to traverse to the current simplex, and

* the last three matrices encode local structural information about the simplicial complex seen during the walk.

More formally, given a walk $W_j$ of length $\ell$ on $k$-simplices, embedding functions $f_{k'}: X_{k'} \to \mathbb{R}^{d_{k'}}$ for $i \in \{k-1, k, k+1\}$, and a local window size $s > 0$, we define the walk feature matrix $\mat{F}_{W_j}$ as

\begin{align}
\setlength{\arraycolsep}{4pt}
\mat{F}_{W_j} = \begin{bmatrix}
\mat{F}_{W_j}^\rightarrow & \mat{F}_{W_j}^\downarrow & \mat{F}_{W_j}^\uparrow & \mat{I}_{W_j,s} & \mat{A}_{W_j,s}^\downarrow & \mat{A}_{W_j,s}^\uparrow
\end{bmatrix} \in \mathbb{R}^{\ell \times (d_k + d_{k-1} + d_{k+1} + 3s - 2)}.
\end{align}

The first sub-matrix $\mat{F}_{W_j} \in \mathbb{R}^{\ell \times d_k}$ records the features of the $k$-simplices appearing on the walk $W$:

\begin{align}
{(\mat{F}_W^\rightarrow)}_{i,-} & = f(v_i) \; \text{ for } \; i \in [\ell].
\end{align}

The second and third sub-matrices $\mat{F}_{W_j}^\downarrow$ and $\mat{F}_{W_j}^\uparrow$ record the features of the faces and cofaces, respectively, that have been used in the random walk to get from the last simplex to the current one:

\begin{align}
{(\mat{F}_{W_j}^\downarrow)}_{i,-} & = \begin{dcases}
f(e_{i-1}), & \text{if } i > 0 \text{ and } |e_{i-1}| < k \\
\mat{0}_{d_{k-1}}, & \text{otherwise.}
\end{dcases} \\
{(\mat{F}_{W_j}^\uparrow)}_{i,-} & = \begin{dcases}
f(e_{i-1}), & \text{if } i > 0 \text{ and } |e_{i-1}| > k \\
\mat{0}_{d_{k+1}}, & \text{otherwise.}
\end{dcases}
\end{align}

Note that whenever a face or a coface feature is present, the entries of the other matrix are zero.
Further, assume that $d_{-1} = d_{k_\text{max} + 1} = 0$ and hence $\mat{F}_{W_j}^\downarrow$ and $\mat{F}_{W_j}^\uparrow$ are empty for walks on $0$-simplices and on maximum order simplices, respectively, where $k_\text{max}$ is the maximum order in the SC.

The remaining three sub-matrices $\mat{I}_{W_j} \in {\{0,1\}}^{\ell \times s}, \mat{A}_{W_j}^\downarrow \in {{\{0,1\}}}^{\ell \times (s-1)}$ and $\mat{A}_{W_j}^\uparrow \in {\{0,1\}}^{\ell \times (s-1)}$ record structural information in the form of local identity and local adjacency features:

\begin{align}
{(\mat{I}_{W_j,s})}_{i,i'} & = \begin{dcases}
1, & \text{if } i \geq i' \text{ and } v_i = v_{i'} \\
0, & \text{otherwise.}
\end{dcases} \\
{(\mat{A}_{W_j,s}^\downarrow)}_{i,i'} & = \begin{dcases}
1, & \text{if } i > i' \text{ and } v_{i'-i} \in N^\downarrow(v_i) \\
0, & \text{otherwise.}
\end{dcases} \\
{(\mat{A}_{W_j,s}^\uparrow)}_{i,i'} & = \begin{dcases}
1, & \text{if } i > i' \text{ and } v_{i'-i} \in N^\uparrow(v_i) \\
0, & \text{otherwise.}
\end{dcases}
\end{align}

In other words, $\mat{I}_{W_j,s}$ keeps track of whether the current simplex has been visited before within the local window $s$, while $\mat{A}_{W_j,s}^\downarrow$ and $\mat{A}_{W_j,s}^\uparrow$ keep track of whether the current simplex is lower or upper adjacent to a previously visited simplex within the local window $s$, respectively.

<Ref> (left) shows an example of the walk feature matrices for the walks on the vertices, the edges, and the triangles shown in <ref>.
Note in particular the last three sub-matrices for each walk:
For $W_0$, vertex $v_4$ is visited twice, once in the third and once in the fifth step, thus we have ${(\mat{I}_{W_0,4})}_{5,2} = 1$.
In the same walk, we go from $v_1$ via the edges $(v_1, v_2)$ and $(v_2, v_4)$ to $v_4$, which is also directly connected with $v_1$.
Thus we have ${(\mat{A}^\uparrow_{W_0,4})}_{2,0} = 1$.
Similarly, for $W_1$, we go from $(v_2, v_4)$ over two steps to $(v_4, v_6)$, which share a common vertex $v_4$.
Thus we have ${(\mat{A}^\downarrow_{W_1,4})}_{3,0} = 1$.

§.§ SCRaWl Module

A SCRaWl module on layer $t$ operating on order $k$-simplices takes as input the random walks on $k$-simplices and the current hidden states $\mat{H}_k^{t-1}, \mat{H}_{k-1}^{t-1}$, and $\mat{H}_{k+1}^{t-1}$ of the previous layer $t-1$.
These matrices are then used in lieu of the simplex, face, and coface features in the walk feature matrices, respectively.
The module first computes the feature matrices $\mat{F}_{W_j}$ for each walk ${W_j}$ as described in the previous subsection and stacks them into a tensor $\mat{F}_k$.

We process the walk features $\mat{F}_k$ by a 1D convolutional network $\operatorname{CNN}_k^t$ (without padding) performed over the walk steps to obtain the convolved matrix $\mat{C}_k^t \in \mathbb{R}^{m_k \times (l - s) \times d}$, where $m_k$ is the number of walks on $k$-simplices.
The CNNs can be configured flexibly depending on the task.

For each row of the convolved matrix $\mat{C}_k^t$, we separately keep track of the simplex in the center row of the convolution window, as is illustrated in <ref> with the indices next to the convolution matrices.
For each $k$-simplex, we pool all entries of $\mat{C}_k^t$ that correspond to that simplex using mean or sum pooling.
This operation results in a vector $\vec{p}^t_k(v)$ for each $k$-simplex $v$, which is then fed into a trainable multilayer perceptron $\operatorname{Update}_k^t$ to obtain the updated hidden state $\mat{H}_k^t(v)$.

The process of a SCRaWl module is illustrated in <ref>.
As can be seen, the new hidden state of $k$-simplices depends on the hidden state of upper and lower adjacent simplices, and hence, over multiple layers, information is propagated through different simplex orders.

§.§ Complete SCRaWl Architecture

0.9 !

Architecture of a SCRaWl model operating on simplex orders $0$, $1$, and $2$ with $L$ SCRaWl layers and an user-defined output layer.
We reuse the same collection of random walks ${\{W_j\}}_{j \in [m]}$ in each SCRaWl module for performance reasons.

Based on the aforementioned building blocks, we can assemble the SCRaWl model, as illustrated in <ref>:
Given a range of simplex orders $[K]$ and the number of layers $L$ as hyperparameters, the model consists of $L$ SCRaWl layers, each consisting of $K$ SCRaWl modules, one for each simplex order.
Within one layer, SCRaWl modules run in parallel for each simplex order $k \in [K]$.
We add skip connections for each module.
Information flows between different simplex orders in between each layer in the form of the hidden states of faces and cofaces appearing in the random walks.

At the end of the last layer, we obtain hidden states $\mat{H}_k^L$ for each simplex order $k \in [K]$.
These can then be processed by user-defined output layers, e.g., by adding multilayer perceptrons that transform the hidden states into a classification or regression output, to obtain the final output of the model.

§.§ Performance Considerations

Most of the computational complexity of SCRaWl is due to the sampling of random walks.
Hence, we only sample the random walks once in each epoch and reuse this collection in all layers, i.e., each module transforms them into module-specific feature matrices with different feature values (first three sub-matrices) but reuses the same random walks for that.
This is also depicted in <ref> with the same collection ${\{W_j\}}_{j \in [m]}$ of random walks passed to each layer.
We further describe how to sample the random walks more efficiently using boundary maps.

The structure of a simplicial complex can be encoded by boundary maps $\mat{B}_k$, which record the incidence relations between $(k-1)$-simplices and $k$-simplices [Hatcher, 2002].
Rows and columns of $\mat{B}_k$ are indexed by $(k-1)$-simplices and by $k$-simplices, respectively.
The entry $(i, j)$ is set to $1$ or $0$ depending on whether the $i$-th $(k-1)$-simplex is incident to the $j$-th $k$-simplex.
Thus, the matrix $\mat{B}_1$ is the unsigned vertex-to-edge incidence matrix, and $\mat{B}_2$ is the edge-to-triangle incidence matrix.
For ease of notation, we define $\mat{B}_0 = \mat{B}_{k+1} \mat{0}$.
Further define $\mat{S}_k = \operatorname{hstack}(\mat{B}_k^\top, \mat{B}_{k+1})$.

For uniform connection sampling on the $k$-simplices, we can efficiently sample the next connection $j$ from $\mat{S}_k$, as the connections for the $i$-th $k$-simplex are precisely those with value one in ${(\mat{B}_k)}_{-,i}$ and ${(\mat{B}_{k+1})}_{i,-}$.
If the sampled index $j \leq n_{k-1}$, the walk transitions over the face $j$, otherwise over the coface $j - n_{k-1}$.
The next walk simplex can then be sampled uniformly from $\mat{S}_k^\top$.

For uniform neighbor sampling, compute $\mat{A}_k = \mat{B}_k^\top \mat{B}_k + \mat{B}_{k+1} \mat{B}_{k+1}^\top$ and $\hat{\mat{A}}_k = \operatorname{sign}(\mat{A}_k)$.
The neighbor of the $i$-th $k$-simplex can be sampled uniformly at random from ${(\hat{\mat{A}}_k)}_{i,-}$.

§ EXPRESSIVE POWER

It is known that the expressive power of message-passing GNNs with anonymous inputs is limited by the WL graph isomorphism test [Morris et al., 2019, Xu et al., 2019, Morris et al., 2023] and thus such GNNs cannot distinguish between non-isomorphic graphs that cannot be separated by color refinement.
In <cit.> a simplicial version of the WL test (SWL) was introduced and it was shown that the expressiveness of message-passing simplicial networks (MPSN) is limited by this test.
The following theorem asserts that for SCRaWl this analysis does not apply.

The expressiveness of SCRaWl is incomparable to the expressiveness of an MPSN, i.e., there are simplicial complexes that can be distinguished by SCRaWl but not by MPSN and vice versa.

Note that when applied to 1-simplicial complexes and considering vertex colorings, SWL is equivalent to WL and by extension, the expressiveness of MPSNs corresponds to the expressiveness of GNNs [Bodnar et al., 2021].
Note also that SCRaWl is equivalent to CRaWl when applied to 1-simplicial complexes.
Thus the theorem follows from <cit.>.

The expressiveness inclusions of the various architectures are shown in <ref>.
This implies that the expressiveness of SCRaWl is incomparable to the WL test and its simplicial counterpart.

valign=b, minipage=0.47

Expressiveness relations of SCRaWl to other neural network models.
SCRaWl and MPSN are strict extensions of their graph counterparts CRaWl and message-passing GNNs, respectively, and their expressive power is incomparable.

valign=b, minipage=0.47

Example of a citation simplicial complex.
Each vertex represents an author and connections (edges and triangle) indicate co-authorship.
The value of the connection is the total number of citations of the papers the authors have written together.

§ EXPERIMENTAL RESULTS

We evaluate SCRaWl on a variety of datasets and compare it to other simplicial neural networks[ Source code and datasets are available at <https://git.rwth-aachen.de/netsci/scrawl>. Source code is available in the supplementary material. ].
Following the training procedure of Dwivedi et al., 2023, we use the Adam optimizer with an initial learning rate $10^{-3}$.
The learning rate decays with a factor of $0.5$ if the validation loss does not improve for $10$ epochs.
Training is stopped once the learning rate drops below $10^{-6}$.
For a full list of hyperparameters and a more detailed description of the training setup, see <ref>.

§.§ Semantic Scholar Co-Authorship Network

Following Ebli et al., 2020, we use SCRaWl to impute missing citation counts for a subset of the Semantic Scholar co-authorship network.
Vertices represent distinct authors, and a paper with $k+1$ authors implies a $k$-simplex connecting these authors in the network (see <ref>).
The value on a $k$-simplex is the sum of citations of the papers written together by the respective $k+1$ authors (including papers with more than these authors).
Citation values are randomly deleted at a rate of $10\%$, $20\%$, $30\%$, $40\%$, and $50\%$ and replaced by the median of the remaining known citations.

We report the imputation accuracies of SCRaWl in <ref> (left) and compare them with $k$-simplex2vec [Hacker, 2020], SNN [Ebli et al., 2020], SCNN [Yang et al., 2022], SAT [Goh et al., 2022], SAN [Giusti et al., 2022], and the MPSN architecture presented in [Bodnar et al., 2021].
For readability, we omitted the results for $k$-simplex2vec, SNN, and SAT in <ref> (right).
Full results are available in <ref> in the appendix.
Following the previous work, a citation value is correctly imputed if the prediction is within $\pm 5\%$[The SNN paper used a threshold of $\pm 10\%$ which was later changed to $\pm 5\%$.] of the true value.
We repeat each experiment $10$ times and report the mean and standard deviation for missing rates of $10\%$, $20\%$, $30\%$, $40\%$, and $50\%$.
The left plot shows the overall accuracies of our model on the different missing rates while the right plot compares the accuracies of the different models on individual simplex orders.

We see that SCRaWl reaches accuracies of $95\%$ to $99\%$ for all missing rates after about $175$ epochs.
Across all experiments, we consistently see a drop in accuracy around epoch $30$ to $50$, which is then overcome.
We forced a minimum of $100$ epochs to avoid early stoppings within this drop.
We also see that SCRaWl is on par with MPSN and outperforms all other models for most missing rates and simplex orders.
Especially for missing rates above $10\%$, the improvement is substantial.
In addition, imputation accuracies are more consistent across different simplex orders for SCRaWl, whereas, for the other models, the accuracy declines faster for smaller simplices than for higher-order simplices.

valign=b, minipage=0.45

valign=b, minipage=0.45

Imputation accuracies on the Semantic Scholar dataset with different percentages of missing citations.
We report the mean and standard deviation over $10$ runs.
Left: Training progress of SCRaWl.
Right: Accuracies of SCRaWl, SIN (MPSN), SAN, and SCNN for each simplex order as more and more data is missing.
Increasingly lighter lines indicate higher simplex orders.

§.§ Social Contact Networks

In a second set of experiments, we perform vertex classification on the primary-school and high-school social contact datasets [Stehlé et al., 2011, Chodrow et al., 2021].
The datasets are constructed based on students' interactions recorded by wearable sensors.
Vertices represent students, and a simplex corresponds to a group of students that were in proximity to each other at a given time.
Vertices are labeled with the classroom to which the student belongs, resulting in a total of $12$ classes for the primary school dataset and $10$ classes for the high school dataset.

We use a cross-entropy loss to train the network on the dataset with $40\%$ of the vertex classes missing and report the validation accuracy of the imputed classes in <ref> (left).
For the primary school dataset, SCRaWl reaches an average accuracy of $0.927 \pm 0.026$ after about $160$ epochs, and for the high school dataset a perfect accuracy of $1.0 \pm 0.0$ after about $150$ epochs.

We compare this result with the performance of MPSN on the same datasets (middle).
Although MPSN generally converges faster in our experiment, it only achieves an average accuracy of $0.727 \pm 0.033$ and $0.943 \pm 0.033$ for the primary and high school dataset, respectively, and is thus outperformed by SCRaWl, especially on the primary school dataset.
Comparisons with SCNN, $k$-simplex2vec and some hypergraph models can be found in <ref> in the appendix.
The other simplicial approaches did not perform well on the datasets, with SCNN reaching an average accuracy of $0.166 \pm 0.056$ for the primary school dataset and $0.228 \pm 0.075$ for the high school dataset.
$k$-simplex2vec reached an average accuracy of $0.398 \pm 0.011$ and $0.308 \pm 0.019$ for the primary and high school dataset, respectively.
The hypergraph approaches performed consistently at most around $0.9$ accuracy on the high school dataset, $10$ percent below our result.
For the primary school dataset, CAt-Walk marginally outperformed SCRaWl with an average accuracy of $0.932 \pm 0.025$.

Vertex classification accuracies of SCRaWl, MPSN, and CRaWl on the social contact datasets.
We report the mean and standard deviation of the validation accuracies over $5$ runs.

Importance of Higher-Order Interactions
We repeated the same experiment using the graph-based CRaWl architecture on the underlying graph skeleton of the datasets, i.e., the graph only records pairwise proximities between students but does not record whether interactions between more than two people occurred.
Keeping all other parameters the same, we found that CRaWl performance is clearly worse than SCRaWl on both datasets (see <ref>).
This highlights the importance of higher-order interactions in the dataset and the ability of SCRaWl to capture them.

§ CONCLUSION

We proposed a simplicial neural network architecture SCRaWl based on random walks and 1D convolutions.
The architecture builds on the idea of [Tönshoff et al., 2023] to use sampled random walks for graph learning tasks.
We demonstrated that these ideas can be extended to work on simplicial complexes if we carefully adjust the notion of random walks.
We showed that the resulting architecture outperforms existing simplicial neural network architectures on a co-authorship network.

Other Complexes
While we focused on simplicial complexes in this work, the architecture can naturally be applied to other types of complexes such as (regular) cell complexes without any changes.
Informally, a cell complex is a collection of cells, where each cell's boundary is also in the complex.
Starting with points as $0$-cells, $1$-cells as lines with their two endpoints as their boundary, a $k$-cell is defined such that its boundary consists of a collection of $(k - 1)$-cells.
A polygon for instance is a $2$-cell with a boundary consisting of a set of lines.
This highlights the difference to a simplicial complex, as a $2$-simplex must form a triangle, while a $2$-cell may be any polygon.

Future Work
The empirical results demonstrate the potential of the proposed architecture.
However, there are several avenues for future work worth exploring.
First, we plan to investigate the performance of SCRaWl in a broader set of applications.
Simplicial complexes have gained a lot of interest for edge flow problems, and it remains to be seen how well SCRaWl performs on these tasks.
While we have paid special attention to the computational complexity in our implementation, due to the many simplices and hence walk feature matrices in a simplicial complex compared to a graph, we expect that using SCRaWl with one random walk per simplex will not scale sufficiently.
Approximation guarantees for different sampling schemes are thus an interesting direction for future work.

§.§.§ Acknowledgments

We acknowledge funding by the Ministry of Culture and Science (MKW) of the German State of North Rhine-Westphalia (“NRW Rückkehrprogramm”) and the European Union (ERC, HIGH-HOPeS, 101039827). Views and opinions expressed are however those of the authors only and do not necessarily reflect those of the European Union or the European Research Council Executive Agency. Neither the European Union nor the granting authority can be held responsible for them.

[Alev & Lau, 2020]
Vedat Levi Alev and Lap Chi Lau.
Improved analysis of higher order random walks and applications.
In Proceedings of the 52nd Annual ACM SIGACT Symposium on
Theory of Computing, pp.  1198–1211, 2020.

[Battiston & Petri, 2022]
Federico Battiston and Giovanni Petri.
Higher-Order Systems.
Springer International Publishing, 2022.

[Behrouz et al., 2023]
Ali Behrouz, Farnoosh Hashemi, Sadaf Sadeghian, and Margo Seltzer.
CAT-Walk: Inductive Hypergraph Learning via Set Walks.
In Thirty-seventh Conference on Neural Information Processing
Systems, 2023.
URL <https://openreview.net/forum?id=QG4nJBNEar>.

[Benson et al., 2018]
Austin R Benson, Rediet Abebe, Michael T. Schaub, Ali Jadbabaie, and Jon
Simplicial closure and higher-order link prediction.
Proceedings of the National Academy of Sciences, 1150
(48):0 E11221–E11230, 2018.

[Billings et al., 2019]
Jacob Charles Wright Billings, Mirko Hu, Giulia Lerda, Alexey N Medvedev,
Francesco Mottes, Adrian Onicas, Andrea Santoro, and Giovanni Petri.
Simplex2vec embeddings for community detection in simplicial
arXiv preprint arXiv:1906.09068, 2019.

[Bodnar et al., 2021]
Cristian Bodnar, Fabrizio Frasca, Nina Otter, Yuguang Wang, Pietro Liò,
Guido F Montufar, and Michael Bronstein.
Weisfeiler and Lehman Go Cellular: CW Networks.
In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wortman
Vaughan (eds.), Advances in Neural Information Processing Systems,
volume 34, pp.  2625–2640. Curran Associates, Inc., 2021a.

[Bodnar et al., 2021]
Cristian Bodnar, Fabrizio Frasca, Yuguang Wang, Nina Otter, Guido F Montufar,
Pietro Lió, and Michael Bronstein.
Weisfeiler and Lehman Go Topological: Message Passing Simplicial
In Marina Meila and Tong Zhang (eds.), Proceedings of the 38th
International Conference on Machine Learning, volume 139 of
Proceedings of Machine Learning Research, pp.  1026–1037. PMLR,
18–24 Jul 2021b.
URL <https://proceedings.mlr.press/v139/bodnar21a.html>.

[Bredon, 1993]
Glen E. Bredon.
Topology and Geometry.
Graduate texts in mathematics. Springer-Verlag, New York, 1993.

[Chen et al., 2020]
Chaofan Chen, Zelei Cheng, Zuotian Li, and Manyi Wang.
Hypergraph Attention Networks.
In 2020 IEEE 19th International Conference on Trust, Security
and Privacy in Computing and Communications (TrustCom), pp.  1560–1565,

[Chien et al., 2022]
Eli Chien, Chao Pan, Jianhao Peng, and Olgica Milenkovic.
You are AllSet: A Multiset Function Framework for Hypergraph Neural
In International Conference on Learning Representations, 2022.
URL <https://openreview.net/forum?id=hpBTIv2uy_E>.

[Chodrow et al., 2021]
Philip S Chodrow, Nate Veldt, and Austin R Benson.
Generative hypergraph clustering: From blockmodels to modularity.
Science Advances, 2021.

[Dwivedi et al., 2023]
Vijay Prakash Dwivedi, Chaitanya K. Joshi, Anh Tuan Luu, Thomas Laurent, Yoshua
Bengio, and Xavier Bresson.
Benchmarking Graph Neural Networks.
Journal of Machine Learning Research, 240
(43):0 1–48, 2023.
URL <http://jmlr.org/papers/v24/22-0567.html>.

[Ebli et al., 2020]
Stefania Ebli, Michaël Defferrard, and Gard Spreemann.
Simplicial Neural Networks.
In TDA & Beyond, 2020.

[Eliasof et al., 2022]
Moshe Eliasof, Eldad Haber, and Eran Treister.
pathGCN: Learning General Graph Spatial Operators from Paths.
In Kamalika Chaudhuri, Stefanie Jegelka, Le Song, Csaba Szepesvari,
Gang Niu, and Sivan Sabato (eds.), Proceedings of the 39th
International Conference on Machine Learning, volume 162 of
Proceedings of Machine Learning Research, pp.  5878–5891. PMLR,
17–23 Jul 2022.
URL <https://proceedings.mlr.press/v162/eliasof22a.html>.

[Giusti et al., 2022]
Lorenzo Giusti, Claudio Battiloro, Paolo Di Lorenzo, Stefania Sardellitti, and
Sergio Barbarossa.
Simplicial Attention Neural Networks.
arXiv preprint arXiv:2203.07485, 2022.

[Goh et al., 2022]
Christopher Wei Jin Goh, Cristian Bodnar, and Pietro Lio.
Simplicial Attention Networks.
In ICLR 2022 Workshop on Geometrical and Topological
Representation Learning, 2022.

[Hacker, 2020]
Celia Hacker.
k-simplex2vec: a simplicial extension of node2vec.
In TDA & Beyond, 2020.

[Hajij et al., 2024]
Mustafa Hajij, Mathilde Papillon, Florian Frantzen, Jens Agerberg, Ibrahem
AlJabea, Ruben Ballester, Claudio Battiloro, Guillermo Bernárdez, Tolga
Birdal, Aiden Brent, Peter Chin, Sergio Escalera, Odin Hoff Gardaa,
Gurusankar Gopalakrishnan, Devendra Govil, Josef Hoppe, Maneel Reddy Karri,
Jude Khouja, Manuel Lecha, Neal Livesay, Jan Meißner, Soham Mukherjee,
Alexander Nikitin, Theodore Papamarkou, Jaro Pr'ilepok, Karthikeyan Natesan
Ramamurthy, Paul Rosen, Aldo Guzm'an-S'aenz, Alessandro Salatiello,
Shreyas N. Samaga, Michael T. Schaub, Luca Scofano, Indro Spinelli, Lev
Telyatnikov, Quang Truong, Robin Walters, Maosheng Yang, Olga Zaghen, Ghada
Zamzmi, Ali Zia, and Nina Miolane.
Topox: A suite of python packages for machine learning on topological
domains, 2024.

[Hatcher, 2002]
Allen Hatcher.
Algebraic Topology.
Cambridge University Press, 2002.

[Heydari & Livi, 2022]
Sajjad Heydari and Lorenzo Livi.
Message Passing Neural Networks for Hypergraphs.
In Elias Pimenidis, Plamen Angelov, Chrisina Jayne, Antonios
Papaleonidas, and Mehmet Aydin (eds.), Artificial Neural Networks and
Machine Learning – ICANN 2022, pp.  583–592, Cham, 2022. Springer Nature
ISBN 978-3-031-15931-2.

[Jin et al., 2022]
Di Jin, Rui Wang, Meng Ge, Dongxiao He, Xiang Li, Wei Lin, and Weixiong Zhang.
RAW-GNN: RAndom Walk Aggregation based Graph Neural Network.
In Luc De Raedt and Luc De Raedt (eds.), Proceedings of the
31st International Joint Conference on Artificial Intelligence, IJCAI 2022,
IJCAI International Joint Conference on Artificial Intelligence, pp.  2108–2114. International Joint Conferences on Artificial Intelligence, 2022.

[Keros et al., 2022]
Alexandros D Keros, Vidit Nanda, and Kartic Subr.
Dist2Cycle: A Simplicial Neural Network for Homology Localization.
Proceedings of the AAAI Conference on Artificial Intelligence,
360 (7):0 7133–7142, Jun. 2022.

[Morris et al., 2019]
Christopher Morris, Martin Ritzert, Matthias Fey, William L. Hamilton, Jan Eric
Lenssen, Gaurav Rattan, and Martin Grohe.
Weisfeiler and Leman Go Neural: Higher-Order Graph Neural Networks.
In Proceedings of the AAAI Conference on Artificial
Intelligence, volume 33, pp.  4602–4609, 2019.

[Morris et al., 2023]
Christopher Morris, Yaron Lipman, Haggai Maron, Bastian Rieck, Nils M. Kriege,
Martin Grohe, Matthias Fey, and Karsten Borgwardt.
Weisfeiler and leman go machine learning: The story so far, 2023.

[Mukherjee & Steenbergen, 2016]
Sayan Mukherjee and John Steenbergen.
Random walks on simplicial complexes and harmonics.
Random Structures & Algorithms, 490 (2):0
379–405, 2016.

[Nikolentzos & Vazirgiannis, 2020]
Giannis Nikolentzos and Michalis Vazirgiannis.
Random Walk Graph Neural Networks.
In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin
(eds.), Advances in Neural Information Processing Systems, volume 33,
pp.  16211–16222. Curran Associates, Inc., 2020.

[Parzanchevski & Rosenthal, 2017]
Ori Parzanchevski and Ron Rosenthal.
Simplicial complexes: spectrum, homology and random walks.
Random Structures & Algorithms, 500 (2):0
225–261, 2017.

[Paszke et al., 2019]
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory
Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban
Desmaison, Andreas Kopf, Edward Yang, Zachary DeVito, Martin Raison, Alykhan
Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith
PyTorch: An Imperative Style, High-Performance Deep Learning
In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché Buc,
E. Fox, and R. Garnett (eds.), Advances in Neural Information
Processing Systems 32, pp.  8024–8035. Curran Associates, Inc., 2019.

[Roddenberry et al., 2021]
T. Mitchell Roddenberry, Nicholas Glaze, and Santiago Segarra.
Principled Simplicial Neural Networks for Trajectory Prediction.
In Marina Meila and Tong Zhang (eds.), Proceedings of the 38th
International Conference on Machine Learning, volume 139 of
Proceedings of Machine Learning Research, pp.  9020–9029. PMLR,
18–24 Jul 2021.
URL <https://proceedings.mlr.press/v139/roddenberry21a.html>.

[Roddenberry et al., 2022]
T. Mitchell Roddenberry, Michael T. Schaub, and Mustafa Hajij.
Signal Processing on Cell Complexes.
In ICASSP 2022 - 2022 IEEE International Conference on
Acoustics, Speech and Signal Processing (ICASSP), pp.  8852–8856, 2022.

[Schaub et al., 2020]
Michael T. Schaub, Austin R Benson, Paul Horn, Gabor Lippner, and Ali
Random Walks on Simplicial Complexes and the Normalized Hodge
SIAM Review, 620 (2):0 353–391, 2020.

[Schaub et al., 2021]
Michael T. Schaub, Yu Zhu, Jean-Baptiste Seby, T. Mitchell Roddenberry, and
Santiago Segarra.
Signal processing on higher-order networks: Livin’ on the edge...
and beyond.
Signal Processing, 187:0 108149, 2021.
ISSN 0165-1684.

[Stehlé et al., 2011]
Juliette Stehlé, Nicolas Voirin, Alain Barrat, Ciro Cattuto, Lorenzo
Isella, Jean-François Pinton, Marco Quaggiotto, Wouter Van den Broeck,
Corinne Régis, Bruno Lina, and Philippe Vanhems.
High-Resolution Measurements of Face-to-Face Contact Patterns in a
Primary School.
PLoS ONE, 60 (8):0 e23176, 2011.

[Tönshoff et al., 2023]
Jan Tönshoff, Martin Ritzert, Hinrikus Wolf, and Martin Grohe.
Walking Out of the Weisfeiler Leman Hierarchy: Graph Learning Beyond
Message Passing.
Transactions on Machine Learning Research, 2023.
ISSN 2835-8856.
URL <https://openreview.net/forum?id=vgXnEyeWVY>.

[Xu et al., 2019]
Keyulu Xu, Weihua Hu, Jure Leskovec, and Stefanie Jegelka.
How Powerful are Graph Neural Networks?
In International Conference on Learning Representations, 2019.
URL <https://openreview.net/forum?id=ryGs6iA5Km>.

[Yadati et al., 2019]
Naganand Yadati, Madhav Nimishakavi, Prateek Yadav, Vikram Nitin, Anand Louis,
and Partha Talukdar.
HyperGCN: A New Method For Training Graph Convolutional Networks on
In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett (eds.), Advances in Neural
Information Processing Systems, volume 32. Curran Associates, Inc., 2019.

[Yang et al., 2022]
Maosheng Yang, Elvin Isufi, and Geert Leus.
Simplicial Convolutional Neural Networks.
In ICASSP 2022 - 2022 IEEE International Conference on
Acoustics, Speech and Signal Processing (ICASSP), pp.  8847–8851. IEEE,

[Yang et al., 2022]
Ruochen Yang, Frederic Sala, and Paul Bogdan.
Efficient Representation Learning for Higher-Order Data With
Simplicial Complexes.
In Bastian Rieck and Razvan Pascanu (eds.), Proceedings of the
First Learning on Graphs Conference, volume 198 of Proceedings of
Machine Learning Research, pp.  13:1–13:21. PMLR, 09–12 Dec
URL <https://proceedings.mlr.press/v198/yang22a.html>.

[Zhou et al., 2023]
Cai Zhou, Xiyuan Wang, and Muhan Zhang.
Facilitating Graph Neural Networks with Random Walk on Simplicial
In Thirty-seventh Conference on Neural Information Processing
Systems, 2023.
URL <https://openreview.net/forum?id=H57w5EOj6O>.

§ TRAINING DETAILS AND RESULTS

Imputation accuracies on the Semantic Scholar dataset with varying percentages of missing citations.
We report the mean and standard deviation of the validation accuracies over $10$ runs.
We compare SCRaWl to the simplicial neural networks SNN, SCNN, SAT, SAN, and MPSN as well as to the representation learning algorithm $k$-simplex2vec.

Vertex classification accuracies on the social contact datasets.

For all experiments, we use the Adam optimizer with an initial learning rate of $10^{-3}$.
The learning rate is reduced by a factor of $0.5$ if the validation loss does not improve for $10$ epochs.
Training is stopped once the learning rate drops below $10^{-6}$.

The walk length $\ell$ is the primary hyperparameter that determines the receptive field of the model.
As the citation counts in the Semantic Scholar dataset can be imputed well with only local information, we choose a walk length of $5$ for this dataset.
On other datasets, the model is trained with a walk length of $50$.
The walk length can be increased independently in prediction, but we did not make use of this option in our experiments.

For the Semantic Scholar dataset, we ran SCRaWl modules on simplex orders $k \in \{ 0, \dots, 5 \}$, i.e., simplices of order $6$ appear only as static coface features in the random walks and simplices of order $7$ do not influence the model.
Using uniform connection sampling, we computed $m = \sum_{i=0}^n n_i$ random walks, i.e., one random walk for each simplex of order $k \in \{ 0, \dots, 5 \}$.
While SCRaWl modules can be configured individually for each simplex order and each layer, we found that the same configuration for all modules works well for this dataset:
Each module is configured with a local window size of $s = 4$, a kernel size of $d_\text{kern} = 8$, a hidden feature size of $d = 32$, and a mean pooling operation.
Complete results for the semantic scholar experiments are given in <ref>.

For the social contact datasets, we use $4$ layers and run SCRaWl modules on simplex orders $k \in \{ 0, \dots, 3 \}$.
We compute one random walk for each simplex.
Each SCRaWl module is again configured identically with a local window size of $s = 8$, a kernel size of $d_\text{kern} = 8$, a hidden feature size of $d = 128$, and a mean pooling operation.
The final vertex embeddings are fed into a $2$-layer MLP with ReLU activation.
Complete results for the social contact experiments are given in <ref>.

A full list of hyperparameters used for each experiment is given in <ref>.
The model has been implemented using PyTorch [Paszke et al., 2019] and TopoX [Hajij et al., 2024].
We adapted the code from Tönshoff et al., 2023 to implement SCRaWl.

List of hyperparameters and training setup used for each experiment.

Parameter     Semantic Scholar     Social Contacts

Walk Length $\ell$     $5$     $50$

Sampling Method     uniform connection     uniform connection

Number of Layers $L$     $3$     $4$

Max. Simplex Order $K$     $5$     -

Local Window Size $s$     $4$     $8$

Kernel Size $d_\text{kern}$     $8$     $8$

Hidden Feature Size $d$     $32$     $32$

Pooling     mean     mean

1-1Total number of
Trainable Parameters     171k     132k

Optimizer     1-2Adam; $\text{LR} = 10^{-3}$

LR Scheduler     1-2Reduce on plateau: factor $0.5$; patience $10$

Stopping Criterion     1-2$\text{LR} < 10^{-6}$

We used Bodnar et al., 2021 reference implementation for MPSN and conducted the experiments using the author's default parameters with $4$ layers, $128$ hidden features, ReLU nonlinearities, and sum aggregations.
Experiments for SCNN were conducted using the implementation provided by TopoX with $3$ layers, $64$ hidden features, and ReLU nonlinearities.
We used our own implementation of $k$-simplex2vec.

§ ABLATION STUDY

Ablation study on the social contact datasets.
Left: Influence of upper and lower adjacency on the performance of SCRaWl.
Middle: Influence of the local window size $s$.
Right: Influence of the walk length $\ell$.
We report the mean and standard deviation of the validation accuracies over $5$ runs.

In this ablation study, we investigate the influence of different hyperparameters and other aspects on the performance of SCRaWl.
For all experiments, the training setup and hyperparameters not under investigation are set to the same values as in the main experiments.
<Ref> shows the results of the ablation study on the social contact datasets.

Lower vs. Upper Adjacency

During a random walk, the next simplex can either be reached by traversing over a face or a coface of the current simplex.
We are interested in the impact of either form of adjacency, i.e., we compare the performance of SCRaWl when using only lower adjacency or only upper adjacency.
For simplices that don't have a face or don't have a coface, respectively, the random walk stays at the current simplex.
We make the exception that for node-level walks, upper adjacency over edges is always allowed, as otherwise, no meaningful walks would be possible.

For the primary school dataset, we see that only using upper adjacency yields a significantly better performance of $0.889 \pm 0.063$ than only using lower adjacency with $0.454 \pm 0.099$.
In this case, only using upper adjacency is thus only marginally below the baseline performance of $0.927 \pm 0.026$.
On the easier-to-learn high school dataset, the difference between using only upper or lower adjacency is smaller with $0.994 \pm 0.009$ vs. $0.975 \pm 0.009$, compared to the baseline with $1.0 \pm 0.0$.

These results support our claim that — at least for these datasets — higher-order interactions give additional insight into the structure of the dataset and methods that capture these (higher-order) structures can perform better.

Influence of Local Window Size $s$

The local window $s$ is a parameter that influences which structural features of the input simplicial complex are visible to the model.
We thus expect that larger windows improve the prediction performance while also increasing the computational complexity.

We validate this hypothesis by training SCRaWl on the social contact datasets with different local window sizes $s \in \{ 1, \dots, 8 \}$.
As expected, we see that the prediction performance of SCRaWl improves with larger windows:
For the primary school dataset, the performance increases almost monotonically from $0.571 \pm 0.079$ with $s = 1$ to $0.927 \pm 0.026$ with $s = 8$.
For the easier-to-learn high school dataset the performance is only marginally worse with $0.977 \pm 0.006$ for $s = 1$ vs. $1.0 \pm 0.0$ for $s = 8$.

Influence of Walk Length $\ell$

Similar to the local window size $s$, we expect that larger walk lengths $\ell$ positively influence the prediction performance of SCRaWl, as longer walks ensure that important structural features of the simplicial complex are captured and processed by the model.

To that end, we trained SCRaWl on the social contact datasets with different walk lengths $\ell \in \{ 5, 10, \dots, 50 \}$.
The prediction performance again improved almost monotonically with increasing walk lengths for both datasets, with the high school dataset reaching accuracies above $0.98$ for $\ell \geq 25$ already.
For the primary school dataset, the performance reached its best accuracy of $0.927 \pm 0.026$ only for $\ell = 50$.

We note that the choice of walk length $\ell$ and window size $s$ has a considerable impact on the computational cost of the model and its value should be chosen carefully to balance predictive performance and computational complexity.


# Supervised deep learning of elastic SRV distances on the shape space of
curves

Emmanuel Hartman1 Yashil Sukurdeep2 Nicolas Charon2 Eric Klassen1 Martin
Bauer1
Department of Mathematics M. Bauer and E. Hartman were supported by NSF grants
1953244 and 1912037. Y. Sukurdeep and N. Charon were supported by NSF grants
1945224 and 1953267.  Florida State University  Tallahassee  USA1
Department of Applied Mathematics  Johns Hopkins University  Baltimore  USA2
<EMAIL_ADDRESS>yashil.sukurdeep<EMAIL_ADDRESS>
###### Abstract

Motivated by applications from computer vision to bioinformatics, the field of
shape analysis deals with problems where one wants to analyze geometric
objects, such as curves, while ignoring actions that preserve their shape,
such as translations, rotations, scalings, or reparametrizations. Mathematical
tools have been developed to define notions of distances, averages, and
optimal deformations for geometric objects. One such framework, which has
proven to be successful in many applications, is based on the square root
velocity (SRV) transform, which allows one to define a computable distance
between spatial curves regardless of how they are parametrized. This paper
introduces a supervised deep learning framework for the direct computation of
SRV distances between curves, which usually requires an optimization over the
group of reparametrizations that act on the curves. The benefits of our
approach in terms of computational speed and accuracy are illustrated via
several numerical experiments on both synthetic and real data.

##  1 Introduction

Motivated by applications from computer vision to bioinformatics, the field of
elastic shape analysis deals with problems where one needs to analyze the
variability of geometric objects [18], [1], [19], [13], [5]. In this article,
we address the computation of elastic geodesic distances between geometric
curves in one and higher dimensions. By geometric curves, we mean curves
modulo shape-preserving transformations, i.e., curves whose images are equal
up to translations, rotations, scalings and reparametrizations.
Mathematically, we model the space of all geometric curves as a quotient space
of the set of absolutely continuous curves. In Section 2, we discuss the
construction of this space, and how the square root velocity (SRV) transform
allows us to define the so-called SRV distance between geometric curves. For
additional details, we refer interested readers to the the vast literature on
elastic metrics [18], [14], [20], [17], [2].

###  1.1 Background & Related Work

The SRV distance [18] quantifies dissimilarity between geometric curves, and
can be used for averaging, classifying and clustering datasets of functions or
curves, which are prevalent tasks in fields such as computer vision and
medical imaging. However, computing SRV distances between pairs of geometric
curves is a nontrivial endeavor, as it typically necessitates solving an
infinite dimensional optimization problem over the set of reparametrizations
that act on the curves. These optimal reparametrizations are guaranteed to
exist for certain classes of curves, such as $C^{1}$-curves [4], and piecewise
linear curves [11]. The beauty of the result in [11] lies in the fact that it
not only provides an existence result, but also describes an algorithm to
explicitly construct optimal reparametrizations for piecewise linear curves.
Although this algorithm allows one to compute exact SRV distances, it has a
high polynomial complexity, rendering it impractical for large datasets that
are typically encountered in applications.

Consequently, faster approaches have emerged to compute SRV distances in
practical contexts [8], [18]. Several such algorithms rely on dynamic
programming (DP) to approximate optimal reparametrizations, operating by
searching over subsets of all possible reparametrizations. The different DP-
based approaches achieve runtimes of $O(n)$ to $O(n^{3})$, where $n$ is the
number of samples used to discretize the curves, thus providing fast over-
estimates of the true SRV distance [3], [6], [18]. Nevertheless, DP also
incurs a significant computational cost when working with very large datasets,
implying that there is a need to develop increasingly efficient approaches for
handling modern datasets of shapes. Towards that end, deep learning (DL)
approaches have recently been introduced to estimate optimal
reparametrizations for functions and curves, including supervised [12], [15]
and unsupervised [10], [16] methods. In particular, Nunez and Joshi train a
convolutional neural network (CNN) on approximate reparametrizations obtained
from DP in order to estimate optimal reparametrizations for functions and
curves [15].

###  1.2 Contributions

In this paper, we propose a supervised DL framework for directly estimating
SRV distances between functions and between spatial curves, without the need
to estimate optimal reparametrizations. More specifically, we train a deep CNN
to learn SRV distances, using training data consisting of pairs of discretized
functions or curves, together with the SRV distance between them as labels.

As a theoretical contribution that is of interest on its own, we extend the
existence results for optimal reparametrizations of [11, 4] to the space of
closed curves, and to the spaces of (open or closed) curves modulo rotations.
These results were previously only known for open curves, and also did not
consider the action of the rotation group. This in turn allows us to directly
generalize the algorithm of [11] for calculating exact SRV distances for open
or closed curves modulo rotations.

Consequently, and in contrast to e.g. [15], we use these exact SRV distances
as training labels for our network, rather than SRV distance over-estimates
computed via DP. This reduces bias in the network’s predictions. Another
distinct feature of our framework is that unlike the aforementioned DP and DL
approaches, we bypass the need to estimate optimal reparametrizations, instead
directly estimating SRV distances. This is especially convenient for certain
unsupervised learning tasks involving datasets of shapes, such as clustering
applications with curves, where one only needs rapidly-computed pairwise
distances rather than optimal reparametrization maps between the shapes.

Moreover, using a neural network to predict a single number (the SRV distance)
rather than a full reparametrization map is an obviously less complex learning
problem, and can thus be achieved with a smaller training set, which is
advantageous in the context of shape analysis where vast amounts of publicly
available data are not as readily available as in e.g., the imaging sciences.
Furthermore, a fundamental property of the SRV distance is its invariance to
parametrizations, which we leverage to introduce a shape-preserving data
augmentation training strategy, outlined in Section 3. This training strategy
allows us to augment the size of our training set while also improving the
variability within training samples, which ultimately allows the trained
network to produce robust, parametrization-invariant SRV distance estimates.

We are also providing an open-source version of the code for our DL framework,
which is publicly available on github at https://github.com/emmanuel-
hartman/supervisedDL-SRVFdistances.

To illustrate our DL framework’s benefits, we show that our trained CNN’s SRV
distance estimates are comparable to or even more accurate than DP distances,
while also being orders of magnitude faster in terms of computation time.

##  2 The Shape Space and SRV Distance

We begin with a brief overview of the square root velocity (SRV) framework for
defining a computable elastic distance between geometric curves. In this
framework, we start with parametrized curves, modelled as elements of the
space of vector-valued absolutely continuous functions, denoted
$\operatorname{AC}(M,\mathbb{R}^{d})$. We have open curves if the parameter
space $M$ is the unit interval $[0,1]\subset{\mathbb{R}}$, and closed curves
if $M$ is the unit circle $S^{1}$. When $d=1$, we have one-dimensional curves,
which we call functions.

In shape analysis for curves, one is interested in the space of all geometric
curves, i.e., parametrized curves whose images are equal up to translations,
rotations and reparametrizations. Note that the SRV framework, and thus our DL
approach, can easily be extended to handle curves modulo scalings as well. We
now briefly outline the construction of the space of all geometric curves.

To identify parametrized curves that only differ by a translation, we work
with the space of absolutely continuous curves such that $c(0)=0$, denoted by
$\operatorname{AC}_{0}(M,\mathbb{R}^{d})$. We will later see that the SRV
distance is naturally defined on this linear subspace of all absolutely
continuous curves. Identifying curves that only differ by a reparametrization
or rotation is a more delicate matter and the main source of complication in
shape analysis. This is accomplished by defining the following equivalence
relation for parametrized curves
$c_{1},c_{2}\in\operatorname{AC}_{0}(M,\mathbb{R}^{d})$. For open curves, we
define $c_{1}\thicksim c_{2}$ if and only if they have the same unit speed
parametrization after the application of an appropriate rotation. For closed
curves, we define $c_{1}\thicksim c_{2}$ if and only if they have the same
unit speed parametrization after an appropriate choice of “starting point” and
an appropriate rotation.111For functions, i.e., $d=1$, the rotation group is
trivial and thus the equivalence relation reduces to factoring out the
reparametrization action only. We denote the equivalence class of a curve $c$
under this relation by $[c]$. We then define the space of all geometric curves
as the quotient space:

$\displaystyle\mathcal{S}(M,{\mathbb{R}}^{d})=\operatorname{AC}_{0}(M,\mathbb{R}^{d})/\thicksim~{},$

and will refer to it as the shape space of curves, or simply as the shape
space for brevity when there is no ambiguity.

We now outline how the SRV transform allows us to define a distance function
on this shape space. The SRV transform is the mapping
$Q:\operatorname{AC}_{0}(M,\mathbb{R}^{d})\rightarrow
L^{2}(M,\mathbb{R}^{d})$, defined by:

$\displaystyle c(\cdot)\mapsto
Q(c)(\cdot):=\begin{cases}\frac{c^{\prime}(\cdot)}{\sqrt{|c^{\prime}(\cdot)|}}&\text{if
}|c^{\prime}(\cdot)|>0,\\\ 0&\text{otherwise.}\end{cases}$

Here, $c^{\prime}$ denotes the first derivative of the parametrized curve
$c\in\operatorname{AC}_{0}(M,\mathbb{R}^{d})$. This transform allows us to
define the SRV distance between parametrized curves
$c_{1},c_{2}\in\operatorname{AC}_{0}(M,\mathbb{R}^{d})$ by pulling back the
$L^{2}$ metric on $L^{2}(M,\mathbb{R}^{d})$ as follows:

$\displaystyle d_{Q}(c_{1},c_{2})^{2}$
$\displaystyle:=\|Q(c_{1})-Q(c_{2})\|_{L^{2}}^{2}$
$\displaystyle=\int_{M}\Bigg{|}\frac{c_{1}^{\prime}(t)}{\sqrt{|c_{1}^{\prime}(t)|}}-\frac{c_{2}^{\prime}(t)}{\sqrt{|c_{2}^{\prime}(t)|}}\Bigg{|}^{2}dt.$

It is worth noting that in the case of open curves, this distance can be
interpreted as the geodesic distance induced by a Riemannian metric. For
closed curves, it is only a first order approximation of a geodesic distance.
The key property of this distance is its invariance under both the action of
the group of rotations $\operatorname{SO}(d)$, and that of the group of
diffeomorphisms of the parameter space ${\operatorname{Diff}}(M)$. The latter
can be seen by a simple change of variables in the above integral. Thus, the
SRV distance descends to a distance on the quotient shape space
$\mathcal{S}(M,{\mathbb{R}}^{d})$, given by:

$\displaystyle
d_{\mathcal{S}}([c_{1}],[c_{2}])=\inf\limits_{\begin{subarray}{c}\gamma\in{\operatorname{Diff}}(M)\\\
\operatorname{O}\in\operatorname{SO}(d)\end{subarray}}d_{Q}\left(c_{1},\operatorname{O}\star\big{(}c_{2}\circ\gamma\big{)}\right).$
(1)

With a slight abuse of terminology, we henceforth refer to the quotient space
distance, namely $d_{\mathcal{S}}(\cdot,\cdot)$, as the SRV distance. It
follows that computing the SRV distance between geometric curves involves
solving a joint optimization problem over the finite dimensional group
$\operatorname{SO}(d)$ and the infinite dimensional reparametrization group
${\operatorname{Diff}}(M)$. The main challenge is the minimization over
${\operatorname{Diff}}(M)$, which is usually accomplished by discretizing the
group into a finite dimensional approximation space, and solving the
discretized problem via a dynamic programming approach.

It is important to note that in general, the existence of a reparametrization
$\gamma\in{\operatorname{Diff}}(M)$ attaining the infimum in (1) is not
guaranteed. However, under some additional regularity assumptions on the
curves $c_{1},c_{2}$, one can recover such existence results. In the
following, we discuss the existence of optimal reparametrizations and
rotations in (1), both for the case of open curves (i.e., $M=[0,1]$) as well
as closed curves (i.e., $M=S^{1}$).

We first introduce the semi-group of generalized reparametrizations for open
curves:

$\bar{\Gamma}([0,1])=\left\\{\gamma\in\operatorname{AC}([0,1],[0,1]):\gamma\text{
is onto; }\gamma^{\prime}\geq 0\text{ a.e.}\right\\}.$

To introduce the analogous construction for the case of closed curves, we view
$S^{1}$ as $\mathbb{R}/\mathbb{Z}$. We then define the shift operator on
$S^{1}$ via:

$S_{\theta}:S^{1}\to S^{1};\qquad\lambda\mapsto\lambda+\theta.$

This allows us to define the semi-group of generalized reparametrizations on
$S^{1}$ via:

$\bar{\Gamma}(S^{1})=\left\\{S_{\theta}\circ\gamma^{*}:\theta\in S^{1}\text{
and }\gamma^{*}\in\bar{\Gamma}([0,1])\right\\}.$

This allows us to formulate the following existence result, which is the main
theoretical contribution of the present article:

######  Theorem 1.

Let $c_{1},c_{2}\in\operatorname{AC}(M,\mathbb{R}^{d})$ such that either both
are of class $C^{1}$, or at least one of them is piecewise linear. Assume also
that $c_{1}^{\prime}$ and $c_{2}^{\prime}$ are both nonzero a.e. on $M$. Then
there exists a pair of generalized reparametrization functions
$(\gamma_{1},\gamma_{2})\in\bar{\Gamma}(M)\times\bar{\Gamma}(M)$ and a
rotation $O\in\operatorname{SO}(d)$ such that:

$d_{\mathcal{S}}([c_{1}],[c_{2}])=d_{Q}(c_{1}\circ\gamma_{1},O\star(c_{2}\circ\gamma_{2})).$

Previously this result was only known for the space of open curves and did not
consider the action of the rotation group, see [11, 4]. The proof of Theorem
1, which builds up on these results, is postponed to the appendix.

For piecewise linear curves
$c_{1},c_{2}\in\operatorname{AC}_{0}(M,\mathbb{R}^{d})$, the results of [11]
and Theorem 1 even lead to an algorithm that allows us to explicitly construct
these optimal reparametrizations for calculating the precise quotient space
distance, see [11]. This algorithm plays a fundamental role in our proposed DL
framework, as we use it to calculate exact quotient SRV distances in order to
generate labels for our training data, as will be outlined in the next
section.

##  3 Deep Learning of SRV Distances

While the algorithm in [11] allows us to compute exact SRV distances, it is
computationally expensive, making it impractical for working with large
datasets of shapes. Consequently, there is a need to develop approaches that
are more computationally efficient in order to calculate SRV distances. We
address this need by introducing a supervised DL framework that provides fast,
accurate and robust SRV distance estimates.

###  3.1 Network Architecture

Figure 1: Training step and network structure diagram for Shape Preserving
Data Augmentation based training: Weights contained in the red blocks are
trainable and the Siamese convolutional nodes have shared weights. Specific
parameter details of the network architecture can be found in Section 1. The
green blocks perform shape preserving data augmentation as described in
Section 3.2.

We train a Siamese convolutional neural network (CNN) to learn the SRV
distance between geometric curves. We use training data consisting of pairs of
discretized ${\mathbb{R}}^{d}$ valued curves, together with their SRV distance
as labels. Each individual curve is sampled at $n$ vertices and represented as
a flattened vector of length $n\times d$, before being fed as input to the
network. Our Siamese CNN has a twin structure, consisting of two components
which have identical architectures and use the same weights. To be more
specific, each component of the CNN operates on an individual discretized
curve, which is passed through a series of convolutional layers with kernels
of size 5, followed each time by batch normalization, a rectified linear unit
(ReLU) activation, and a max-pooling layer with pool size $2$. This produces
two outputs, which are concatenated and passed through four dense layers whose
widths are proportional to $d$, with ReLU activations being used in each dense
layer. The network then outputs a single real number: the SRV distance between
the two curves. We provide a schematic description of the network architecture
in Fig. 1.

###  3.2 Training Method

Figure 2: Example of Shape Preserving Data Augmentation: The top curve is an
example of a parameterization of a curve from the Swedish Leaf II dataset, see
Section 4 for a description of this dataset. The three curves on the bottom
represent parameterizations and rotations of this curve as produced by the
shape preserving data augmentation described in Section 3.2.

We create training and testing sets for our network by randomly generating
pairs of functions or curves, or by picking them from an existing dataset, and
labelling them with their SRV distance. We use exact distances computed with
the algorithm of [11] as labels for functions and 2D curves, but due to this
algorithm’s high complexity, we instead use DP distances as labels for 3D
curves.

Computing these SRV distance labels using the exact algorithm or DP may be
very time consuming, which could limit the size of our training set in
practice. Thankfully, from a base training set, one can easily generate more
training samples at no extra cost by applying shape-preserving
transformations, such as resampling and rotations, to both curves. Indeed, the
quotient SRV distance is invariant to reparametrizations (i.e., to resampling
in the discrete situation) and to rotations, implying that the distance
between resampled and/or rotated curves remains unchanged and need not be
recomputed. This data augmentation strategy allows the network to see a wider
variety of sampling patterns and rotations for the same curve during training,
which helps it to learn and predict distances that are truly invariant to
reparametrizations and rotations. Moreover, since resampling and rotating
curves is computationally inexpensive, this procedure can be performed at each
iteration of the training step, without incurring any additional storage for
new training samples. We empirically observed that this shape-preserving data
augmentation-based training method reduced overfitting in the distance learned
by the network.

The training itself is performed using an Adam optimization procedure [9]. We
observed relatively fast convergence in all cases, with convergence curves
shown in Fig. 3. We refer readers to the code documentation on
https://github.com/emmanuel-hartman/supervisedDL-SRVFdistances for further
training details, including information on the exact training parameters such
as the batchsize for each epoch of training, and the parameters involved in
the shape-preserving data augmentation training step.

Figure 3: On both figures, the $x$-axis represents epochs, and on the
$y$-axis, we plot the mean squared error of the network on the training data
(blue), as well as on unseen testing data (red). Convergence curves for
network trained on open, real-valued functions discretized at 90 points from
our Synthetic I data set, trained for 500 epochs (left figure). Convergence
curves for network trained on closed, 2-dimensional curves discretized at 100
points from the Kimia dataset, trained for 50 epochs (right figure).
Descriptions of the datasets are given in Section 4.

##  4 Numerical Experiments

We now present empirical results demonstrating the performance of our DL
approach for estimating SRV distances on real-valued functions, and on curves
in ${\mathbb{R}}^{2}$ and ${\mathbb{R}}^{3}$. As we shall see, the experiments
show that when compared to DP, our approach produces SRV distance estimates at
a significantly lower numerical cost, while being comparable, and sometimes
superior, in terms of accuracy.

###  4.1 Computation Method

To compute SRV distances using the exact algorithm and DP, we used Martins
Bruveris’ package222https://github.com/martinsbruveris/libsrvf, which builds
on the DP code of FSU’s Statistical Shape Analysis and Modeling Group. Our
network was implemented on TensorFlow. All computation times using the
different algorithms were recorded on an Intel Xeon X5650 2.66 GHz CPU with a
Gigabyte GeForce GTX 1060 1582 MHz GPU. A comparision of the computation times
for the various algorithms can be found in Table 1. One can clearly see that
the trained network is several orders of magnitude faster than both DP and the
exact algorithm.

| Exact | DP | DL (CPU) | DL (GPU)
---|---|---|---|---
1D | $2\times 10^{3}$ | $5\times 10^{2}$ | $5\times 10^{-2}$ | $2\times 10^{-2}$
2D | $2\times 10^{5}$ | $8\times 10^{2}$ | $2\times 10^{-1}$ | $3\times 10^{-2}$
Table 1: Computation time for one SRV distance using several different
algorithms, in milliseconds.

###  4.2 Evaluation Method

For functions and 2D curves: To evaluate the trained network’s accuracy for
functions and 2D curves, we use the mean relative error (MRE) between its
output and the true SRV distances on a test set, computed via the exact
algorithm. As a secondary measure of estimation quality, we use the Pearson
correlation coefficient

$\rho_{y\widehat{y}}=\frac{\sum\limits_{i=1}^{N}(y_{i}-\overline{y})(\widehat{y_{i}}-\overline{\widehat{y}})}{\sqrt{\sum\limits_{i=1}^{N}(y_{i}-\overline{y})^{2}}\sqrt{\sum\limits_{i=1}^{N}(\widehat{y_{i}}-\overline{\widehat{y}})^{2}}}$

between the network’s output and exact distances on a test set of functions or
2D curves. Here $N$ is the number of training samples, $\\{y_{i}\\}_{i=1}^{N}$
are the exact distances, $\\{\widehat{y_{i}}\\}_{i=1}^{N}$ are the outputs of
the network, with $\overline{y}=\frac{1}{N}\sum_{i=1}^{N}y_{i}$ and
$\overline{\widehat{y}}=\frac{1}{N}\sum_{i=1}^{N}\widehat{y_{i}}$ being their
respective sample means. A low MRE and a strong positive correlation
coefficient indicates a good performance of the network.

For 3D curves: In addition to the results for functions and planar curves, we
present preliminary results for curves in $\mathbb{R}^{3}$. As the
computational complexity of the exact algorithm is orders of magnitude higher
for curves in $\mathbb{R}^{3}$ when compared to the case of functions and
curves in ${\mathbb{R}}^{2}$, we only label 3D curves with DP distances.
Consequently we can only evaluate the CNN’s performance for 3D curves via the
correlation coefficient between its output and DP distances.

To avoid bias in our results, elements of the test set are never contained in
the training set, which is used solely for the purpose of calibrating the
network. Furthermore, to assess our network’s generalization capabilities, we
make sure that the trained network is tested on data that is significantly
different compared to the data used for training, c.f. Figures 4 and 5.

###  4.3 Experiments with functions

Datasets: First, we tested our network’s ability to predict SRV distances for
functions, using both synthetic and real data. The synthetic data was created
by generating functions sampled at $90$ evenly spaced points on the unit
interval with random arc length. We note that the shape class of a function
modulo reparametrizations is entirely determined by its local maxima and
minima, as it is determined by its constant speed parameterization which is a
linear interpolation between the local maxima and minima. Thus, the synthetic
function data is generated by drawing the number of extrema for our function
from a $\mathcal{N}(\mu,\sigma)$ distribution, randomly assigning values for
these extrema, and then randomly choosing a function with $90$ breakpoints
from the shape class determined by the generated extrema. We created two
different synthetic datasets: the first one with parameters
$(\mu,\sigma)=(18,6)$, and the second with $(\mu,\sigma)=(30,10)$. These
datasets, dubbed Synthetic I and Synthetic II respectively, each contain
100,000 pairs of functions labelled with their exact distances, partitioned
into 99,000 training cases and 1,000 testing cases. For the real dataset, we
use CPC Global Unified Precipitation data from the NOAA/OAR/ESRL PSL, Boulder,
Colorado, USA333https://psl.noaa.gov/, from which we extracted 90 days of
precipitation data across several locations and years. We randomly selected
400 samples from this database, computed exact pairwise distances, and
partitioned them into a set of 89700 distances for training, and 9,900 for
testing. See Fig. 4 for examples from the different datasets.

Figure 4: Five examples from Synthetic I (top-left) and the CPC Precipitation
dataset (top-right). Third and fourth figure: Comparison of DP (red) and
Trained Network (blue). Scatter plot of relative errors for 1000 testing cases
from the CPC Precipitation dataset, using a network trained on Synthetic I
(bottom-left). Corresponding correlation plot for both methods, with exact
distances on the $y$-axis, and estimated distances on the $x$-axis, and the
line $y=x$ in green (bottom-right).

Results: First, we highlight the difference in performance between our trained
network and DP, see Table 2 and Fig. 4. When trained and validated on the same
type of data, the network significantly outperforms DP across all three
datasets, both in terms of the MRE and correlation coefficient with respect to
the exact distances, see Table 2.

To demonstrate our network’s generalization capabilities, we trained it on one
type of data and tested it on a different dataset, e.g., we trained on
synthetic data but tested on CPC precipitation data. While this leads to a
slight increase in prediction error, the network still outperforms DP on both
measures by a large margin, see Table 3.

Dataset | MRE | Corr.
---|---|---
DP | DL | DP | DL
Synthetic I | 0.44551 | 0.04346 | 0.84307 | 0.96770
Synthetic II | 0.45949 | 0.03806 | 0.87123 | 0.97053
CPC Precip. | 0.45853 | 0.03722 | 0.85452 | 0.96090
Table 2: Comparison between DP and DL Training Set | Testing Set | MRE | Corr.
---|---|---|---
Synthetic I | Synthetic II | 0.05520 | 0.96782
Synthetic I | CPC Precip. | 0.08088 | 0.95264
Synthetic II | Synthetic I | 0.05206 | 0.96110
Synthetic II | CPC Precip. | 0.07093 | 0.94888
Table 3: Generalization results across several testing sets

###  4.4 Experiments with curves in ${\mathbb{R}}^{2}$

Datasets: We used data from the
MPEG-7444https://dabi.temple.edu/external/shape/MPEG7/dataset.html and Swedish
leaf datasets555https://www.cvl.isy.liu.se/en/research/datasets/swedish-leaf/,
which contain images of objects whose boundaries were extracted and treated as
2D curves, discretized with $100$ points, see Fig. 5. To extract discretized
boundary curves from these datasets, we binarized each image via Otsu’s
algorithm, then extracted vertices on the boundary using the Moore-Neighbor
tracing algorithm, before downsampling to $100$ points.

We trained the network on 229162 distinct pairs of curves labelled with exact
distances from the MPEG-7 dataset, which contains a diverse array of 2D shapes
from many different shape classes, see Fig. 5. We tested the network on two
versions of the Swedish leaf dataset, called Swedish Leaf I and II
respectively. Swedish Leaf I contains curves with arc length parametrizations,
i.e., discretized with $n$ points that are uniformly distributed across the
curve. Swedish Leaf II contains “adversarial parametrizations”, i.e., curves
with $n$ points that are far from uniformly distributed across the curve, with
many points concentrated on a small portion of the curve, see Fig. 6.

Figure 5: Five examples from the MPEG-7 dataset (left). Five examples from the
Swedish leaf dataset (right).

Figure 6: Example of a curve from the Swedish Leaf I dataset, where curves
have arc length parametrization (left), and from the Swedish Leaf II dataset
with adversarial parametrizations (right).

Results: Due to the higher dimensionality and complexity of the data, the
network’s performance drops compared to the case of functions. However, we
still obtain a high correlation coefficient with the exact distance across
both datasets, namely 0.924 for Swedish Leaf I, and 0.917 for Swedish Leaf II.
Meanwhile, the corresponding correlation coefficient for DP distances is 0.996
for Swedish Leaf I, but drops significantly to 0.899 for Swedish Leaf II.
These observations show comparable performance between our DL framework and
DP, with DP being more accurate for curves that are already well-aligned
(e.g., for those in Swedish Leaf I), and DL being superior in terms of
accuracy for data requiring larger reparametrizations (e.g., for curves in
Swedish Leaf II).

As yet another proof of concept for our DL approach, we perform an
unsupervised clustering experiment using $40$ curves taken from the Swedish
Leaf I dataset. These curves are evenly distributed across four categories of
leaves. To cluster the curves, we compute all pairwise SRV distances using
both our DL framework and the exact algorithm, and apply classical
multidimensional scaling (CMDS) to the resulting pairwise distance matrices in
order to obtain a 2D projection of the dataset, see Fig. 7. While the
resulting 2D visualizations are slightly different, the clusters produced are
comparable.

Figure 7: CMDS clusters of 40 curves selected from the Swedish leaf dataset
using exact distances (left) and DL distances (right).

###  4.5 Preliminary experiments for curves in ${\mathbb{R}}^{3}$

Datasets: Finally we present preliminary results for curves in
${\mathbb{R}}^{3}$. We use two distinct datasets of open 3D curves for the
experiments: hurricane paths from the National Hurricane Center Data
Archive666https://www.nhc.noaa.gov/data/ and plant
roots777https://github.com/RSA-benchmarks/collaborative-comparison from which
we only keep the taproots (i.e. main stem of the roots). All curves were
discretized with $100$ points. We trained the network on 284622 distinct pairs
of hurricane paths, labelled with distances computed with the DP algorithm. We
validated the network on either a different testing set of hurricane paths, or
on the dataset of taproots. The reason for choosing DP distances instead of
exact distances for training the network is the high computational cost of
calculating exact distances in the situation of 3D curves.

Results: The correlation coefficient between the predicted DL and the DP
distances on the test set of hurricane paths is 0.977, but drops significantly
to 0.823 when tested on taproots. This drop-off in prediction quality can be
most likely explained by the lack of sufficiently diverse samples in the
training set, which limits the network’s generalization capabilities. We
expect that enriching the training set with 3D curves displaying more varied
geometries will help to improve the network’s performance. However, due to the
scarcity of publicly available datasets of 3D curves, we leave it as future
work to build a better training set for 3D curves and invest the required
computational resources to train our network using SRV distance labels
computed via the exact algorithm, as was done for functions and 2D curves.

##  5 Conclusion

We have introduced a supervised DL framework to compute SRV distances for
curves in ${\mathbb{R}}^{d}$. The main advantage of our approach is that our
trained CNN provides fast and accurate estimates of the SRV distance between
pairs of geometric curves, without the need to find optimal
reparametrizations. Moreover, we exploited the invariance of the SRV distance
to shape-preserving actions in order to propose a shape-preserving data-
augmentation based training strategy, which is a flexible and efficient
procedure for creating and augmenting our training set. Empirical observations
show that this training strategy allows our network to estimate SRV distances
which are truly invariant to rotations and reparametrizations, while also
reducing overfitting. Moreover, our experiments show that when compared to DP,
our approach produces SRV distance estimates at a significantly lower
numerical cost, while also being comparable, and sometimes superior, in terms
of accuracy.

### Acknowledgements

We want to thank Martins Bruveris for helpful discussions during the
preparation of this manuscript.

## Appendix

###### Proof of Theorem 1.

Not taking into account the action of the group of rotations and considering
only open curves, this result was shown in [11], assuming that one of the
curves is piecewise linear, and in [4] under the assumptions that both curves
are of class $C^{1}$.

In the case where $M=S^{1}$, by the definition of $\bar{\Gamma}(S^{1})$, the
existence of a pair of optimal
$(\gamma_{1},\gamma_{2})\in\bar{\Gamma}(S^{1})\times\bar{\Gamma}(S^{1})$ is
equivalent to the existence of an optimal $\tau\in S^{1}$, an optimal
$O\in\operatorname{SO}(d)$ and a pair of optimal
$\gamma^{*}_{1},\gamma^{*}_{2}\in\overline{\Gamma}([0,1])$. Consider the
function $F:S^{1}\times\operatorname{SO}(d)\to{\mathbb{R}}$ given by

$F(\lambda,O)=\inf\limits_{\gamma^{*}_{1},\gamma^{*}_{2}\in\overline{\Gamma}([0,1])}d_{Q}(c_{1}\circ\gamma^{*}_{1},O\star(c_{2}\circ
S_{\lambda}\circ\gamma^{*}_{2})).$

We will first show that $F$ is continuous. Let $\tau\in S^{1}$,
$O\in\operatorname{SO}(d)$ and $\epsilon>0$. Since $C(S^{1},{\mathbb{R}}^{d})$
is dense in $L^{2}(S^{1},{\mathbb{R}}^{d})$, let $g\in
C(S^{1},{\mathbb{R}}^{d})$ such that $||Q(c_{2})-g||_{L^{2}}<\epsilon/4$. As
$g$ is continuous on a compact domain it follows, by the Heine-Cantor theorem,
that it is uniformly continuous. Thus, there exists $\delta>0$ such that for
each $\theta\in S^{1}$ and each $\lambda\in S^{1}$ such that
$|\lambda|\leq\delta,$ we have
$|g(\theta)-g(S_{\lambda}(\theta))|<\epsilon/4.$ Thus, for each $\lambda$ such
that $|\lambda|\leq\delta$, we have

$\displaystyle||g-g\circ S_{\lambda}||^{2}_{L^{2}}$
$\displaystyle=\int_{S^{1}}|g(\theta)-g(S_{\lambda}(\theta))|^{2}d\theta$
$\displaystyle<\int_{S^{1}}(\epsilon/4)^{2}d\theta=(\epsilon/4)^{2}.$

Pick this $\delta$ and let $\lambda\in S^{1}$ such that
$|\tau-\lambda|<\delta.$ Thus,

$\displaystyle||g\circ S_{\tau}-g\circ S_{\lambda}||_{L^{2}}=||g-g\circ
S_{\lambda-\tau}||_{L^{2}}<\epsilon/4.$

By a change of variable argument, we can show for any $g_{1},g_{2}\in
C(S^{1},{\mathbb{R}})$ and $\theta\in S^{1}$, we have

$||g_{1}\circ S_{\theta}-g_{2}\circ
S_{\theta}||_{L^{2}}=||g_{1}-g_{2}||_{L^{2}}.$

Furthermore, it is easy to show that for any
$c\in\operatorname{AC}(S^{1},\mathbb{R}^{d})$ and any $\theta\in S^{1}$, we
have $Q(c\circ S_{\theta})=Q(c)\circ S_{\theta}$. On the other hand, the
action of the rotation group on curves induces a corresponding action on their
SRV transform which we write for any $O\in\operatorname{SO}(d)$ as $Q(O\star
c)=O\cdot Q(c)$, where we have specifically that $O\cdot
Q(c)(\cdot)=\frac{1}{\sqrt{|c^{\prime}(\cdot)|}}Oc^{\prime}(\cdot)$. Note that
this action on SRV transforms is by isometry for $\|\cdot\|_{L^{2}}$.

Now, for any $(\lambda,O^{\prime})\in S^{1}\times\operatorname{SO}(d)$ with
$|\tau-\lambda|<\delta$ and
$||O-O^{\prime}||<\epsilon/(4||Q(c_{2})||_{L^{2}})$ (for the operator norm on
matrices), we can write

$\displaystyle|F(\lambda$ $\displaystyle,O)-F(\tau,O^{\prime})|$
$\displaystyle\leq$ $\displaystyle||O\cdot Q(c_{2}\circ
S_{\lambda})-O^{\prime}\cdot Q(c_{2}\circ S_{\tau})||_{L^{2}}$
$\displaystyle=$ $\displaystyle||O\cdot Q(c_{2})\circ
S_{\lambda}-O^{\prime}\cdot Q(c_{2})\circ S_{\tau}||_{L^{2}}$
$\displaystyle\leq$ $\displaystyle||O\cdot Q(c_{2})\circ
S_{\lambda}-O^{\prime}\cdot Q(c_{2})\circ S_{\lambda}||_{L^{2}}$
$\displaystyle+||O^{\prime}\cdot Q(c_{2})\circ S_{\lambda}-O^{\prime}\cdot
Q(c_{2})\circ S_{\tau}||_{L^{2}}.$

For the first term on the right hand side, we can see that

$||O\cdot Q(c_{2})\circ S_{\lambda}-O^{\prime}\cdot Q(c_{2})\circ
S_{\lambda}||_{L^{2}}\leq||O-O^{\prime}||.||Q(c_{2})\circ
S_{\lambda}||_{L^{2}},$

and since $||Q(c_{2})\circ S_{\lambda}||_{L^{2}}=||Q(c_{2})||_{L^{2}}$, we can
bound this term by $\epsilon/4$. On the other hand, we have

$\displaystyle||O^{\prime}\cdot Q(c_{2})\circ S_{\lambda}-O^{\prime}\cdot
Q(c_{2})\circ S_{\tau}||_{L^{2}}$ $\displaystyle=||Q(c_{2})\circ
S_{\lambda}-Q(c_{2})\circ S_{\tau}||_{L^{2}}$
$\displaystyle\leq||Q(c_{2})\circ S_{\lambda}-g\circ
S_{\lambda}||_{L^{2}}+||g\circ S_{\lambda}-g\circ S_{\tau}||_{L^{2}}$
$\displaystyle\phantom{||}+||g\circ S_{\tau}-Q(c_{2})\circ S_{\tau}||_{L^{2}}$
$\displaystyle=||Q(c_{2})-g||_{L^{2}}+||g\circ S_{\lambda}-g\circ
S_{\tau}||_{L^{2}}+||g-Q(c_{2})||_{L^{2}}$
$\displaystyle<\epsilon/4+\epsilon/4+\epsilon/4=3\epsilon/4,$

which finally leads to $|F(\lambda,O)-F(\tau,O^{\prime})|<\epsilon$. Now,
since $F$ is continuous on the compact set $S^{1}\times\operatorname{SO}(d)$,
there exists an optimal $\tau\in S^{1}$ and an optimal
$O\in\operatorname{SO}(d)$ such that
$F(\tau,O)=\inf_{S^{1}\times\operatorname{SO}(d)}F$. Note that the curves
$c_{1}$ and $O\star(c_{2}\circ S_{\tau})$ belong to
$\operatorname{AC}(S^{1},\mathbb{R}^{d})$ and thus in particular to
$\operatorname{AC}([0,1],\mathbb{R}^{d})$, and that by assumption, they are
either both of class $C^{1}$ or one of them is piecewise linear. By the
results of [4, 11], there exist optimal
$\gamma^{*}_{1},\gamma^{*}_{2}\in\overline{\Gamma}([0,1])$ such that

$d_{\mathcal{S}}([c_{1}],[c_{2}])=d_{Q}(c_{1}\circ\gamma^{*}_{1},O\star(c_{2}\circ
S_{\tau}\circ\gamma^{*}_{2})),$

which concludes the proof for the case of closed curves. The proof for open
curves modulo reparametrizations and rotations can be done exactly as above,
by considering a function $\tilde{F}$ that only depends on rotations. ∎

## References

  * [1] Boulbaba Ben Amor, Jingyong Su, and Anuj Srivastava.  Action recognition using rate-invariant analysis of skeletal shape trajectories.  IEEE transactions on pattern analysis and machine intelligence, 38(1):1–13, 2015.
  * [2] Martin Bauer, Martins Bruveris, and Peter W Michor.  Overview of the geometries of shape spaces and diffeomorphism groups.  Journal of Mathematical Imaging and Vision, 50(1-2):60–97, 2014\.
  * [3] Javier Bernal, Gunay Dogan, and Charles R Hagwood.  Fast Dynamic Programming for Elastic Registration of Curves.  In The IEEE Conference on Computer Vision and Pattern Recognition (CVPR) Workshops, June 2016.
  * [4] Martins Bruveris.  Optimal Reparametrizations in the Square Root Velocity Framework.  SIAM Journal on Mathematical Analysis, 48, 07 2015.
  * [5] Mengyu Dai, Zhengwu Zhang, and Anuj Srivastava.  Analyzing dynamical brain functional connectivity as trajectories on space of covariance matrices.  IEEE transactions on medical imaging, 39(3):611–620, 2019.
  * [6] Gunay Dogan, Javier Bernal, and Charles R Hagwood.  A Fast Algorithm for Elastic Shape Distances Between Closed Planar Curves.  In The IEEE Conference on Computer Vision and Pattern Recognition (CVPR), June 2015.
  * [7] E. Hartman, Y. Sukurdeep, E. Klassen, N. Charon, and M. Bauer.  SupervisedDL-SRVFdistances: https://github.com/emmanuel-hartman/supervisedDL-SRVFdistances, 2021.
  * [8] Wen Huang, Kyle A Gallivan, Anuj Srivastava, and Pierre-Antoine Absil.  Riemannian optimization for registration of curves in elastic shape analysis.  Journal of Mathematical Imaging and Vision, 54(3):320–343, 2016\.
  * [9] Diederik P Kingma and Jimmy Ba.  Adam: A Method for Stochastic Optimization, 2014.
  * [10] K. Koneripalli, S. Lohit, R. Anirudh, and P. Turaga.  Rate-Invariant Autoencoding of Time-Series.  In ICASSP 2020 - 2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), pages 3732–3736, 2020.
  * [11] Sayani Lahiri, Daniel Robinson, and Eric Klassen.  Precise Matching of PL Curves in $R^{N}$ in the Square Root Velocity Framework.  Geometry, Imaging and Computing, 2, 01 2015.
  * [12] Suhas Lohit, Qiao Wang, and Pavan Turaga.  Temporal transformer networks: Joint learning of invariant and discriminative time warping.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 12426–12435, 2019.
  * [13] James Stephen Marron, James O Ramsay, Laura M Sangalli, and Anuj Srivastava.  Functional data analysis of amplitude and phase variation.  Statistical Science, pages 468–484, 2015.
  * [14] Washington Mio, Anuj Srivastava, and Shantanu Joshi.  On shape of plane elastic curves.  International Journal of Computer Vision, 73(3):307–324, 2007.
  * [15] Elvis Nunez and Shantanu H Joshi.  Deep learning of warping functions for shape analysis.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR) Workshops, June 2020.
  * [16] Kaushik Koneripalli Seetharam.  Structured Disentangling Networks for Learning Deformation Invariant Latent Spaces.  PhD thesis, Arizona State University, 2019.
  * [17] Anuj Srivastava and Eric Klassen.  Functional and Shape Data Analysis.  Springer Series in Statistics. Springer New York, 2016.
  * [18] Anuj Srivastava, Eric Klassen, Shantanu H Joshi, and Ian H Jermyn.  Shape analysis of elastic curves in euclidean spaces.  IEEE Transactions on Pattern Analysis and Machine Intelligence, 33(7):1415–1428, 2010.
  * [19] Jingyong Su, Sebastian Kurtek, Eric Klassen, Anuj Srivastava, et al.  Statistical analysis of trajectories on riemannian manifolds: bird migration, hurricane tracking and video surveillance.  The Annals of Applied Statistics, 8(1):530–552, 2014.
  * [20] Laurent Younes.  Computable elastic distances between shapes.  SIAM Journal on Applied Mathematics, 58(2):565–586, 1998.