# Let There Be Order: Rethinking Ordering in Autoregressive Graph Generation

Jie Bu
Department of Computer Science
Virginia Tech
Blacksburg, VA 24060
<EMAIL_ADDRESS>
&Kazi Sajeed Mehrab
Department of Computer Science
Virginia Tech
Blacksburg, VA 24060
<EMAIL_ADDRESS>
&Anuj Karpatne
Department of Computer Science
Virginia Tech
Blacksburg, VA 24060
<EMAIL_ADDRESS>

###### Abstract

Conditional graph generation tasks involve training a model to generate a
graph given a set of input conditions. Many previous studies employ
autoregressive models to incrementally generate graph components such as nodes
and edges. However, as graphs typically lack a natural ordering among their
components, converting a graph into a sequence of tokens is not
straightforward. While prior works mostly rely on conventional heuristics or
graph traversal methods like breadth-first search (BFS) or depth-first search
(DFS) to convert graphs to sequences, the impact of ordering on graph
generation has largely been unexplored. This paper contributes to this problem
by: (1) highlighting the crucial role of ordering in autoregressive graph
generation models, (2) proposing a novel theoretical framework that perceives
ordering as a dimensionality reduction problem, thereby facilitating a deeper
understanding of the relationship between orderings and generated graph
accuracy, and (3) introducing "latent sort," a learning-based ordering scheme
to perform dimensionality reduction of graph tokens. Our experimental results
showcase the effectiveness of latent sort across a wide range of graph
generation tasks, encouraging future works to further explore and develop
learning-based ordering schemes for autoregressive graph generation.

##  1 Introduction

We consider the problem of generating graphs given conditional inputs, where
there exists precisely one target graph we are interested in generating for
every conditional input. We call this problem as paired graph generation,
since every input condition is paired with exactly one target graph. This
setting is highly relevant in many real-world applications such as predicting
the road network of a region given its satellite image [1], generating scene
graphs with semantic relationships among objects given an image of the objects
[2, 3, 4], and inferring the underlying circuit graph of an electrical
instrument given its sensor observations [5]. This is different from the
unpaired setting of graph generation problems studied by several previous
works, where we are interested in generating distributions of graphs rather
than a specific target graph.

The problem of graph generation has some similarity to natural language
generation [6, 7, 8], since they both involve generating outputs with variable
sizes. However, a key distinction is that there is no natural ordering of
graph components (e.g., nodes and edges) as opposed to the sequential nature
of words in a sentence. Graph generation is thus a particularly challenging
problem due to the inherent diversity of graphs involving varying sizes and
different topologies. Many previous approaches for graph generation [9, 10,
11, 1, 12] have adopted autoregressive models such as recurrent neural
networks [13, 14] and auto-regressive transformers [15, 16, 17] to handle the
variable sizes and topologies in graphs. These models often convert graph
components into vector representations or tokens, and incrementally generate
tokens representing graph components in an ordered sequence.

However, existing autoregressive models for graph generation face unique
challenges in sorting multi-dimensional graph tokens into ordered sequences.
As graph components often have no natural ordering [18], there can be $n!$
possible orderings for $n$ graph tokens. It has been shown by [19, 20, 11]
that finding an optimal ordering of unsorted elements in a set is key to the
success of autoregressive models, since certain orderings may be favored by
the models leading to better performance. Despite this, existing works on
autoregressive graph generation [10, 11, 1, 12] make arbitrary choices for
ordering graph components based on graph traversals such as breadth-first-
search (BFS) or depth-first-search (DFS) without fully justifying their
choices empirically or theoretically.

In this work, we present a novel perspective to address this issue by
reframing the problem of ordering a set of graph tokens as that of performing
dimensionality reduction (DR). Specifically, by learning a generalizable
mapping of every possible graph token to a 1-D representation, we
theoretically and empirically show that a set of graph tokens can be
conveniently ordered by sorting their 1-D values. Our contributions can be
summarized as follows:

  * •

We develop novel theoretical frameworks to study the accuracy of DR techniques
for sorting unordered tokens in autoregressive models, opening a new line of
research in developing learning-based ordering schemes for paired graph
generation.

  * •

We propose the latent sort algorithm, a novel learning-based ordering scheme
employing DR for autoregressive graph generation. We provide theoretical
bounds of the errors of latent sort and highlight an interesting connection
between performing latent sort and finding approximate solutions to the
shortest path problem.

  * •

We propose a strong backbone network for autoregressive graph generation named
Graph Auto-Regressive Transformer (GraphART) based on a modified GPT-2
architecture [17]. GraphART is purely autoregressive in nature and does not
use any graph-specific architectural components such as graph neural networks
(GNNs), in contrast to prior works.

  * •

We empirically show that different token ordering schemes result in nontrivial
performance gaps in autoregressive models such as GraphART. We also show that
Latent Sort is versatile and achieves competitive performance over a wide
range of graph generation tasks.

##  2 Related Works

A number of approaches have been developed for graph generation in a
probabilistic setting where the goal is to generate a distribution of graphs
that resemble a target distribution given conditional inputs [9, 10, 21, 22,
23, 11, 24, 25, 26]. These works fall in the “unpaired” graph generation
category, since there is no unique ground-truth target graph that is paired
with every input condition. Hence, the evaluation of performance of these
works typically focuses on measuring the validity of the generated graphs in
terms of statistical similarity with the target distribution, rather than
measuring deviation from its corresponding ground-truth graph in a “paired”
fashion, which is the focus of our work.

One line of work for paired graph generation involves using non-autoregressive
approaches to generate graphs [27, 28, 29, 3, 4]. These approaches generate
entire graphs rather than individual components, which are then matched to the
target graphs using bipartite matching algorithms (e.g., Hungarian matching in
DETR [30] for object detection). However, bipartite matching loss functions
can slow down training due to their high computational costs [31] and can
potentially suffer from slow convergence [32, 33]. Moreover, these models
cannot direclty handle output graphs with varying sizes and have to resort to
padding, making them inefficient for real-world applications.

Another line of work involves autoregressive models to generate variable-
length sequences of graph tokens. As a necessary preprocessing step, these
methods require an approach to sequentialize target graphs into an ordered
sequence of tokens, where finding the ideal scheme for token ordering can be
non-trivial. Vinyals et al. [19] were the first to reveal the importance of
ordering when using sequence-to-sequence models to process sets. They proposed
a training algorithm that searches for the optimal order, but scalability
issues would arise with large sets, and inexact searches or sampling methods
could negatively impact performance. More recently, Chen et al. [20] discuss
the importance of ordering for unpaired autoregressive graph generation, where
they derive the joint probability over the graph for sorting the nodes. To
sequentialize graphs, most of the existing works favor traditional graph
traversal methods like BFS [10, 11] or DFS [34, 12, 11] . However, to the best
of our knowledge, no prior work has theoretically and empirically compared the
performance of autoregressive models for graph generation with varying token
ordering schemes, which is one of the contributions of our work.

##  3 Rethinking Sorting As a Dimensionality Reduction Problem

Given a set of $M$ unordered points (or vector tokens)
$\mathcal{X}\subseteq\mathbb{R}^{N}$, $|\mathcal{X}|=M$ in an $N$-dimensional
space ($N>1$), we are interested in finding the “optimal” ordering of points
in $\mathcal{X}$, denoted by $Y^{*}\in\mathbb{R}^{M\times N}$, which when used
as the target sequence to supervise an autoregressive model yields optimal
performance in generating $\mathcal{X}$. Formally, let us denote the ordered
sequence $Y^{*}$ as the matrix
$\left[\bm{x_{1}}^{*},\bm{x_{2}}^{*},...,\bm{x_{M}}^{*}\right]^{\intercal}$,
where every row $\bm{x_{i}}^{*}$ of $Y^{*}$ is a point in $\mathcal{X}$ and
$i$ denotes its sorted index in the optimal ordering. For now, we assume that
such an ordering exists for every set $\mathcal{X}$. Later, in Section 3.3, we
will describe some of the ideal properties of $Y^{*}$.

We are interested in finding the optimal ordering $Y^{*}$ by sorting high-
dimensional vector tokens in $\mathcal{X}$. This is challenging due to the
absence of a well-defined comparison operation in the $N$-dimensional space of
vector tokens, resulting in a factorial number of possible orderings to be
evaluated that is computationally prohibitive. To tackle this issue, we
consider sorting points in $\mathcal{X}$ as a dimensionality reduction (DR)
problem, as described in the following.

###  3.1 DR-based Sorting Algorithms

Figure 1: A general pipeline for sorting algorithms using dimensionality
reduction (DR) mapping $f$.

Figure 1 provide an overview of the generic DR-based sorting pipeline for any
DR mapping $f$. The key insight that we leverage here is that while sorting
points is ill-defined in $N$-dimensional space, it is well-defined on a 1-D
space (see Remark E.1 in Appendix). Hence, if we can learn a generalizable DR
mapping from every possible token $\bm{x}$ to a corresponding 1D
representation $h$, we can order tokens in any arbitrary set $\mathcal{X}$ by
sorting their 1D representations. We can thus define DR-based sorting
algorithms as follows.

######  Definition 3.1 (DR-based Sorting Algorithm).

Given a sorting algorithm $s:\mathcal{X}\rightarrow Y$, let
$Y=\left[\bm{x_{1}},\bm{x_{2}},...,\bm{x_{M}}\right]^{\intercal}$ be the
resulting ordered sequence for the unordered set
$\mathcal{X}\subseteq\mathbb{R}^{N}$. We can then represent the sorting
algorithm $s$ by a DR mapping $f:\mathbb{R}^{N}\rightarrow\mathbb{R}$, such
that the ordering in the original $N$-D space is given by sorting their 1-D
values, $h_{i}=f(\bm{x_{i}})\in\mathbb{R}$, and $\forall i\leq j\leq
M,f(\bm{x_{i}})\leq f(\bm{x_{j}})$.

Note that the DR mapping $f$ maps every token $\bm{x_{i}}$ to a universal
position $h_{i}$ in the latent space, regardless of the combination of tokens
in $\mathcal{X}$ that it appears with. Ordering $\mathcal{X}$ thus reduces to
finding a sorted path traversal in the 1D latent space. Thus, we formulate the
problem of finding the optimal sorting $s^{*}:\mathcal{X}\rightarrow Y^{*}$ as
finding an optimal DR mapping $f^{*}:\mathbb{R}^{N}\rightarrow\mathbb{R}$,
$h_{i}^{*}=f^{*}(\bm{x_{i}}^{*})$. While sorting itself is challenging to
perform in a differentiable way, finding differentiable approaches for DR is
relatively easier. This opens up new possibilities for developing learning-
based methods to order points in high-dimensional spaces. For simplicity, we
assume the 1-D representations of all the points in $\mathcal{X}$, denoted by
$\mathcal{H}$, are normalized to $[0,1]$.

Latent Sort Algorithm: We propose latent sort as our DR-based sorting
algorithm, where the DR mapping $f$ is represented by an auto-encoder model
consisting of an MLP encoder $f_{e}$ and an MLP decoder $f_{d}$. We first
train $f_{e}$ and $f_{d}$ to reconstruct all tokens in the dataset. We then
freeze $f_{e}$ and plug it in the pipeline shown in Figure 1 as the DR mapping
$f$ for sorting. The property of latent sort will be analyzed in following
sections.

###  3.2 Analyzing The Errors of DR-based Sorting Algorithms

To quantify the errors in a sorting algorithm $s$, we introduce a probability
matrix $P\in\mathbb{R}^{M\times M}$ where each entry $p_{ij}$ of $P$
represents the probability of $\bm{x_{i}}$ being $\bm{x_{j}}^{*}$. Thus, we
have $\mathbb{E}[Y]=PY^{*}$, and $P$ can be viewed as a soft permutation
matrix that permutes $Y^{*}$ into $Y$. We use the Frobenius norm of the
difference between $\mathbb{E}[Y]$ and $Y^{*}$ as the expected error in a
sorting algorithm,

$\vspace{-0.05cm}\mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}=\big{\|}\mathbb{E}[Y]-Y^{*}\big{\|}_{F}^{2}=\big{\|}PY^{*}-Y^{*}\big{\|}_{F}^{2},$
(1)

where $\|\cdot\|_{F}$ denotes the Frobenius norm. We can see that two factors
contribute to the error term $\mathcal{E}$: $P$ and $Y^{*}$. It is obvious how
$P$ affects the error, as $P$ measures the deviation from the optimal ordering
$Y^{*}$, which reflects the property of an ordering method. We show that
minimizing $\mathcal{E}\big{(}\mathbb{E}[Y],Y^{*}\big{)}$ is equivalent to
minimizing $\|P-I_{M}\|_{F}^{2}$ (see Lemma D.1 in the Appendix). On the other
hand, it is less straight-forward to see how the choice of $Y^{*}$ affects the
error, which we will discuss in Section 3.3. Here, we give two sources of
errors that can lead to suboptimal $P$.

Figure 2: Illustration of the ordering ambiguity on 2D points within $[0,1]$.
Colored contour lines represent points with same 1D latent value, i.e., they
belong to the same ordering ambiguity set.

Errors from Ordering Ambiguity: The first source of errors is the ordering
ambiguity. An ideal mapping $f^{*}$ is one that does not result in any
ambiguity in sorting, i.e., $f^{*}(\bm{x_{i}})=f^{*}(\bm{x_{j}})$ iff $i=j$.
However, such an ideal mapping only exists for certain conditions (see Remark
E.3 in Appendix), otherwise $f$ is surjective (non-bijective) even if it has
no reconstruction errors. Thus, it is possible that a DR mapping $f$ assigns
the same 1-D representation to different points in $\mathcal{X}$, which can
cause ambiguity in sorting these points. To describe this ambiguity, we
introduce the concept of ordering ambiguity set.

######  Definition 3.2 (Ordering Ambiguity Set).

Let $\mathcal{A}_{i}$ be the ordering ambiguity set for a point
$\bm{x_{i}}\in\mathcal{X}$, where
$\mathcal{A}_{i}=\left\\{\bm{x_{j}}|\forall\bm{x_{j}}\in\mathcal{X},f(\bm{x_{i}})=f(\bm{x_{j}})\right\\}$.

Ordering ambiguity occurs when there exists more than one element in the
ordering ambiguity set, i.e., $|\mathcal{A}_{i}|>1$. Different DR mappings ,
can have varying ordering ambiguity patterns. Figure 2 shows ordering
ambiguity patterns in 2D space for two DR mappings: (i) Mean Squared Sort:
$h_{i}=(x_{i1}^{2}+x_{i2}^{2})/2$, (ii) Summation Sort: $h_{i}=x_{i1}+x_{i2}$.

To quantify the error introduced by ordering ambiguity, we consider the
expectation of the ambiguity set $\mathcal{A}_{i}$ for each point
$\bm{x_{i}}\in\mathcal{X}$. It is reasonable to assume that each row $i$ of
$P$ is a uniform distribution over the corresponding ambiguity set
$\mathcal{A}_{i}$ (i.e., $\forall\bm{x_{j}}^{*}\in\mathcal{A}_{i}$,
$p_{ij}={1}/{|\mathcal{A}_{i}|}$ and
$\forall\bm{x_{j}}^{*}\not\in\mathcal{A}_{i}$, $p_{ij}=0$. Putting this in
Equation (1), we get (see Section C.3 in Appendix for derivation):

$\displaystyle\vspace{-0.08cm}\mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}$
$\displaystyle=\sum_{i=1}^{M}\Bigg{\|}\left(\frac{1}{|\mathcal{A}_{i}|}\sum_{\bm{x_{j}}\in\mathcal{A}_{i}}\bm{x_{j}}\right)-\bm{x_{i}}\
\Bigg{\|}^{2}$ (2)

Errors From Imperfect Reconstruction In Latent Sort: Since latent sort trains
an autoencoder using reconstruction loss for DR mapping, it is beneficial to
discuss the potential errors from imperfect reconstructions (See Section C.4
and C.5 in the Appendix for details). Imperfect reconstruction leads to a
discrepancy between the reconstructed token $\hat{x}$ and the original token
$x$. Assuming the reconstruction follows normal distribution, the
reconstruction error will also result in a normal distributed error with
bounded mean and variance in the latent space (See Theorem C.3 in Appendix).
This error can affect the ordering results when two points swap their
positions in the sequence due to their 1-D latent values crossing each other.
The probability of swap can be quantified by examining the overlap between
their 1-D latent distributions. Since the 1-D latent values follow Gaussian
distribution, and the overlaps of normal distributions are very small if the
mean values of the distributions are far from each other, we only consider the
swapping between neighboring points since the non-neighboring points will have
larger mean value differences and thus are less likely to be swapped (See
Section C.6 in Appendix for details).

###  3.3 The Shortest Path Property of Desirable Ordering Schemes

In previous theories, we have discussed the error of sorting algorithms using
the probability matrix $P$. In this section, we discuss the desired properties
of $Y^{*}$ that minimize the errors induced by given $P$ matrices.

From Section C.7 in the Appendix, we have the simplified equation for error
$\mathcal{E}$ as follows:

$\displaystyle\vspace{-0.05cm}\mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}=\sum_{i=2}^{M-1}\left\|p_{i(i-1)}(\bm{x_{i-1}}^{*}-\bm{x_{i}}^{*})+p_{i(i+1)}(\bm{x_{i+1}}^{*}-\bm{x_{i}}^{*})\right\|_{F}^{2}\vspace{-0.05cm}$
(3)

Recall that $f$ maps each point in $\mathcal{X}$ to its 1D latent
representation independently of other points. Therefore, it is unaware of the
global information of the entire sequence. That being said, we expect (3) to
be minimized for any three points
$\\{\bm{x_{i-1}}^{*},\bm{x_{i}}^{*},\bm{x_{i+1}}^{*}\\}\subset\mathcal{X}$. To
make $\mathcal{E}\big{(}\mathbb{E}[Y],Y^{*}\big{)}$ minimized for all possible
$\mathcal{X}\subseteq\mathbb{R}^{N}$, we need to minimize the upper bound of
(3), giving by:

$\displaystyle\ \ \ \
\sum_{i=2}^{M-1}\left\|p_{i(i-1)}(\bm{x_{i-1}^{*}}-\bm{x_{i}}^{*})+p_{i(i+1)}(\bm{x_{i+1}}^{*}-\bm{x_{i}}^{*})\right\|_{F}^{2}$
(4)
$\displaystyle\leq\sum_{i=2}^{M-1}\left\|p_{i(i-1)}(\bm{x_{i-1}}^{*}-\bm{x_{i}}^{*})\right\|_{F}^{2}+\sum_{i=2}^{M-1}\left\|p_{i(i+1)}(\bm{x_{i+1}}^{*}-\bm{x_{i}}^{*})\right\|_{F}^{2},\vspace{-0.06cm}$
(5)

assuming both $p_{i(i-1)}$ and $p_{i(i+1)}$ are nonzero for each $i$. In order
to minimize errors from imperfect $P$, the ideal sorting $Y^{*}$ should
minimize $(\bm{x_{i-1}}^{*}-\bm{x_{i}}^{*})$ and
$(\bm{x_{i+1}}^{*}-\bm{x_{i}}^{*})$, i.e., minimize the distances between
neighboring pairs. This property connects to the shortest path problem or
traveling salesman problem (TSP), with the target sorting being the TSP
solution on the input point set $\mathcal{X}$. However, TSP is NP-hard and
solving it is computationally expensive. Using TSP solutions also induces
difficulties for the autoregressive models to capture the ordering rule
represented by the TSP solving algorithms. Thus, an approximation to the TSP
solution should provide a balanced option between simplicity for learning and
error tolerance.

Theorem C.1 and C.2 in the Appendix state the correlation between the distance
between two points and their 1-D latent distance. Since neighboring points in
the sorted sequence have the closest latent values, the distance between
neighboring points in the sorted sequence are expected to be small when the
decoder is properly regularized. However, our experiments suggest that
reconstruction loss and regularization techniques alone may not sufficiently
approximate shortest path solutions. To further enhance the shortest path
property in the latent sort algorithm, we introduce the latent gradient
penalty (LGP) loss given by:

$\displaystyle\text{LGP}(\mathcal{X},f_{e})=\min_{f_{e}}\sum_{i=2}^{M-1}\sum_{j=2}^{M-1}{g_{ij}^{2}},$
$\displaystyle\text{where}\
g_{ij}=\begin{cases}\frac{\|\bm{x_{i}}-\bm{x_{j}}\|}{|h_{i}-h_{j}|+\beta}-\alpha&\text{if}\
|i-j|=1,\\\ 0&\text{otherwise}\end{cases}$ (6)

where $f_{e}$ is the latent sort encoder, and $\alpha$ and $\beta$ are
positive constants. The LGP aims to make the 1-D latent value distance
$|h_{i}-h_{j}|$ between two neighboring points proportional to their distance
$\|\bm{x_{i}}-\bm{x_{j}}\|$ in the original space, with proportionality
constant constrained by $\alpha$ (we chose $\alpha=1$). Without extra
justification, the latent sort algorithm in the paper are by default trained
with LGP.

##  4 Experiment Setups

We briefly summarize the experiment setups in this section (see Appendix B for
complete details of the processing pipelines for full reproducibility). All
codes are publicly on Github 111https://github.com/jayroxis/ordering-in-graph-
generation.

###  4.1 Application Tasks and Datasets

(a) Paired Graph Generation Applications and Pipelines.

(b) GraphART Decoder Block

Figure 3:  Left: schematic views of the three paired graph generation tasks
(from top to bottom): scene graph generation [3], topological graph extraction
[1], and circuit graph prediction [5]. Right: The decoder block of the
proposed GraphART model using multi-head causal attention architecture.

We consider the following three application tasks for paired graph generation
(see Figure 3(a)).

Scene Graph Generation: The goal of this task is to take an image as input and
generate a scene graph as output. We use the OpenPSG dataset [3], which
consists of semantic-level triplets including two objects (nodes) and
predicates (edges) between them.

Topological Graph Extraction: The goal here is to predict the spatial location
of nodes in a graph and the undirected edges between them, given an image
representation of the graph. We use the Toulouse Road Network (TRN) dataset
[1] and our synthetic Planar Graph dataset, which can randomly generate image-
graph pairs that are more challenging than the ones in the TRN dataset. We
average the results over 10 random runs for every method on the Planar Graph
dataset.

Circuit Graph Prediction: For this task, we consider the Terahertz channelizer
dataset [5], which comprises 347k circuit graphs. The goal is to predict the
ground-truth circuit graph comprises of 3 to 6 resonators, given desirable
electromagnetic (EM) properties as inputs. Each resonator is represented by a
9-dimensional vector. Since the scale of some dimensions in the resonator
vectors is much larger than others, which may dominate the prediction errors,
we normalize the data using min-max normalization per dimension. We report the
errors on the normalized data.

###  4.2 Models and Pipelines

Additional details for the baseline models can be found in Section B in the
Supplementary.

Backbone Models: We propose two backbone models for paired graph generation:
GraphART and GraphLSTM. These models rely solely on their autoregressive
nature to sequentially predict graph tokens, without involving graph-specific
architecture components such as graph neural networks (GNNs). Both models
consist of an encoder to encode the input conditions (e.g., images) and an
autoregressive decoder for graph generation. The decoder of GraphART is based
on our modified version of the GPT-2 architecture, as shown in Figure 3(b). It
incorporates pre-normalization [35] and LayerScale [36] designs. In GraphLSTM,
we simply replace the autoregressive transformer of GraphART with a long
short-term memory (LSTM) model [13]. These models can be coupled with any
sorting algorithm for ordering graph tokens during training.

Baseline Sorting Algorithms: (1) Mean Squared Sort: $n$-dimensional tokens are
sorted based on the mean squared value of their $n$ dimensions. This heuristic
sorts the graph tokens based on their overall magnitudes, from larger to
smaller $L^{2}$ norms. (2) Lexicographical Sort: We start by sorting tokens
based on their first dimension, then the second dimension if there are ties,
and so on until all dimensions have been considered. Note that this heuristic
assumes a meaningful ordering of dimensions in the data, which may not be
true. (3) SVD Low-Rank Sort: We perform a Singular Value Decomposition (SVD)
of the tokens in a graph and sort the tokens by projecting them onto the
direction of maximum variance, identified using the largest principal
component. This serves as a linear baseline for dimensionality reduction, as
opposed to the non-linear Latent Sort. (4) BFS or DFS Sort: We use the
breadth-first search (BFS) or depth-first search (DFS) traversal on edges,
implemented by NetworkX [37].

Graph Tokenization Pipeline: To address the potential disadvantages of
transformers attending different token distributions for each category of
graph components (e.g., nodes and edges), as advocated in [36], we implemented
a novel edge-based tokenization approach. This approach is similar to, but
distinct from, the method used in Graph Transformers [38, 39]. Specifically,
whenever possible, we concatenated the features of each edge with graph-level
features and the features of the two connecting nodes to create edge-based
tokens. For example, in the topological graph extraction task where only node
coordinates are available as features, we constructed edge tokens by combining
the coordinates of the two connecting nodes for each edge. In the scene graph
generation task, we followed the convention of using one-hot encoded triplets
as edge tokens. In the circuit graph prediction task, we utilized a
9-dimensional vector as a token to represent each node (resonator), as
employed in [5].

Loss Functions: We trained all autoregressive models using a teacher forcing
setup. For directed graphs, we calculated the loss between each element in the
ground-truth sequence $Y_{gt}$ and its corresponding element in the predicted
sequence $Y_{pred}$ using an appropriate choice of loss function
$\mathcal{L}(Y_{pred},Y_{gt})$. For circuit extraction, we employed the
elastic loss function, which is the sum of L1 and MSE loss, as $\mathcal{L}$.
In scene graph generation, we used the binary cross-entropy (BCE) loss instead
of cross-entropy (CE) loss, following previous works [40, 41]. For the
topological graph extraction task involving undirected graphs, we utilized the
elastic loss as $\mathcal{L}$, but we considered its undirected variant
$\mathcal{L}_{u}$. This variant accounts for swapping the node features in
every edge token. The undirected graph loss is defined as
follows:$\mathcal{L}_{u}(Y_{pred},Y_{gt})=\mathcal{L}(Y_{pred},Y_{gt})+\mathcal{L}(Y_{pred},\overline{Y}_{gt})$,
where $\overline{Y}_{gt}$ is the sequence where the two node features in every
edge token are swapped.

Evaluation Metrics: During testing, we employed the Edge Mover distance (EMD),
Edge Hausdorff distance (EHD), and StreetMover Distance (SMD) from [1] as
evaluation metrics to assess the similarity between predicted graphs and
ground-truth graphs. The EMD calculates the average of the minimum distances
between each point in the sequence $Y_{pred}$ and its nearest neighbor in the
set $Y_{gt}$. The EHD measures the maximum distance between the closest points
of the two sequences. On the other hand, the SMD is based on point cloud
distances. Further details on these metrics can be found in Supplementary B.
For scene graph generation, we consider precision, recall, and F-1 scores as
the predicted sequences consist of multi-class labels. Additionally, we report
the differences in sizes between the predicted graphs and the target graphs as
another metric in each application.

##  5 Experimental Results

###  5.1 Comparing Performance on Paired Graph Generation Tasks

Table 1: Toulouse Road Network Generation. Model | Sorting / Matching |  SMD [1]
---|---|---
GraphART | Latent Sort | 0.0075
GraphART | Mean Squared Sort | 0.0081
GraphART | Lexicographical Sort | 0.0422
GraphART | SVD Low-Rank Sort | 0.0751
GraphART | BFS Sort | 0.0893
GraphART | DFS Sort | 0.0761
GGT [1] | BFS Sort | 0.0158
GraphRNN [1, 10] | BFS Sort | 0.0245

Tables 1 and 2 compare the performance of the Latent Sort algorithm with
baseline sorting algorithms and models on the road network generation and
planar graph generation tasks, respectively. The top-2 models are highlighted
in bold in all tables in this paper. We observe that GraphART coupled with
Latent Sort consistently demonstrates superior performance for topological
graph extraction compared to baseline methods such as GGT [1], GraphRNN [1,
10], and GraphTR. This is noteworthy as GGT and GraphRNN utilize graph-
specific architectures to predict adjacency matrices in addition to graph
generation, unlike our approach. These results showcase the effectiveness of
using purely autoregressive models for graph generation when combined with an
appropriate graph token ordering scheme like Latent Sort.

Among the various sorting algorithms used with GraphART, Mean Squared Sort
stands out for its impressive performance on both datasets, despite being a
simple heuristic. This could potentially be attributed to its simplicity in
sorting based on token magnitudes, which the autoregressive models can
effectively capture. It highlights the importance of ordering schemes being
"easy to learn" which is commonly neglected. On the other hand, BFS and DFS
sort consistently exhibit inferior performance on both datasets for paired
graph generation, despite being the de facto convention used in existing
methods for probabilistic graph generation. We hypothesize that although BFS
and DFS can generate viable distributions of graphs (e.g., molecules) in a
probabilistic setting, they perform poorly when an exact match with a target
graph is required. This is because they do not consider the values of graph
tokens, but only their relative positions in the graph, which can be
arbitrarily defined. Moreover, the

Table 2: Planar Graph Generation. Model | Sorting / Matching | EMD | EHD | SMD | Size. Diff.
---|---|---|---|---|---
GraphART | Latent Sort | 0.0038 | 0.0166 | 0.0002 | 0.0021
GraphART | Mean Squared Sort | 0.0031 | 0.0151 | 0.0002 | -0.0360
GraphART | Lexicographical Sort | 0.0141 | 0.0564 | 0.0013 | -0.2434
GraphART | SVD Low-Rank Sort | 0.0381 | 0.0951 | 0.0068 | 4.8280
GraphART | BFS Sort | 0.0217 | 0.0579 | 0.0017 | 1.7296
GraphART | DFS Sort | 0.0213 | 0.0572 | 0.0017 | 1.6399
GraphLSTM | Latent Sort | 0.0100 | 0.0293 | 0.0006 | 0.1639
GraphLSTM | Mean Squared Sort | 0.0063 | 0.0211 | 0.0004 | 0.0060
GraphLSTM | Lexicographical Sort | 0.0248 | 0.0609 | 0.0030 | 3.9318
GraphLSTM | SVD Low-Rank Sort | 0.0480 | 0.0940 | 0.0069 | 4.5328
GraphLSTM | BFS Sort | 0.0446 | 0.0888 | 0.0059 | 3.6330
GraphLSTM | DFS Sort | 0.0391 | 0.0738 | 0.0045 | 1.7268
GraphTR | Hungarian Matcher | 0.0069 | 0.0348 | 0.0003 | -0.0130

autoregressive model needs to understand the global structure of the graph to
determine the sequence ordering, adding significant complexity in learning the
ordering. Additionally, we observe that GraphLSTM generally does not perform
as well as GraphART. This could be due to the advantage of transformers, which
are capable of attending to longer context that help it gain better
understanding on global information in graphs compared to LSTMs.

Table 3: Semantic Scene Graph Generation on OpenPSG [3] Dataset . Model | Sorting / Matching | Additional Labels | Precision | Recall | F1-Score
---|---|---|---|---|---
GraphART | Latent Sort | None | 0.2537 | 0.1724 | 0.2053
GraphART | Lexicographical Sort | None | 0.2425 | 0.1867 | 0.2110
GraphART | SVD Low-Rank Sort | None | 0.2496 | 0.1696 | 0.2020
GraphLSTM | Latent Sort | None | 0.0115 | 0.0094 | 0.0104
GraphLSTM | Lexicographical Sort | None | 0.0110 | 0.0093 | 0.0101
GraphLSTM | SVD Low-Rank Sort | None | 0.0099 | 0.0079 | 0.0088
PSGTR [3] | Hungarian Matcher | Panoptic Segmentation | 0.1810 | 0.2116 | 0.1920
PSGFormer [3] | Query Matching Block | Panoptic Segmentation | 0.0438 | 0.0520 | 0.0467
Table 4: Results on Circuit Graph Prediction [5].  Model | Token Sorting / Matching | EMD | EHD | Size. Diff.
---|---|---|---|---
GraphART | Latent Sort | 0.0509 | 0.0856 | 0.6320
GraphART | Mean Squared Sort | 0.0558 | 0.0860 | 0.0794
GraphART | Lexicographical Sort | 0.0623 | 0.1080 | 1.4757
GraphART | SVD Low-Rank Sort | 0.0590 | 0.0863 | 1.7485
GraphLSTM | Latent Sort | 0.0496 | 0.0757 | 0.4249
GraphLSTM | Mean Squared Sort | 0.0576 | 0.0836 | -0.1305
GraphLSTM | Lexicographical Sort | 0.0528 | 0.0803 | -0.1205
GraphLSTM | SVD Low-Rank Sort | 0.0516 | 0.0777 | 0.6112
Circuit-GNN [5] | N/A | 0.0738 | 0.0943 | (Set to GT)

Tables 3 and 4 present the results for the scene graph generation and circuit
graph prediction tasks, respectively. Both tasks involve directed graphs. It
is important to note that the mean squared sort cannot be applied to scene
graph generation due to the one-hot encoding of graph tokens, which leads to
identical mean squared values for all tokens. Additionally, the BFS and DFS
sorting algorithms could hardly be applied to both datasets as they mostly
consist of disconnected small graphs. In the scene graph generation task, we
observe that GraphART exhibits higher precision compared to the state-of-the-
art (SOTA) baselines PSGTR and PSGFormer for all sorting algorithm choices,
where PSGTR [3] and PSGFormer [3] use additional panoptic segmentation labels
for training while we only rely on the scene graph labels. However, it has
lower recall than PSGTR. It is worth mentioning that PSGTR has an advantage
over our method as it utilizes additional supervision from panoptic
segmentation during training and knows the ground-truth sizes of target graphs
as user-specified hyperparameters during testing. Despite this, GraphLSTM
shows better F-1 scores than PSGTR for multiple sorting algorithms. Moreover,
we observe that lexicographical sort performs the best on this dataset,
followed closely by Latent Sort. Using lexicographical sort as a domain-
specific heuristic for this problem aligns well with the categorical nature of
graph tokens. We include additional results in Section A in the Supplementary,
including the discussion of label noises in the dataset that lead to size
differences between GraphART generated scene graph vs the ground-truth labels.
For the circuit graph prediction task, both GraphART and GraphLSTM models
outperform Circuit-GNN. Among the ordering schemes, Latent Sort exhibits a
slight advantage over others.

In summary, while specific heuristics may excel in certain applications (e.g.,
mean squared sort in topological graph generation and lexicographical sort in
scene graph generation), Latent Sort demonstrates versatility and consistently
delivers competitive performance across all four datasets.

###  5.2 Analyzing Sorting Ambiguity and Shortest Path Property of Latent Sort

Figure 4: 1-D representation of 2-D points learned by latent sort with (right)
and without (left) using latent gradient penalty (LGP). White lines are
contours showing regions with same latent values.

To analyze the importance of using latent gradient penalty (LGP) in the
training of 1-D representations by latent sort, Figure 4 shows the contour
lines in the 1-D space learned by latent sort for a toy 2-D dataset, with and
without LGP. The latent space learned without LGP exhibits large regions with
identical values, which, according to the theory of ordering ambiguity, can
lead to significant errors. In contrast, the latent space learned with LGP
demonstrates finer partitioning, resulting in smaller regions sharing the same
values. This leads to reduced ambiguity in ordering and, consequently, lower
errors due to ordering ambiguity.

Figure 5: Approximating shortest path solutions using latent sort (with LGP).

Moreover, to empirically support the theoretical claim that latent sort finds
an approximate solution to the shortest path problem, Figure 5 showcases the
paths identified by latent sort over different collections of points in 2-D
space. These paths are determined by sorting the points based on their 1-D
values. The examples demonstrate that latent sort reasonably approximates the
shortest path solution, even for a large number of points. For additional
quantitative results, please refer to Appendix A. These findings highlight an
intriguing connection between Latent Sort and shortest path problems,
motivating a new use-case for Latent Sort as a fast, low-cost method for
approximating shortest path solutions through GPU-based parallelization of
neural networks. This is particularly beneficial for applications with strict
runtime efficiency requirements.

##  6 Conclusions

We presented a novel autoregressive framework for paired graph generation by
reframing the task of sorting unordered tokens as a dimensionality reduction
problem from a theoretical standpoint. We showed the efficacy of our proposed
latent sort algorithm on various graph generation tasks. However, we also
acknowledge certain limitations of our work. While latent sort performs
competitively across a wide range of datasets, it does not consistently
outperform application-specific heuristics. This suggests the need for future
works to explore more sophisticated learning-based ordering schemes.
Additionally, our experiments are currently limited to small graphs,
emphasizing the importance of investigating the scalability of our approach to
larger graphs in future research.

## References

  * [1] Davide Belli and Thomas Kipf.  Image-conditioned graph generation for road network extraction.  arXiv preprint arXiv:1910.14388, 2019.
  * [2] Yichao Lu, Himanshu Rai, Jason Chang, Boris Knyazev, Guangwei Yu, Shashank Shekhar, Graham W Taylor, and Maksims Volkovs.  Context-aware scene graph generation with seq2seq transformers.  In Proceedings of the IEEE/CVF international conference on computer vision, pages 15931–15941, 2021.
  * [3] Jingkang Yang, Yi Zhe Ang, Zujin Guo, Kaiyang Zhou, Wayne Zhang, and Ziwei Liu.  Panoptic scene graph generation.  In ECCV, 2022.
  * [4] Rongjie Li, Songyang Zhang, and Xuming He.  Sgtr: End-to-end scene graph generation with transformer.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 19486–19496, 2022.
  * [5] Guo Zhang, Hao He, and Dina Katabi.  Circuit-gnn: Graph neural networks for distributed circuit design.  In International Conference on Machine Learning, pages 7364–7373, 2019.
  * [6] Tom Brown, Benjamin Mann, Nick Ryder, Melanie Subbiah, Jared D Kaplan, Prafulla Dhariwal, Arvind Neelakantan, Pranav Shyam, Girish Sastry, Amanda Askell, et al.  Language models are few-shot learners.  Advances in neural information processing systems, 33:1877–1901, 2020.
  * [7] Long Ouyang, Jeffrey Wu, Xu Jiang, Diogo Almeida, Carroll Wainwright, Pamela Mishkin, Chong Zhang, Sandhini Agarwal, Katarina Slama, Alex Ray, et al.  Training language models to follow instructions with human feedback.  Advances in Neural Information Processing Systems, 35:27730–27744, 2022.
  * [8] Sébastien Bubeck, Varun Chandrasekaran, Ronen Eldan, Johannes Gehrke, Eric Horvitz, Ece Kamar, Peter Lee, Yin Tat Lee, Yuanzhi Li, Scott Lundberg, et al.  Sparks of artificial general intelligence: Early experiments with gpt-4.  arXiv preprint arXiv:2303.12712, 2023.
  * [9] Yujia Li, Oriol Vinyals, Chris Dyer, Razvan Pascanu, and Peter Battaglia.  Learning deep generative models of graphs.  arXiv preprint arXiv:1803.03324, 2018.
  * [10] Jiaxuan You, Rex Ying, Xiang Ren, William Hamilton, and Jure Leskovec.  Graphrnn: Generating realistic graphs with deep auto-regressive models.  In International conference on machine learning, pages 5708–5717. PMLR, 2018.
  * [11] Renjie Liao, Yujia Li, Yang Song, Shenlong Wang, Will Hamilton, David K Duvenaud, Raquel Urtasun, and Richard Zemel.  Efficient graph generation with graph recurrent attention networks.  Advances in neural information processing systems, 32, 2019.
  * [12] Chia-Cheng Liu, Harris Chan, Kevin Luk, and AI Borealis.  Auto-regressive graph generation modeling with improved evaluation methods.  In 33rd Conference on Neural Information Processing Systems. Vancouver, Canada, 2019.
  * [13] Sepp Hochreiter and Jürgen Schmidhuber.  Long short-term memory.  Neural computation, 9(8):1735–1780, 1997.
  * [14] Kyunghyun Cho, Bart Van Merriënboer, Caglar Gulcehre, Dzmitry Bahdanau, Fethi Bougares, Holger Schwenk, and Yoshua Bengio.  Learning phrase representations using rnn encoder-decoder for statistical machine translation.  arXiv preprint arXiv:1406.1078, 2014.
  * [15] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin.  Attention is all you need.  Advances in neural information processing systems, 30, 2017.
  * [16] Alec Radford, Karthik Narasimhan, Tim Salimans, Ilya Sutskever, et al.  Improving language understanding by generative pre-training.  OpenAI blog, 2018.
  * [17] Alec Radford, Jeffrey Wu, Rewon Child, David Luan, Dario Amodei, Ilya Sutskever, et al.  Language models are unsupervised multitask learners.  OpenAI blog, 1(8):9, 2019.
  * [18] Will Hamilton, Zhitao Ying, and Jure Leskovec.  Inductive representation learning on large graphs.  Advances in neural information processing systems, 30, 2017.
  * [19] Oriol Vinyals, Samy Bengio, and Manjunath Kudlur.  Order matters: Sequence to sequence for sets.  arXiv preprint arXiv:1511.06391, 2015.
  * [20] Xiaohui Chen, Xu Han, Jiajing Hu, Francisco JR Ruiz, and Liping Liu.  Order matters: Probabilistic modeling of node sequence for graph generation.  arXiv preprint arXiv:2106.06189, 2021.
  * [21] Aleksandar Bojchevski, Oleksandr Shchur, Daniel Zügner, and Stephan Günnemann.  Netgan: Generating graphs via random walks.  In International conference on machine learning, pages 610–619. PMLR, 2018.
  * [22] Qi Liu, Miltiadis Allamanis, Marc Brockschmidt, and Alexander Gaunt.  Constrained graph variational autoencoders for molecule design.  Advances in neural information processing systems, 31, 2018.
  * [23] Tengfei Ma, Jie Chen, and Cao Xiao.  Constrained generation of semantically valid graphs via regularizing variational autoencoders.  Advances in Neural Information Processing Systems, 31, 2018.
  * [24] Carl Yang, Peiye Zhuang, Wenhan Shi, Alan Luu, and Pan Li.  Conditional structure generation through graph variational generative adversarial nets.  Advances in neural information processing systems, 32, 2019.
  * [25] Jaehyeong Jo, Seul Lee, and Sung Ju Hwang.  Score-based generative modeling of graphs via the system of stochastic differential equations.  In International Conference on Machine Learning, pages 10362–10383. PMLR, 2022.
  * [26] Clement Vignac, Igor Krawczuk, Antoine Siraudin, Bohan Wang, Volkan Cevher, and Pascal Frossard.  Digress: Discrete denoising diffusion for graph generation.  International Conference on Learning Representations (ICLR 2023), 2022.
  * [27] Manzil Zaheer, Satwik Kottur, Siamak Ravanbakhsh, Barnabas Poczos, Russ R Salakhutdinov, and Alexander J Smola.  Deep sets.  Advances in neural information processing systems, 30, 2017.
  * [28] Yan Zhang, Jonathon Hare, and Adam Prugel-Bennett.  Deep set prediction networks.  Advances in Neural Information Processing Systems, 32, 2019.
  * [29] Adam R Kosiorek, Hyunjik Kim, and Danilo J Rezende.  Conditional set generation with transformers.  arXiv preprint arXiv:2006.16841, 2020.
  * [30] Nicolas Carion, Francisco Massa, Gabriel Synnaeve, Nicolas Usunier, Alexander Kirillov, and Sergey Zagoruyko.  End-to-end object detection with transformers.  In Computer Vision–ECCV 2020: 16th European Conference, Glasgow, UK, August 23–28, 2020, Proceedings, Part I 16, pages 213–229. Springer, 2020.
  * [31] Roy Jonker and Ton Volgenant.  A shortest augmenting path algorithm for dense and sparse linear assignment problems.  In DGOR/NSOR: Papers of the 16th Annual Meeting of DGOR in Cooperation with NSOR/Vorträge der 16. Jahrestagung der DGOR zusammen mit der NSOR, pages 622–622. Springer, 1988.
  * [32] Zhiqing Sun, Shengcao Cao, Yiming Yang, and Kris M Kitani.  Rethinking transformer-based set prediction for object detection.  In Proceedings of the IEEE/CVF international conference on computer vision, pages 3611–3620, 2021.
  * [33] Gongjie Zhang, Zhipeng Luo, Yingchen Yu, Kaiwen Cui, and Shijian Lu.  Accelerating detr convergence via semantic-aligned matching.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 949–958, 2022.
  * [34] Wengong Jin, Regina Barzilay, and Tommi Jaakkola.  Junction tree variational autoencoder for molecular graph generation.  In International conference on machine learning, pages 2323–2332. PMLR, 2018.
  * [35] Ruibin Xiong, Yunchang Yang, Di He, Kai Zheng, Shuxin Zheng, Chen Xing, Huishuai Zhang, Yanyan Lan, Liwei Wang, and Tieyan Liu.  On layer normalization in the transformer architecture.  In International Conference on Machine Learning, pages 10524–10533. PMLR, 2020.
  * [36] Hugo Touvron, Matthieu Cord, Alexandre Sablayrolles, Gabriel Synnaeve, and Hervé Jégou.  Going deeper with image transformers.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 32–42, 2021.
  * [37] Aric A. Hagberg, Daniel A. Schult, and Pieter J. Swart.  Exploring network structure, dynamics, and function using networkx.  In Gaël Varoquaux, Travis Vaught, and Jarrod Millman, editors, Proceedings of the 7th Python in Science Conference, pages 11 – 15, Pasadena, CA USA, 2008.
  * [38] Chengxuan Ying, Tianle Cai, Shengjie Luo, Shuxin Zheng, Guolin Ke, Di He, Yanming Shen, and Tie-Yan Liu.  Do transformers really perform badly for graph representation?  Advances in Neural Information Processing Systems, 34:28877–28888, 2021.
  * [39] Jinwoo Kim, Dat Nguyen, Seonwoo Min, Sungjun Cho, Moontae Lee, Honglak Lee, and Seunghoon Hong.  Pure transformers are powerful graph learners.  Advances in Neural Information Processing Systems, 35:14582–14595, 2022.
  * [40] Lucas Beyer, Olivier J Hénaff, Alexander Kolesnikov, Xiaohua Zhai, and Aäron van den Oord.  Are we done with imagenet?  arXiv preprint arXiv:2006.07159, 2020.
  * [41] Ross Wightman, Hugo Touvron, and Hervé Jégou.  Resnet strikes back: An improved training procedure in timm.  arXiv preprint arXiv:2110.00476, 2021.
  * [42] Petar Veličković, Guillem Cucurull, Arantxa Casanova, Adriana Romero, Pietro Lio, and Yoshua Bengio.  Graph attention networks.  arXiv preprint arXiv:1710.10903, 2017.
  * [43] Mingxing Tan and Quoc V. Le.  Efficientnet: Rethinking model scaling for convolutional neural networks, 2020.
  * [44] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition, 2015.
  * [45] Olga Russakovsky, Jia Deng, Hao Su, Jonathan Krause, Sanjeev Satheesh, Sean Ma, Zhiheng Huang, Andrej Karpathy, Aditya Khosla, Michael Bernstein, Alexander C. Berg, and Li Fei-Fei.  ImageNet Large Scale Visual Recognition Challenge.  International Journal of Computer Vision (IJCV), 115(3):211–252, 2015.
  * [46] Alexey Dosovitskiy, Lucas Beyer, Alexander Kolesnikov, Dirk Weissenborn, Xiaohua Zhai, Thomas Unterthiner, Mostafa Dehghani, Matthias Minderer, Georg Heigold, Sylvain Gelly, Jakob Uszkoreit, and Neil Houlsby.  An image is worth 16x16 words: Transformers for image recognition at scale, 2021.
  * [47] Alessandro Achille, Matteo Rovere, and Stefano Soatto.  Critical learning periods in deep neural networks, 2019.

##  Appendix A Additional Experiment Results

In this section, we present additional experiment results from our study,
offering further insights into the aspects and factors contributing to our
proposed framework.

###  A.1 Approximating TSP Solutions Using Latent Sort

To investigate the performance of latent sort algorithms in approximating the
Traveling Salesman Problem (TSP), we present additional results using a
synthetic 2-D dataset in $\mathbb{R}^{2}$. All points in the dataset are
generated using a uniform distribution, and the objective is to find the
shortest path that traverses all the given points. Table 5 provides a summary
of the experiment outcomes.

Table 5: Approximating TSP solutions using latent sort algorithm. We report the percentage of all possible paths traversing all points (generated using brute-force enumeration) that are longer than the solution predicted by the latent sort algorithm. The results are reported in “mean $\pm$ standard deviation” from 10 random runs of generating $N$ random points and running latent sort. | N = 5 | N = 6 | N = 7 | N = 8
---|---|---|---|---
LS w. LGP | $93.17\%\pm 11.77\%$ | $92.00\%\pm 8.69\%$ | $97.01\%\pm 6.10\%$ | $99.59\%\pm 0.34\%$
LS w.o. LGP | $82.17\%\pm 13.12\%$ | $90.44\%\pm 11.44\%$ | $94.67\%\pm 6.64\%$ | $96.61\%\pm 5.14\%$

Figure 6: Examples solutions of of TSP problems ($N=8$) with randomly
initialized points. It includes approximated solutions using the latent sort
algorithm and the exact solutions using brute force search.

The results consistently demonstrate that employing latent sort with the
latent gradient penalty (LGP) achieves higher approximation rates compared to
using latent sort without LGP. This finding supports the conclusions discussed
in Section 5.2 regarding the effectiveness of LGP in the latent sort algorithm
for approximating TSP solutions. Note that the methods in general appear to
perform better for a higher number of points, as there are more possible paths
in total, leading to a higher percentage of paths being longer than the one
predicted by latent sort.

###  A.2 Semantic Scene Graph Generation on OpenPSG

Here we provide some examples of the results (see Figure 7) for the semantic
scene graph generation task on the OpenPSG dataset.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 7: Randomly sampled results for semantic scene graph generation on the
OpenPSG dataset. Each image consists of the input image (left), ground-truth
labels (middle), and predictions generated by GraphART using the latent sort
algorithm (right). The samples show the label noises in the dataset.

We can see that GraphART generally tends to predict slightly smaller scene
graphs compared to the ground-truth labels. This model behavior can be
attributed to the presence of label noise in the OpenPSG dataset, although the
dataset is already known to have less noise compared to other scene graph
datasets [3]. Due to this issue in the training data, GraphART becomes
cautious in predicting triplets with low confidence. Instead, it shows a
preference for the most common objects and predicates, as well as the
prevalent patterns present in the dataset. Consequently, the model tends to
predict the stop token earlier than expected, resulting in predicted graph
sequences that are likely to include only the most common objects and
predicates. As a result, the prediction sequences conclude before reaching the
size of the ground-truth labels (see Figure 7).

We believe that this issue does not seem to be specific to GraphART alone, but
rather a common problem encountered by autoregressive models when dealing with
label noise in the dataset. Interestingly, we observed that the model’s
predictions are occasionally more reasonable and logically consistent than the
ground-truth labels. For example, in Figure 7(a), the model predicts person -
holding - teddy bear, which is more accurate than person - attached to - teddy
bear in the ground-truth labels. In Figure 7(d), the model predicts window-
other - in - building-other-merged, which is missing from the ground-truth
labels, and the labels have a lot of duplicated relations for person. In
Figure 7(e), the model predicts book - on - shelf, which is more accurate than
book - hanging from - wall-other-merged in the ground-truth labels. We believe
that for future research, it is necessary to develop scene graph datasets with
cleaner labels in order to enable fair benchmarking and evaluate the
performance of different models more accurately.

##  Appendix B Additional Information For Experiment Setups

In this section, we will provide additional details regarding the experiment
setups, as well as provide details of procedures for conducting experiments
involving the models and tasks discussed in the paper for full reproducibility
of our results. Note that, the amount of details provided here are limited
comparing to the released code 222https://anonymous.4open.science/r/ordering-
in-graph-gen/.

###  B.1 Additional Information For Datasets

Scene Graph Generation: We use the OpenPSG dataset [3] which contains 49k
annotated images, and we follow the same train-test split as [3] , which has
around 45k for training and 4k for testing. Our problem setup for scene graph
generation is slightly different from most of the existing works on scene
graph generation [2, 4, 3]. In particular, we consider semantic-level scene
graphs where every node represents a generic object and is not associated to a
particular localized (or segmented) object. For example, a “person1 - beside -
person2” triplet will be converted into “person - beside - person” in our
problem setup. We chose this problem setup because our main focus is to
investigate the effect of different ordering methods on graph generation
rather than striving for state-of-the-art (SOTA) performances. As such, we
avoid complex segmentation settings that could potentially divert our focus.
However, this doesn’t imply an easier task. In fact, our model must infer the
number of objects implicitly without resorting to object-level granularity.
This contrasts with methods like PSGTR and PSGFormer used in conventional
scene graph generation tasks, where the notion of objects is explicitly
provided for training via segmentation masks.

We adopt the same data preprocessing pipelines used for training PSGTR and
PSGFormer, except that we do not use the panoptic segmentation labels for
training. We also have a simpler cropping function instead of the
sophisticated multi-size cropping used in PSGTR and PSGFormer, and our
cropping generates smaller images that result in faster training. More
specifically, we pad the input images to $384\times 384$, where the shorter
side will be padded by 0 values to make the image squared (equal aspect
ratio). In contrast, the PSGTR and PSGFormer uses $640\times 384$ or
$1333\times 640$. The difference in input size may be the factor that GraphART
performs worse in identifying small objects. Then, we parse the object-level
scene graphs into semantic-level scene graphs involving 133 object classes and
56 predicate classes. While state-of-the-art (SOTA) methods for this problem,
such as PSGTR [3] and PSGFormer [3], use additional panoptic segmentation
labels for training, we do not use this additional supervision in our
framework. This is done to demonstrate the capability of autoregressive models
to directly generate scene graphs from images.

Note that the central goal of this paper is not to find a superior vision
model that achieves state-of-the-art performances in scene graph generation.
Thus, we did not fully optimize the pipeline and there are still potentials in
GraphART to reach better performances for this particular task (scene graph
generation). We welcome future works to explore the improvement that can be
done using GraphART.

Topological Graph Extraction: The task is based on the Toulouse Road Network
dataset from [1], which consists of 99k image-graph pairs. The dataset
comprises $64\times 64$ grayscale images paired with corresponding topological
graphs. We follow the same train-test split as described in [1], where
approximately 80k image-graph pairs are used for training, and 19k pairs are
used for testing. The dataset is derived from OpenStreetMap and involves the
processing and labeling of raw road network data from the city of Toulouse,
France, which was obtained from GeoFabrik and timestamped in June 2017. Each
image represents a small cropped patch of the city map, and the coordinates
are normalized to [-1, 1]. Upon inspecting sample images from the dataset, we
observed that the graph structures in most patches are generally simple,
reflecting the inherent simplicity of real-world road networks.

To assess the models’ performance under more challenging conditions, we
decided to create a synthetic Planar Graph dataset from scratch. This dataset
is generated by utilizing Delaunay triangulation and the NetworkX package to
generate random planar graphs. The node coordinates are sampled from a uniform
distribution $\mathcal{U}(0,1)^{2}$. Initially, the graphs contain 15 nodes.
Preprocessing involves collapsing nodes that are too close to each other, with
a threshold of 0.1. Similarly, edges are collapsed if they form an angle less
than 30 degrees. After collapsing nodes and edges, self-loops are removed to
obtain the target graph. Subsequently, we utilize the OpenCV package to render
the target graph into $256\times 256$ binary images, with edges represented by
lines of random widths. To introduce some variability, random artifacts are
added to the images. The graph structures in this dataset are generally more
challenging when compared to the Road Network dataset. During training,
minibatches of images and ground-truth graphs are generated on the fly,
resulting in non-repeating dataset (of 40k samples every epoch) and zero
generalization gap in the trained model. During testing, we test the
performance of the results using 10 random seeds on 1000 generated samples.

Circuit Graph Prediction: We use the Terahertz channelizer dataset333Please
make sure to carefully review and comply with the license permissions and
restrictions (if there are any) before using the dataset [5]:
https://github.com/hehaodele/circuit-gnn consisting of 347k circuit graphs.
This dataset was introduced in the work by He et al. [5], which presented a
novel approach to represent circuits using graphs, where resonators are
represented as nodes and electromagnetic couplings as edges. The authors
mapped the geometric configuration of a circuit to the graph structure, using
node attributes for resonator properties, while edge attributes captured the
relative positions, gap lengths, and shifts between resonators. The presence
of edges in the graph is determined by a distance threshold based on prior
research. We recommend readers to go through the original paper for further
details on the dataset, which includes visual examples of the circuits and
their graph representations.

In our paper, we consider the goal of solving the “inverse” problem, where we
are expected to predict the ground-truth circuit graph that consists of 3 to 6
resonators (or nodes), given desirable electromagnetic (EM) properties (e.g.,
transfer function) as inputs. Each resonator (node) is represented by a
9-dimensional vector. We follow the problem statement of the forward and
inverse problems as the original work [5]. Since the scale of some dimensions
in the resonator vectors is much larger than others, which may dominate the
prediction errors, we normalize the data using min-max normalization per
dimension. We report the errors on the normalized data.

###  B.2 Model Configurations And Baseline Methods

Baselines in Scene Graph Generation: For the scene graph generation task
(Table 3), we consider PSGTR [3] and PSGFormer [3] as state-of-the-art (SOTA)
baselines. We ran the evaluation scripts and trained checkpoints as provided
by Yang et al 444OpenPSG Github Repo: https://github.com/Jingkang50/OpenPSG to
obtain the predictions on the test set, where each node in the predicted scene
graph is grounded by its pixel-accurate segmentation mask in the image. Then,
we convert the segmentation masks based scene graph into semantic scene graph.
Since the output graph size of these methods depends on a confidence
thresholding hyper-parameter, it is challenging to tune it to match the
ground-truth graph sizes. Moreover, these methods are tuned to be biased
towards achieving a good recall score, generally producing much larger graphs
than the ground-truth (GT) labels. Therefore, we downsampled the output
triplets to match the GT numbers and report their average performance over 10
random runs. Note that the PSGTR and PSGFormer utilize panoptic segmentation
labels as additional supervision during training, which we do not use in our
framework. Additionally, we do not use the GT size of the scene graph as an
input parameter in our approach, in contrast to the implementation of PSGTR
and PSGFormer.

Baselines in Topological Graph Extraction: For the road network generation
task (Table 1), we use the GGT [1] and GraphRNN [10] algorithms as baselines.
These models are supervised on adjacency matrices using BCELoss, in contrast
to our purely autoregressive framework that does not generate graph-specific
components. For the planar graph generation task (Table 2), we designed a non-
autoregressive architecture called GraphTR, inspired by the DETR architecture
[30]. GraphTR utilizes bipartite matching loss and does not rely on
autoregressive generation.

Baselines in Topological Graph Prediction: Regarding the circuit prediction
task (Table 4), the original work [5] does not directly solve the inverse
problem from circuit specifications to circuit graphs. Therefore, we
implemented a model adapted from the original design using a graph attention
network (GAT) [42] and fed it with the ground-truth size of output graphs. We
model each circuit as a Graph Attention Network that is initialized to a
projection of the input signals. We also initialize the edge attributes to the
electromagnetic coupling, which can be seen as a distance measure between the
nodes. Our task is to predict the nine “raw” node features for each node of
the circuit graphs. We focus on only predicting the node features as the
original work [5] derives both node attributes and edge attributes from these
node features. Note that the CircuitGNN in the original work does not have the
capability to predict the target graph size (the number of resonators in the
circuit) for a given input, therefore we train different models for different
graph sizes, specifically one model for one size, as handled by the original
work [5].

GraphART And GraphLSTM: We have different model configurations for each task
included in the paper, and GraphLSTM use the same configurations as GraphART
for all the tasks except the decoders are LSTM models instead of
autoregressive transformers. Unless specified, the dropout or drop-path
probabilities are set to 0.

Figure 8: The encoder architecture of the GraphART and GraphLSTM for different
tasks. Figure 9: The decoder architecture of the GraphART for different tasks.
The autoregressive transformer (ART) is based on GPT-2 architecture [17] with
our modified decoder block architecture (see Figure 3(b)). Figure 10: The
decoder architecture of the GraphLSTM for different tasks.

The model configurations include the following details:

  * •

Visual Encoders: Figure 8 presents the common encoder architecture for both
GraphART and GraphLSTM. These models utilize either an Efficient-B0 [43] or a
ResNet-50 [44] pretrained on the ImageNet-1k dataset [45] as the visual
backbone network. The input images are encoded into feature vectors after
flattening the spatial dimensions of the output feature map. The positional
embeddings used in the original Vision Transformer (ViT) [46] are adopted.

  * •

Autoregressive Decoders: Figure 9 and 10 illustrate the detailed architectures
for GraphART and GraphLSTM, respectively. The GraphART model employs a
GPT-2-based architecture with stacked transformer decoder blocks, as depicted
in Figure 3(b). During training, teacher-forcing is used by feeding the
ground-truth sequence to the model. During inference, the model performs
iterative next token prediction until a stop criterion is met. For road
network extraction, planar graph generation, and circuit graph prediction,
stop tokens are represented by vectors with all $-1$ values. If the mean
absolute distance of a generated token to the vector of all $-1$ values is
less than 0.8, it is considered a stop token, leading to the termination of
iterative sequence generation.

We try to control the number of parameters to be similar among different
models for the same task. Please refer to the figures for a visual
representation of the mentioned architectures.

###  B.3 Additional Information For Evaluation Metrics

To evaluate the performance of graph generation results in relation to the
ground truth, it is appropriate to employ metrics that quantify the
dissimilarity between two sets. For the topological graph generation tasks, we
employ the StreetMover distance as was done in [1]. In addition, drawing
inspiration from the Earth Mover’s distance and Hausdorff distance, we
introduce the Edge Mover Distance (EMD) and Edge Hausdorff distance (EHD)
metrics specifically tailored for the "edge-based tokenization" scheme used in
our models (refer to tokenization in Section 4.2).

StreetMover distance (SMD): The SMD is proposed by [1] as an evaluation metric
for evaluating generative models in the context of road networks. It addresses
the limitations of existing metrics by jointly capturing the accuracy of the
reconstructed graph while being invariant to changes in graph representation,
transformations, and size. Unlike pixel-based metrics, SMD considers the
global alignment and magnitude of errors in the reconstructed graphs, making
it more suitable for evaluating road network generation. It overcomes the
drawbacks of other metrics, such as Average Path Length Similarity (APLS),
which require post-processing steps and are designed for different purposes.
SMD is easily interpretable and computationally efficient, leveraging Sinkhorn
iterations for its computation.

Edge Mover Distance (EMD) And Edge Hausdorff Distance (EHD):

Given two finite length sequences of N-dimensional points $X$ and $Y$ (where
$X$ and $Y$ are matrices of N columns), we denote $\mathcal{X}$ and
$\mathcal{Y}$ as the corresponding unordered sets of $X$ and $Y$. In other
words, $\mathcal{X}$ and $\mathcal{Y}$ are sets that contain the rows of the
corresponding matrices, where each row represents an N-dimensional point.
Additionally, $\mathcal{X}$ and $\mathcal{Y}$ are sets that contain all the
rows, and every element within these sets is a row from the matrices.

The EMD between two sequences $X$ and $Y$ (represented by matrices) can be
defined using its unordered sets $\mathcal{X}$ and $\mathcal{Y}$:

$\text{EMD}(\mathcal{X},\mathcal{Y})=\frac{1}{|\mathcal{X}|}\sum_{\bm{x}\in\mathcal{X}}\min_{\bm{y}\in\mathcal{Y}}\rho(\bm{x},\bm{y}),$
(7)

and the EHD can be defined as:

$\text{EHD}(\mathcal{X},\mathcal{Y})=\max\left[\sup_{\bm{x}\in\mathcal{X}}\inf_{\bm{y}\in\mathcal{Y}}\rho(\bm{x},\bm{y}),\sup_{\bm{y}\in\mathcal{Y}}\inf_{\bm{x}\in\mathcal{X}}\rho(\bm{x},\bm{y})\right],$
(8)

where $|\mathcal{X}|$ is the cardinality of set $\mathcal{X}$ representing the
number of points in $\mathcal{X}$. The function $\rho(\bm{x},\bm{y})$
represents the distance between points $\bm{x}$ and $\bm{y}$. This distance
can be calculated using various metrics, such as the Euclidean distance or any
other suitable distance metric for the given problem domain. The first term
$\sup_{\bm{x}\in\mathcal{X}}\inf_{\bm{y}\in\mathcal{Y}}\rho(\bm{x},\bm{y})$
computes the maximum distance from any point in $\mathcal{X}$ to its nearest
neighbor in $\mathcal{Y}$, while the second term
$\sup_{\bm{y}\in\mathcal{Y}}\inf_{\bm{x}\in\mathcal{X}}\rho(\bm{x},\bm{y})$
computes the maximum distance from any point in $\mathcal{Y}$ to its nearest
neighbor in $\mathcal{X}$. The EMD calculates the average of the minimum
distances between each point in set $\mathcal{X}$ and its nearest neighbor in
set $\mathcal{Y}$, whereas the EHD measures the maximum distance between the
closest points of two sets $\mathcal{X}$ and $\mathcal{Y}$. The EHD is the
larger of these two values, representing the maximum discrepancy between the
two sets.

Classification Metrics For Sequences And Sets of Different Sizes:

In traditional classification settings, precision, recall, and F1 score are
typically computed on a per-instance basis. Each individual item (or instance)
in the prediction set is compared to its corresponding item in the true set,
and the counts of true positives, false positives, and false negatives are
aggregated over all instances, often using "macro" averaging. This approach
works well when the prediction set and true set have an equal number of
instances. However, in situations where the prediction set and true set have a
different number of instances, the traditional per-instance calculation
becomes inapplicable. It is crucial to establish new definitions for these
metrics that can appropriately handle the discrepancy in set sizes.

We use the same definition of $\mathcal{X}$ and $\mathcal{Y}$ as unordered
sets and $X$ and $Y$ as predicted and ground-truth sequences. To calculate TP
(True Positives), we compute the cardinality of the intersection of
$\mathcal{X}$ and $\mathcal{Y}$, which represents the count of matches or true
positives. Specifically, we express TP, FP, and FN as follows:

$\displaystyle TP$ $\displaystyle=|\mathcal{X}\cap\mathcal{Y}|$ $\displaystyle
FP$ $\displaystyle=|\mathcal{X}|-TP$ $\displaystyle FN$
$\displaystyle=|\mathcal{Y}|-TP$

where $|\cdot|$ represents the cardinality of a set. Using this new definition
of TP, FP and FN, we define Precision, Recall, and F1 Score between two sets
of different sizes (cardinality) as follows:

Precision is the proportion of correct predictions (TP) to the total number of
predictions in the predicted sets. It reflects the percentage of predictions
that are correct. Mathematically, Precision is defined as:

$Precision=\frac{TP}{TP+FP}=\frac{|\mathcal{X}\cap\mathcal{Y}|}{|\mathcal{X}|}$
(9)

Recall is the proportion of correct predictions (TP) to the total number of
elements in the targe set. It reflects the percentage of ground-truth elements
that are included in the predictions. Mathematically, Recall is defined as:

$Recall=\frac{TP}{TP+FN}=\frac{|\mathcal{X}\cap\mathcal{Y}|}{|\mathcal{Y}|}$
(10)

Similar to the traditional definition, the F1 Score for sets is the harmonic
mean of Precision and Recall. The F1 Score gives equal weight to both
Precision and Recall and provides a balanced measure of model performance.
High Precision and high Recall will yield a high F1 Score. The F1 Score is
defined as:

$F1=\frac{2\times Precision\times
Recall}{Precision+Recall}=\frac{2\times\frac{|\mathcal{X}\cap\mathcal{Y}|}{|\mathcal{X}|}\times\frac{|\mathcal{X}\cap\mathcal{Y}|}{|\mathcal{Y}|}}{\frac{|\mathcal{X}\cap\mathcal{Y}|}{|\mathcal{X}|}+\frac{|\mathcal{X}\cap\mathcal{Y}|}{|\mathcal{Y}|}}$
(11)

Using the above definitions, there are precision, recall and F1 scores for
each pair of predicted and ground-truth graphs. We report the averaged values
over the entire test set.

###  B.4 Training Setups For Graph Generation Models

The configurations for the graph generation models are described as follows:

  * •

Planar Graph (Table 2) & Road Network (Table 1): We train the models using the
Adam (planar graph) or AdamW (road network) optimizer with the following
settings: learning rate = $1.0\times 10^{-4}$, weight decay = 0 (planar graph)
or $1.0\times 10^{-4}$ (road network), and an effective batch size = 300 (150
on 2 GPUs). The warm-up period is set to $10\%$ of the total epochs, and we
employ "cosine" annealing. Sequences are padded to the maximum length of
sequences plus one within a batch, using a padding value of -1. The stop token
is also represented by the value -1. Training 1000 epochs on the planar graph
dataset typically takes approximately 16 to 20 hours on two NVIDIA A100 GPUs
with 128 CPU cores, using 32-bit precision555We observed degradation in model
accuracy when using 16-bit precision in training., and it takes around 6 to 8
hours to train for 200 epochs on the road network dataset using a single
NVIDIA A100 GPU with 16 CPU cores.

  * •

Scene Graph (Table 3): We train the models for 20 epochs using the AdamW
optimizer with the following settings: learning rate = $1.0\times 10^{-4}$,
weight decay = $5.0\times 10^{-4}$, and an effective batch size of 100. The
warm-up period is set to $10\%$ of the total epochs, and we employ "cosine"
annealing. Sequences are padded within a batch to the maximum length of
sequences plus one, with the addition of stop tokens. For each object and
their predicate that are one-hot encoded, an additional dimension is added for
the "stop token" class, resulting in a total of three additional dimensions.
Training typically takes approximately 2 hours on a single NVIDIA A100 GPU
with 128 CPU cores, using 32-bit precision.

  * •

Circuit Graph (Table 4): We train the models using the AdamW optimizer with
the following settings: learning rate = $1.0\times 10^{-4}$, weight decay =
$1.0\times 10^{-4}$, and an effective batch size = 300 (150 on 2 GPUs). The
warm-up period is set to $10\%$ of the total epochs, and we employ "cosine"
annealing. Sequences are padded to the maximum length of sequences plus one
within a batch, using a padding value of -1. The stop token is also
represented by the value -1. Training 200 epochs typically takes approximately
20 to 24 hours on two NVIDIA A100 GPUs with 64 CPU cores, using 32-bit
precision.

###  B.5 Training Latent Sort Encoders

We trained 10 latent sort encoders for each experiment and selected the one
that yielded the best validation performance for the paper. We observed that
the difference in graph generation performance due to different random
initializations of the latent sort encoders is less significant than the
impact of the hyperparameters used in training the encoders. Based on our
observations, here are some key takeaways to guide the design of latent sort
encoders:

  * •

We found that different random initializations of the latent sort encoders can
lead to less than a $5\%$ difference in graph generation performance (measured
in terms of EMD values using GraphART).

  * •

Larger latent sort encoders do not necessarily result in better graph
generation performances, even if they have lower reconstruction errors. In our
experiments, we utilized multilayer perceptrons as encoders and decoders, each
consisting of three hidden layers with 512 neurons and tanh activation..

  * •

The most important hyperparameters we identified include: the coefficient of
LGP in the loss function, the strength of $L_{2}$ regularization, and the
number of training epochs. The trade-off parameter controlling the strength of
LGP regularization is the most important hyperparameter, and we used a range
of 0.01 to 0.1 in our experiments. We found that $L_{2}$ regularization was
unnecessary (thus set to 0) when using LGP. Generally, longer training epochs
do not lead to a degradation in graph generation performance.

  * •

Learning rate scheduling plays a crucial role in the performance of latent
sort encoders. In our experiments, we employ Adam optimizers and cosine
learning rate annealing. The model undergoes a warm-up phase at the beginning
of training, lasting $10\%$ of the total training epochs. During this phase,
the learning rate is exponentially increased from an initial value of
$1.0\times 10^{-5}$. Subsequently, the learning rate gradually decreases until
it reaches a final value of $1.0\times 10^{-8}$.

  * •

We observe that using the undirected graph loss as reconstruction does not
produce significant differences when compare to the L1 reconstruction loss for
the topological graph generation tasks. For ciruit graph generation, we use L1
loss for reconstruction. For the scene graph generation task, we use binary
cross entropy (BCE) loss as used in training the graph generation models.

##  Appendix C Supplementary Theoretical Results And Discussions

This section supplements the theoretical results presented in Section 3 by
providing a more detailed and comprehensive explanation of the overall
framework. The structure of this section is as follows:

  1. C.1

includes the definition and motivations of sorting algorithms from the
perspective of dimensionality reduction (DR), along with other important
notations. It also includes a discussion of how sequence orderings affect
model behavior.

  2. C.2

introduces the definition of ordering ambiguity for sorting algorithms.

  3. C.3

discusses the expected prediction errors caused by ordering ambiguity.

  4. C.4

introduces the idea of the latent sort algorithm. It also includes two
theorems to reveal some of its properties.

  5. C.5

discusses prediction errors from imperfection reconstruction in the latent
sort algorithm. First, it connects the reconstruction error (assuming normal
distribution) with the error distribution of the 1-D latent values. Theorem
C.3 shows how to estimate the 1-D latent distribution (in terms of mean and
variance) from the distribution of reconstruction errors. Then, it discusses
the ordering errors that may be caused by imperfection reconstruction.

  6. C.6

provides a way to simplify the ordering errors from imperfection
reconstruction.

  7. C.7

discusses the desirable properties of an optimal ordering to minimize the
errors from imperfection reconstruction using the simplified expression of
errors from C.6.

###  C.1 Sorting As A Dimensionality Reduction Problem

The Ordering Problem in Set Generation: Since we tokenize each graph into a
set of $N$-dimensional tokens ($N>1$), we denote a graph as a set
$\mathcal{X}$ consisting of $M$ unique points (we also refer to these as
vector tokens) within an $N$-dimensional space. These points are unordered,
which means they don’t follow a specific sequence or pattern.

Our goal is to discover the best possible (or optimal) ordering for these
points, which when used as the target sequence to supervise an autoregressive
model yields optimal performance in generating $\mathcal{X}$. We denote the
sequence obtained by applying the optimal ordering as $Y^{*}$. To visualize
this, imagine $Y^{*}$ as a matrix where every row represents a point in
$\mathcal{X}$. Each point is assigned an index which signifies its position in
this optimal order.

Finding $Y^{*}$ is very challenging since we do not know what orderings the
autoregressive models “prefer” that yield optimal results. In fact, different
model architectures or different applications can have different optimal
orderings. For example, in an RNN, each output at time $t$ depends on the
previous outputs and the current input, forming a sequence of dependencies
(i.e., a chain-like structure). This means information has to flow
sequentially from one step to the next. In contrast, the transformer model’s
attention mechanism allows it to have direct dependencies between all pairs of
positions in the sequence. Therefore, information can flow directly between
any two positions in a sequence, which allows transformers to capture longer
dependencies more efficiently. As a result, the optimal ordering of input
elements for an RNN may show temporal or spatial continuity, given its
inherent sequential processing nature, while a Transformer may not strictly
require such continuity due to its capability to process all elements
concurrently and directly. Moreover, the type of application also plays a
vital role in determining the optimal ordering. Therefore, our goal isn’t
simply to find an ordering that fits all models, but rather to discover an
ordering that is best suited to the model’s architecture and the specific use-
case it is being applied to. The optimal orderings could also potentially be
influenced by factors like the data distribution, the complexity of the
learning task, and even the particular optimization algorithms used during the
training process.

How About End-to-End Learning Based Ordering?

One might naturally consider employing an end-to-end learning approach, in
which the process of identifying an optimal ordering is integrated directly
into the training of the model. This approach could, for example, involve a
Transformer model, known for its permutation invariance at the encoder stage
when no positional embeddings are used. However, at the decoder stage, it
could be designed to attend to a learnable positional embedding, essentially
“learning” the optimal sequence ordering as part of the training process. This
can be an appealing idea since the model is potentially capable of dynamically
adapting to the best ordering for the given data and task.

However, there’s a significant challenge when employing such an end-to-end
training scheme. Neural networks, including both the ordering model and the
autoregressive model, tend to have a critical learning phase [47] in the early
epochs, where they learn essential features and structures of the data. In an
end-to-end setup, the ordering model is consistently changing the sequence
ordering to try to best suit the autoregressive model, while the
autoregressive model is simultaneously trying to adapt to the new ordering
provided by the ordering model. This dual adaptation can create issues with
convergence, since we cannot guarantee the two models will ultimately "agree"
with each other. This can result in a chaotic learning process where the two
models are continuously trying to adapt to each other’s changes without
reaching a stable state. Indeed, in our experiments, we have observed that
such end-to-end learning systems often fail to converge, which motivated us to
discontinue them in our work.

Why We Opt For Dimensionality Reduction (DR)?

A natural question to ask is, how is sorting related to Dimension Reduction
(DR)? We argue that, although sorting in a differentiable manner is
challenging, finding differentiable methods for DR is more straightforward.
This offers new opportunities for creating learning-based methods to order
points in high-dimensional spaces.

In a simpler scenario, let’s say we are working with 1-D representations of
all points in $\mathcal{X}$, which we denote as $\mathcal{H}$, and have
normalized them to a range between 0 and 1. This allows us to define the
problem of finding the optimal ordering as discovering the optimal DR mapping.

Consider a sorting algorithm, $s$, which takes $\mathcal{X}$ and transforms it
into $Y$, the ordered sequence of the points. We can represent this sorting
algorithm with a DR mapping, $f$, which takes points in the original high-
dimensional space $\bm{x}_{i}$ and reduces it to their 1-D value, $h_{i}$.

Note that this DR mapping, $f$, assigns each token $\bm{x_{i}}$ a specific
position, $h_{i}$, in the latent space, regardless of the combinations of
tokens it appears with in $\mathcal{X}$. Therefore, the task of ordering any
arbitrary combination of tokens in $\mathcal{X}$ essentially boils down to
finding a sorted path traversal in their 1-D latent space representations.

Why Local DR Instead of Non-Local DR?

We broadly classify DR methods that can be used for sorting into two
categories: local and non-local. Local DR methods compute the DR mapping based
solely on the value of an individual token, independent of other tokens in the
set. In contrast, non-local DR methods involve the interaction between
different tokens in the set. Non-local methods such as Convolutional Neural
Networks (CNNs) or autoregressive Transformer models inherently assume some
form of structure or relationships among tokens when performing DR.

For the problem at hand, we do not have any prior knowledge about the best
ordering of tokens. Therefore, assuming any form of relationship among tokens,
as non-local DR methods do, might not be appropriate. For example, in the case
of 1D-CNNs, the 1-D representation depends on the arrangement of neighboring
tokens. Since we do not know the correct ordering, deciding what tokens should
be considered neighbors can be ambiguous and arbitrary. This ambiguity will
then be reflected in the 1-D representations, leading to potential
inconsistencies and inaccuracies. Similarly, using autoregressive Transformers
for DR involves an even greater level of complexity since every token can
potentially influence every other token. However, unlike standard
Transformers, the autoregressive nature of these models limits the flow of
information to a directional sequence, which complicates the global
interaction of tokens when the correct ordering is unknown. This suggests an
implicit assumption of a global sequence structure among all tokens in the
set. Yet, without knowing the correct ordering of tokens, determining this
structure becomes a chicken-and-egg problem.

There is one notable exception to the aforementioned issue with non-local DR
methods: the self-attention mechanism of Transformers without using positional
encoding. It treats the token set as a whole, without assuming any specific
order or structure, thanks to its full-attention property. This full-attention
mechanism, which allows the model to attend to all tokens in the set, is the
key to the Transformer’s permutation invariance. However, the full-attention
mechanism has one significant drawback when used in the context of
autoregressive models. The full-attention mechanism is inherently
bidirectional, meaning that it allows for information flow from both "left"
and "right" tokens, effectively modeling
$P(\bm{x_{i}}|\bm{x_{1}},\bm{x_{2}},...,\bm{x_{M}})$. This is contrary to the
autoregressive models that we use for graph generation, which only attend to
"left" side tokens and learn
$P(\bm{x_{i}}|\bm{x_{1}},\bm{x_{2}},...,\bm{x_{i-1}})$ [19]. The bidirectional
nature of the full-attention mechanism means that it has access to future
information that is unavailable to the autoregressive models. Consequently,
the autoregressive models cannot learn the ordering rule represented by the
full-attention Transformer, as they are inherently unable to incorporate
“right" side information. This discrepancy between the two types of models
poses a significant challenge in aligning the orderings learned by them, and
thus limits the utility of using a full-attention Transformer for ordering in
the context of training autoregressive models.

In contrast to all the non-local DR methods mentioned above, local DR methods
like MLPs do not require any assumptions about relationships among tokens. The
1-D representations are computed based solely on the properties of individual
tokens, which are the only certain entities in our context. Therefore, in our
scenario where the correct ordering of tokens is unknown, local DR proves to
be a more reliable and less assumptive approach for obtaining 1-D
representations of tokens.

Restated Definition of the Ordering Problem: Given a set of $M$ unordered
points (or vector tokens) $\mathcal{X}\subseteq\mathbb{R}^{N}$,
$|\mathcal{X}|=M$ in an $N$-dimensional space ($N>1$), we are interested in
finding the “optimal” ordering of points in $\mathcal{X}$, denoted by
$Y^{*}\in\mathbb{R}^{M\times N}$, which when used as the target sequence to
supervise an autoregressive model yields optimal performance in generating
$\mathcal{X}$. Formally, let us denote the ordered sequence $Y^{*}$ as the
matrix
$\left[\bm{x_{1}}^{*},\bm{x_{2}}^{*},...,\bm{x_{M}}^{*}\right]^{\intercal}$,
where every row $\bm{x_{i}}^{*}$ of $Y^{*}$ is a point in $\mathcal{X}$ and
$i$ denotes its sorted index in the optimal ordering. For now, we assume that
such an ordering exists for every set $\mathcal{X}$, and we will discuss the
desired property of it in Section C.7. For simplicity, we assume the 1-D
representations of all the points in $\mathcal{X}$, denoted by $\mathcal{H}$,
are normalized to $[0,1]$.

See 3.1

We can thus formulate the problem of finding the optimal sorting
$s^{*}:\mathcal{X}\rightarrow Y^{*}$ as finding an optimal DR mapping
$f^{*}:\mathbb{R}^{N}\rightarrow\mathbb{R}$,
$h_{i}^{*}=f^{*}(\bm{x_{i}}^{*})$, and $\forall i\leq j\leq
M,f^{*}(\bm{x_{i}}^{*})\leq f^{*}(\bm{x_{j}}^{*})$.

Modeling the Errors of a Sorting Algorithm: Assuming we know the target
sequence $Y^{*}$ and we only want to measure the alignment between a sequence
$Y$ and $Y^{*}$, we introduce a probability matrix $P\in\mathbb{R}^{M\times
M}$, where each entry $p_{ij}$ within $P$ is designed to represent the
probability of $\bm{x_{i}}$ ($i$-th element in $Y$) being $\bm{x_{j}}^{*}$
($j$-th element in $Y^{*}$). In this way, the expected value of the ordering
$Y$, denoted by $\mathbb{E}[Y]$, can be computed as $PY^{*}$. The matrix $P$
serves as a “soft" permutation matrix that transforms $Y^{*}$ into $Y$. This
error representation is generic for analyzing different types of errors that
may emerge in a sorting algorithm.

For simplicity, we use Frobenius norm along with the $P$ matrix
representation:

$\mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}=\big{\|}\mathbb{E}[Y]-Y^{*}\big{\|}_{F}^{2}=\big{\|}PY^{*}-Y^{*}\big{\|}_{F}^{2},$
(12)

where $\|\cdot\|_{F}$ denotes the Frobenius norm. We show that the
minimization over $\mathcal{E}\big{(}\mathbb{E}[Y],Y^{*}\big{)}$ can be solved
by an equivalent problem of minimizing $\|P-I_{M}\|_{F}^{2}$ (see Lemma D.1).

###  C.2 Definition of The Ordering Ambiguity

We refer to the problem where multiple points in $\mathcal{X}$ are assigned
the same 1-D representation value as ordering ambiguity. It is problematic
since the ordering of tokens with the same 1-D latent value becomes undefined
if no further tie-breaking schemes are defined.

Why Ordering Ambiguity Is Almost Inevitable to Local DR? As we discussed
earlier, local DR methods are more suitable for sorting tokens used for
training the autoregressive models. Local DR methods are characterized by the
computation of dimensionality reduction solely based on the value of an
individual token, without considering the context or interactions among other
tokens within the set. Since high-to-low-dimensional mappings are mostly
surjective (each element in the target space has a pre-image in the domain)
and not bijective (one-to-one correspondence), it is very likely that multiple
N-D points share the same 1-D representation value (see Remark E.3). This
property is a fundamental source of ordering ambiguity, making it an almost
inevitable occurrence in local DR methods. Non-local DR methods, on the other
hand, can dynamically adapt to alterations in other data points, thereby
circumventing collisions in the low-dimensional space. This adaptability can
potentially prevent the occurrence of ordering ambiguity.

The ordering ambiguity of a sorting algorithm can be characterized by an
ordering ambiguity set, which is defined as follows: See 3.2

###  C.3 The Errors From Ordering Ambiguity

In the presence of ordering ambiguity, the resulting sorted sequence $Y$ may
not be an accurate representation of the target sorted sequence $Y^{*}$. We
can use the above $P$ matrix to quantify the error introduced by ordering
ambiguity.

Assuming that $Y=\mathbb{E}[Y]=Y^{*}$ in the absence of ordering ambiguity, we
assume the autoregressive model will converge to the expectation of the
ambiguity set $\mathcal{A}_{i}$ for each point $\bm{x_{i}}$ in
$\mathcal{A}_{i}$, i.e.,

$\bm{x_{i}}=\mathbb{E}_{\bm{x_{j}}\in\mathcal{A}_{i}}\left[\bm{x_{j}}\right],$
(13)

When ordering ambiguity is present, we assume the probabilities are spread
over the members of the corresponding ambiguity set $\mathcal{A}_{i}$
uniformly, i.e.,

$p_{ij}=\begin{cases}\frac{1}{|\mathcal{A}_{i}|}&\text{if}\
\bm{x_{j}}^{*}\in\mathcal{A}_{i}\\\ 0&\text{otherwise}\end{cases}$ (14)

Therefore we have:

$\mathbb{E}[Y]-Y^{*}=PY^{*}-Y^{*}=\begin{bmatrix}\left(\frac{1}{|\mathcal{A}_{1}|}\sum_{\bm{x_{j}}^{*}\in\mathcal{A}_{1}}\bm{x_{j}}^{*}\right)^{\intercal}-{\bm{x_{1}}^{*}}^{\intercal}\\\
\left(\frac{1}{|\mathcal{A}_{2}|}\sum_{\bm{x_{j}}^{*}\in\mathcal{A}_{2}}\bm{x_{j}}^{*}\right)^{\intercal}-{\bm{x_{2}}^{*}}^{\intercal}\\\
\vdots\\\
\left(\frac{1}{|\mathcal{A}_{M}|}\sum_{\bm{x_{j}}^{*}\in\mathcal{A}_{M}}\bm{x_{j}}^{*}\right)^{\intercal}-{\bm{x_{M}}^{*}}^{\intercal}\end{bmatrix}$
(15)

Now we can calculate the error $\mathcal{E}\big{(}\mathbb{E}[Y],Y^{*}\big{)}$
by taking the squared Frobenius norm of the difference between $PY^{*}$ and
$Y^{*}$:

$\displaystyle\mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}$
$\displaystyle=\|PY^{*}-Y^{*}\|_{F}^{2}$ (16)
$\displaystyle=\sum_{i=1}^{M}\Bigg{\|}\left(\frac{1}{|\mathcal{A}_{i}|}\sum_{\bm{x_{j}}^{*}\in\mathcal{A}_{i}}\bm{x_{j}}^{*}\right)-\bm{x_{i}}^{*}\
\Bigg{\|}^{2}$ (17)
$\displaystyle=\sum_{i=1}^{M}\Bigg{\|}\left(\frac{1}{|\mathcal{A}_{i}|}\sum_{\bm{x_{j}}\in\mathcal{A}_{i}}\bm{x_{j}}\right)-\bm{x_{i}}\
\Bigg{\|}^{2}$ (18)

Here, $\|\cdot\|_{F}$ denotes the Frobenius norm, which is the square root of
the sum of the squared elements of a matrix or vector. Note that we leverage
the summation in the formula to make the error term independent from the
definition of the target sequence $Y^{*}$.

The Impact of Ordering Ambiguity on Autoregressive Models: The ordering of
tokens within an ordering ambiguity set is determined using a uniform
distribution in the above analysis, which might seem confusing, given that the
DR methods we employ are deterministic and, theoretically, should always yield
a consistent order rather than a “random" one. However, the randomness can be
understood from the perspective of the training data input into the
autoregressive models.

In training autoregressive models, our objective is to enable the model to
learn generalizable rules automatically from the training data, including
those determining the ordering of tokens. Suppose the ordering varies
dramatically among different sequences in the training data with similar
patterns. In that case, the model may struggle to identify the factors causing
these discrepancies in ordering. Instead, it may resort to learning an average
of possible orderings (see Figure 11), a phenomenon that is commonly observed
when the models have difficulty fitting the data.

The issue of ordering ambiguity compounds this problem. Without a consistent
and reliable ordering rule for tokens within the ambiguity set, different
samples in the dataset might exhibit varied ordering. Such inconsistency makes
it challenging for autoregressive models to effectively learn the underlying
ordering rule. As such, the presence of ordering ambiguity in the training
data can significantly hamper the performance of autoregressive models in
terms of sequence prediction.

Figure 11: Illustration of Ordering Ambiguity Problem In Autoregressive Graph
Generation.

###  C.4 The Latent Sort Algorithm

The main idea of the latent sort algorithm is to use an auto-encoder to learn
a dimensionality reduction mapping, which enjoys high degree of flexibility
when paired with different loss functions to achieve different goals. In our
design, the auto-encoder consists of two multi-layer perceptrons (MLPs) as the
encoder and decoder, respectively, denoted as $f_{e}$ and $f_{d}$. We can use
the encoder for sorting $\bm{x_{i}}$ and the ordering is given by the encodede
1-D latent representation $h_{i}=f_{e}(\bm{x_{i}})$.

The Ordering Ambiguity of Latent Sort: Since the auto-encoder uses
deterministic neural networks (e.g., MLPs), the encoder and decoder are both
surjective (see Remark E.4). The auto-encoder uses reconstruction loss during
training. If the auto-encoder has perfection reconstruction, its encoder and
decoder are both bijective (see Remark E.5) and there would be no ordering
ambiguity. However, we know it cannot be true for most cases (see Remark E.3).
Therefore, when there exists no bijective mapping between $\mathcal{X}$ and
$\mathcal{H}$, ordering ambiguity in inevitable in latent sort. When multiple
points in $\mathcal{X}$ are mapped to the same point in $\mathcal{H}$, we can
use Equation (2) to estimate the error due to ordering ambiguity.

Some Theoretical Properties of The Latent Sort Autoencoders:

For latent sort, since we adopt auto-encoder for dimensionality reduction,
there are more theoretical properties due to the fact that the neural networks
are often Lipschitz continuous.

Assume that the auto-encoder is initialized using a random distribution with a
small variance, and the model is trained with proper regularization, such that
after training, the encoder and decoder are Lipschitz continuous with
constants $K_{e}$ and $K_{d}$ respectively, and $K_{e}K_{d}\geq 1$. Assume for
all $\bm{x_{i}}\in\mathcal{X}$ and its reconstructed
$\hat{\bm{x_{i}}}=f_{d}(f_{e}(\bm{x_{i}}))$, $\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|$
has a constant upper bound $B$. Assuming $\forall\bm{x_{i}}\in\mathcal{X}$,
there always exist at least one point $h_{i}^{*}\in R$ such that
$f_{d}(h_{i}^{*})=\bm{x_{i}}$, we have the following theorems on the
relationship between the latent representation error and the reconstruction
error:

######  Theorem C.1 (Bounded Original Distance).

Let $\bm{x_{i}}$ and $\bm{x_{j}}$ be two distinct points in
$\mathcal{X}\subseteq\mathbb{R}^{N}$, and let $h_{i}=f_{e}(\bm{x_{i}})$ and
$h_{j}=f_{e}(\bm{x_{j}})$ be their corresponding latent representations, and
$\hat{\bm{x_{i}}}=f_{d}(h_{i})$ and $\hat{\bm{x_{j}}}=f_{d}(h_{j})$ be their
corresponding reconstructions. If $|h_{i}-h_{j}|\leq\epsilon$, where constant
$\epsilon>0$, then the distance $\|\bm{x_{i}}-\bm{x_{j}}\|$ has an upper bound
$2B+{K_{d}}\epsilon$.

###### Proof.

From the assumptions, we know that the decoder $f_{d}$ is $K_{d}$-Lipschitz
continuous. Let $f_{d}(h_{i}^{*})=\bm{x_{i}}$ and
$f_{d}(h_{j}^{*})=\bm{x_{j}}$, where $h_{i}^{*}$ and $h_{j}^{*}$ are the two
1-D representation values that could result in perfect reconstruction. Since
the reconstruction error has a constant upper bound $B$ for all points in
$\mathcal{X}$, we get

$\displaystyle\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|$
$\displaystyle=\|\bm{x_{i}}-f_{d}(h_{i})\|$
$\displaystyle=\|f_{d}(h_{i}^{*})-f_{d}(h_{i})\|$ $\displaystyle\leq
K_{d}|h_{i}^{*}-h_{i}|\leq B$ $\displaystyle\|\bm{x_{j}}-\hat{\bm{x_{j}}}\|$
$\displaystyle=\|\bm{x_{j}}-f_{d}(h_{j})\|$
$\displaystyle=\|f_{d}(h_{j}^{*})-f_{d}(h_{j})\|$ $\displaystyle\leq
K_{d}|h_{j}^{*}-h_{j}|\leq B$

Using $|h_{i}-h_{j}|\leq\epsilon$ and triangle inequality we have:

$\displaystyle\|\bm{x_{i}}-\bm{x_{j}}\|$
$\displaystyle=\|f_{d}(h_{i}^{*})-f_{d}(h_{j}^{*})\|$ $\displaystyle\leq
K_{d}|h_{i}^{*}-h_{j}^{*}|$
$\displaystyle=K_{d}(|h_{i}^{*}-h_{i}+h_{i}-h_{j}+h_{j}-h_{j}^{*}|)$
$\displaystyle\leq K_{d}(|h_{i}^{*}-h_{i}|+|h_{i}-h_{j}|+|h_{j}-h_{j}^{*}|)$
$\displaystyle\leq
K_{d}\left(\frac{1}{K_{d}}\big{\|}f_{d}(h_{i}^{*})-f_{d}(h_{i})\big{\|}+\epsilon+\frac{1}{K_{d}}\big{\|}f_{d}(h_{j}^{*})-f_{d}(h_{j})\big{\|}\right)$
$\displaystyle=\|f_{d}(h_{i}^{*})-f_{d}(h_{i})\|+{K_{d}}\epsilon+\|f_{d}(h_{j}^{*})-f_{d}(h_{j})\|$
$\displaystyle\leq 2B+{K_{d}}\epsilon,$

which gives the upper bound of $\|\bm{x_{i}}-\bm{x_{j}}\|$.

In addition, we need to examine whether the upper bound exists. Since the
encoder is $K_{e}$-Lipschitz continuous everywhere on $\mathcal{X}$ and we
have assumed $K_{e}K_{d}\geq 1$, thus

$\displaystyle|h_{i}-h_{j}|$
$\displaystyle\leq{K_{e}}\|\bm{x_{i}}-\bm{x_{j}}\|$ (19)
$\displaystyle\leq{K_{e}}(2B+{K_{d}}\epsilon)$ (20)
$\displaystyle=2{K_{e}}B+{K_{e}}{K_{d}}\epsilon$ (21)

Since $|h_{i}-h_{j}|\leq\epsilon$ and $2K_{e}B\geq 0$, inequality (21) is
always true for any $\bm{x_{i}}$ and $\bm{x_{j}}$ in $\mathcal{X}$. Therefore,
the upper bound exists as stated in the theorem. ∎

Theorem C.1 establishes a relationship between the distance in the original
space and the distance in the latent representation space learned by an auto-
encoder. It states that if two data points have similar latent
representations, their distance in the original space is also small, with an
upper bound proportional to the Lipschitz constants of the encoder and
decoder.

######  Theorem C.2 (Bounded Latent Distance).

Let $\bm{x_{i}}$ and $\bm{x_{j}}$ be two distinct points in
$\mathcal{X}\subseteq\mathbb{R}^{N}$, and let $h_{i}=f_{e}(\bm{x_{i}})$ and
$h_{j}=f_{e}(\bm{x_{j}})$ be their corresponding latent representations, and
$\hat{\bm{x_{i}}}=f_{d}(h_{i})$ and $\hat{\bm{x_{j}}}=f_{d}(h_{j})$ be their
corresponding reconstructions. If the distance $\|\bm{x_{i}}-\bm{x_{j}}\|$ is
greater than or equal to a constant $D$, and $D\geq 2B$, then the difference
in their latent representations $|h_{i}-h_{j}|$ is larger than a constant
lower bound $(D-2B)/{K_{d}}$.

###### Proof.

Let $f_{d}(h_{i}^{*})=\bm{x_{i}}$ and $f_{d}(h_{j}^{*})=\bm{x_{j}}$. Since
$f_{d}$ is Lipschitz continuous with constant $K_{d}$, we know that

$\displaystyle\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|$ $\displaystyle\leq
K_{d}|h_{i}^{*}-h_{i}|\leq B$ $\displaystyle\|\bm{x_{j}}-\hat{\bm{x_{j}}}\|$
$\displaystyle\leq K_{d}|h_{j}^{*}-h_{j}|\leq B$ $\displaystyle
D\leq\|\bm{x_{i}}-\bm{x_{j}}\|$ $\displaystyle\leq K_{d}|h_{i}^{*}-h_{j}^{*}|$

From the assumptions, we know that the distance $\|\bm{x_{i}}-\bm{x_{j}}\|\geq
D$ and $D\geq 2B$. Using reverse triangle inequality, we get

$\displaystyle|h_{i}-h_{j}|$
$\displaystyle=|(h_{i}^{*}-h_{j}^{*})-(h_{i}^{*}-h_{i})+(h_{j}^{*}-h_{j})|$
(22)
$\displaystyle\geq\Big{|}|h_{i}^{*}-h_{j}^{*}|-|(h_{i}^{*}-h_{i})-(h_{j}^{*}-h_{j})|\Big{|}$
(23)
$\displaystyle\geq\Big{|}|h_{i}^{*}-h_{j}^{*}|-|h_{i}^{*}-h_{i}|-|h_{j}^{*}-h_{j}|\Big{|}$
(24) $\displaystyle\geq\frac{1}{K_{d}}\Big{|}D-B-B\Big{|}$ (25)
$\displaystyle\geq\frac{1}{K_{d}}(D-2B)$ (26)

Thus we obtain the lower bound for $|h_{i}-h_{j}|$.

In addition, we need to examine whether the lower bound exists. Since $f_{e}$
is $K_{e}$-Lipschitz continuous and $K_{e}K_{d}\geq 1$, we have

$\displaystyle|h_{i}-h_{j}|$
$\displaystyle\leq{K_{e}}\|\bm{x_{i}}-\bm{x_{j}}\|$ (27)
$\displaystyle\frac{1}{K_{d}}(D-2B)$
$\displaystyle\leq{K_{e}}\|\bm{x_{i}}-\bm{x_{j}}\|$ (28)
$\displaystyle\frac{D-2B}{K_{e}K_{d}}$ $\displaystyle\leq
D\leq\|\bm{x_{i}}-\bm{x_{j}}\|$ (29)

Inequality (29) is always true therefore the lower bound exists as stated in
the theorem. ∎

Theorem C.1 provides a complementary perspective to Theorem C.2. While the
first theorem bounds the distance in the original space based on the distance
in the latent representation space, the second theorem bounds the difference
in the latent representations based on the distance in the original space.
This is useful to estimate the error from latent ambiguity in latent sort,
which we will discuss in the later sections.

Interestingly, it has a connection to ordering ambiguity. If the ordering
ambiguity problem is serious and the ambiguity sets are large, then the upper
bound $B$ for the reconstruction error is expected to be large as well. In
such cases, $D-2B$ in the original space can shrink to close to zero, which
inversely supports the existence of ordering ambiguity.

###  C.5 The Errors From Imperfect Reconstruction

Previously, when discussing the errors from ambiguity, we assumed that
$Y=Y^{*}$ in the absence of ordering ambiguity, even though it is not
realistic since the auto-encoder used in latent sort is not bijective in most
cases. To quantify the errors from imperfect reconstruction, we assume no
ordering ambiguity in this section.

We model the imperfect reconstruction by assigning a non-zero variance to the
reconstructed $\hat{\bm{x_{i}}}$, and we assume the mean value equals to the
target $\bm{x_{i}}$. We assume $\hat{\bm{x_{i}}}$ follows Gaussian
distribution (see Remark E.8), and larger variances indicate worse
reconstruction, and vice versa.

Figure 12 shows a schematic plot of an autoencoder that has imperfect
reconstruction. The ideal autoencoder that has perfect reconstruction are
shown as the dashed lines in the plot, indicating the bijectivity of the
autoencoder. However, when the reconstruction is not perfect, we assume that
the $\hat{h_{i}}$ is a sample from the distribution that are centered around
$h_{i}$.

Figure 12: A schematic representation of the autoencoder under imperfect
reconstruction.

Now, we want to understand, when observing the reconstruction error in the
original N-D space, can we estimate the distribution of the 1-D latent
distribution? This motivates the following theorem.

######  Theorem C.3 (Latent Distribution Under Imperfect Reconstruction).

Let $\bm{\sigma_{\hat{\bm{x_{i}}}}^{2}}\in\mathbb{R}^{N}$ be the component-
wise variance vector of the reconstructed point $\hat{\bm{x_{i}}}$, and
$\bm{x_{i}}=f_{d}(h_{i})$, $\hat{\bm{x_{i}}}=f_{d}(\hat{h_{i}})$. Assume that
the encoder $f_{e}$ and decoder $f_{d}$ are Lipschitz continuous on
$\mathcal{X}$ with constants $K_{e}$ and $K_{d}$ respectively. If
$\hat{h_{i}}$ also follows a normal distribution, the following bounds hold
for the mean $\mu_{\hat{h_{i}}}$ and variance $\sigma_{\hat{h_{i}}}^{2}$ of
$\hat{h_{i}}$:

$\displaystyle\frac{\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}}{K_{d}}\leq\big{|}\mu_{\hat{h_{i}}}-h_{i}\big{|}\leq
K_{e}\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}$ (30)
$\displaystyle\frac{1}{K_{d}^{2}}\
\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}\leq\sigma_{\hat{h_{i}}}^{2}\leq
4K_{e}^{2}\ \big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}$ (31)

###### Proof.

Consider the following chain of inequalities based on the Lipschitz continuity
of $f_{e}$ and $f_{d}$:

$\displaystyle\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|$
$\displaystyle=\|f_{d}(h_{i})-f_{d}(\hat{h_{i}})\|$ (32) $\displaystyle\leq
K_{d}\big{|}h_{i}-\hat{h_{i}}\big{|}.$ (33)

Now, we can also write the Lipschitz continuity for the encoder $f_{e}$:

$\displaystyle\big{|}h_{i}-\hat{h_{i}}\big{|}$
$\displaystyle=|f_{e}(\bm{x_{i}})-f_{e}(\hat{\bm{x_{i}}})|$ (34)
$\displaystyle\leq K_{e}\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|.$ (35)

Combining the two inequalities, we get:

$\displaystyle\frac{\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|}{K_{d}}\leq\big{|}h_{i}-\hat{h_{i}}\big{|}$
$\displaystyle\leq K_{e}\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|$ (36)
$\displaystyle\text{s.t.\ }\ K_{e}K_{d}$ $\displaystyle\geq 1$ (37)

Taking the squared expectation on both sides:

$\frac{\mathbb{E}\Big{[}\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|^{2}\Big{]}}{K_{d}^{2}}\leq\mathbb{E}\Big{[}\big{|}h_{i}-\hat{h_{i}}\big{|}^{2}\Big{]}\leq
K_{e}^{2}\cdot\mathbb{E}\Big{[}\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|^{2}\Big{]}$
(38)

Considering that $\mathbb{E}[\hat{\bm{x_{i}}}]=\bm{x_{i}}$ we have:

$\displaystyle\bm{\sigma_{\hat{\bm{x_{i}}}}^{2}}$
$\displaystyle=\operatorname{Var}(\hat{\bm{x_{i}}})$ (39)
$\displaystyle=\mathbb{E}\Big{[}(\hat{\bm{x_{i}}}-\mathbb{E}\big{[}\hat{\bm{x_{i}}}\big{]})^{2}\Big{]}$
(40) $\displaystyle=\mathbb{E}\big{[}(\hat{\bm{x_{i}}}-\bm{x_{i}})^{2}\big{]}$
(41)

And we have

$\displaystyle\mathbb{E}\Big{[}\|\bm{x_{i}}-\hat{\bm{x_{i}}}\|^{2}\Big{]}$
$\displaystyle=\mathbb{E}\Big{[}(\hat{\bm{x_{i}}}-\bm{x_{i}})^{\intercal}(\hat{\bm{x_{i}}}-\bm{x_{i}})\Big{]}$
(42) $\displaystyle=\Big{\|}\
\mathbb{E}\big{[}(\hat{\bm{x_{i}}}-\bm{x_{i}})^{2}\big{]}\ \Big{\|}$ (43)
$\displaystyle=\sum_{j=1}^{N}\sigma_{\hat{\bm{x_{i}}},j}^{2}=\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}$
(44)

where $\sigma_{\hat{\bm{x_{i}}},j}^{2}$ is the $j$-th component of the
variance vector $\bm{\sigma_{\hat{\bm{x_{i}}}}^{2}}$. Therefore we have the
upper bound and lower bound for the mean value of $\hat{h_{i}}$:

$\displaystyle\frac{\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}}{K_{d}^{2}}$
$\displaystyle\leq\mathbb{E}\Big{[}\big{|}h_{i}-\hat{h_{i}}\big{|}^{2}\Big{]}\leq
K_{e}^{2}\cdot\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}$ (45)
$\displaystyle\frac{\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}}{K_{d}}$
$\displaystyle\leq\big{|}\mu_{\hat{h_{i}}}-h_{i}\big{|}\leq
K_{e}\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}$ (46)

Now, consider the variance of $\hat{h_{i}}$. From the Lipschitz continuity of
the encoder and decoder, we can use triangle inequality:

$\displaystyle\sigma_{\hat{h_{i}}}^{2}=\operatorname{Var}(\hat{h_{i}})$
$\displaystyle=\mathbb{E}\Big{[}(\hat{h_{i}}-\mathbb{E}\big{[}\hat{h_{i}}\big{]})^{2}\Big{]}=\mathbb{E}\big{[}(\hat{h_{i}}-\mu_{\hat{h_{i}}})^{2}\big{]}$
(47)
$\displaystyle=\mathbb{E}\Big{[}\big{|}\hat{h_{i}}-h_{i}+h_{i}-\mu_{\hat{h_{i}}}\big{|}^{2}\Big{]}$
(48)
$\displaystyle\leq\mathbb{E}\Big{[}\big{|}f_{e}(\hat{\bm{x_{i}}})-f_{e}(\bm{x_{i}})\big{|}^{2}\Big{]}+\big{|}\mu_{\hat{h_{i}}}-h_{i}\big{|}^{2}$
(49) $\displaystyle\ \ \ \
+2\mathbb{E}\Big{[}\big{|}f_{e}(\hat{\bm{x_{i}}})-f_{e}(\bm{x_{i}})\big{|}\Big{]}\cdot\big{|}\mu_{\hat{h_{i}}}-h_{i}\big{|}$
(50) $\displaystyle\leq K_{e}^{2}\
\mathbb{E}\Big{[}|\bm{x_{i}}-\hat{\bm{x_{i}}}|^{2}\Big{]}+K_{e}^{2}\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}$
(51) $\displaystyle\ \ \ \ +2K_{e}^{2}\
\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}\cdot\mathbb{E}\Big{[}|\bm{x_{i}}-\hat{\bm{x_{i}}}|\Big{]}$
(52) $\displaystyle=4K_{e}^{2}\
\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}$ (53)

Conversely, using the Lipschitz continuity of the decoder, we can derive the
lower bound of the variance:

$\displaystyle\sigma_{\hat{h_{i}}}^{2}=\operatorname{Var}(\hat{h_{i}})$
$\displaystyle=\mathbb{E}\Big{[}(\hat{h_{i}}-\mathbb{E}\big{[}\hat{h_{i}}\big{]})^{2}\Big{]}=\mathbb{E}\big{[}(\hat{h_{i}}-\mu_{\hat{h_{i}}})^{2}\big{]}$
(54)
$\displaystyle=\mathbb{E}\Big{[}\big{|}(\hat{h_{i}}-h_{i})-(\mu_{\hat{h_{i}}}-h_{i})\big{|}^{2}\Big{]}$
(55)

We can apply the Cauchy-Schwarz inequality to above, which gives:

$\mathbb{E}\left[\left(\hat{h}_{i}-h_{i}\right)\left(\hat{h}_{i}-\mu_{\hat{h}_{i}}\right)\right]^{2}\leq\mathbb{E}\left[\left(\hat{h}_{i}-h_{i}\right)^{2}\right]\mathbb{E}\left[\left(\hat{h}_{i}-\mu_{\hat{h}_{i}}\right)^{2}\right]$
(56)

We can rearrange the inequality:

$\mathbb{E}\left[\left(\hat{h}_{i}-\mu_{\hat{h}_{i}}\right)^{2}\right]\geq\frac{\mathbb{E}\left[\left(\hat{h}_{i}-h_{i}\right)\left(\hat{h}_{i}-\mu_{\hat{h}_{i}}\right)\right]^{2}}{\mathbb{E}\left[\left(\hat{h}_{i}-h_{i}\right)^{2}\right]}$
(57)

Using the triangle inequality, we have:

$\displaystyle\Big{[}\
\mathbb{E}\big{[}(\hat{h_{i}}-h_{i})(h_{i}-\mu_{\hat{h_{i}}})\big{]}\Big{]}^{2}$
$\displaystyle\leq\Big{[}\
\mathbb{E}\big{[}|\hat{h_{i}}-h_{i}||h_{i}-\mu_{\hat{h_{i}}}|\big{]}\Big{]}^{2}$
(58)
$\displaystyle\leq\mathbb{E}\Big{[}(\hat{h_{i}}-h_{i})^{2}\Big{]}\mathbb{E}\Big{[}(h_{i}-\mu_{\hat{h_{i}}})^{2}\Big{]}$
(59)

Now, we can use this bound in the expression for the lower bound:

$\displaystyle\mathbb{E}\left[\left(\hat{h}_{i}-\mu_{\hat{h}_{i}}\right)^{2}\right]$
$\displaystyle\geq\frac{\mathbb{E}\left[\left(\hat{h}_{i}-h_{i}\right)\left(\hat{h}_{i}-\mu_{\hat{h}_{i}}\right)\right]^{2}}{\mathbb{E}\left[\left(\hat{h}_{i}-h_{i}\right)^{2}\right]}$
(60)
$\displaystyle\geq\frac{\mathbb{E}\Big{[}\big{(}\hat{h_{i}}-h_{i}\big{)}^{2}\Big{]}\mathbb{E}\Big{[}\big{(}h_{i}-\mu_{\hat{h_{i}}}\big{)}^{2}\Big{]}}{\mathbb{E}\left[\left(\hat{h}_{i}-h_{i}\right)^{2}\right]}$
(61)
$\displaystyle=\mathbb{E}\Big{[}\big{(}h_{i}-\mu_{\hat{h_{i}}}\big{)}^{2}\Big{]}$
(62)

Note that
$\mathbb{E}\big{[}(h_{i}-\mu_{\hat{h_{i}}})^{2}\big{]}=\big{|}\mu_{\hat{h_{i}}}-h_{i}\big{|}^{2}\geq\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}/K_{d}^{2}$.
Thus, we have obtained the bounds for the variance of $\hat{h_{i}}$ as
follows:

$\frac{1}{K_{d}^{2}}\cdot\big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}\leq\sigma_{\hat{h_{i}}}^{2}\leq
4K_{e}^{2}\ \big{\|}\bm{\sigma_{\hat{\bm{x_{i}}}}}\big{\|}^{2}$ (63)

∎

Theorem C.3 provides bounds on the mean and variance of the 1-D representation
in terms of the Lipschitz continuity of the encoder and decoder and the
component-wise variance of the reconstruction error. To build theoretical
bounds for the error from imperfect reconstruction, we need to consider how
the latent distribution affect the sorting results.

How Does Imperfect Reconstruction Affect Ordering Results?

We have discussed the distribution of 1-D representation values when the
reconstruction is not perfect. Connecting back to the ordering results, when
two points $\bm{x_{i}}$ and $\bm{x_{j}}$ follow an ordering where $i<j$, we
want the learned 1-D latent to always satisfy $\hat{h_{i}}<\hat{h_{j}}$.
However, since their 1-D latent representations have non-zero variances, it is
possible that $\hat{h_{i}}\geq\hat{h_{j}}$, resulting in a swap between
$\bm{x_{i}}$ and $\bm{x_{j}}$ and leading to an error in the sorted sequence.
Therefore we want to quantitatively measure how much error can we expect from
imperfect reconstruction.

Assuming for all $\bm{x_{i}}^{*}$ in the target sequence $Y^{*}$, we have
$f_{e}(\bm{x_{i}}^{*})=\hat{h_{i}^{*}}\sim\mathcal{N}(\mu_{\hat{h_{i}^{*}}},\sigma_{\hat{h_{i}^{*}}}^{2})$
independent of each other (using local DR). For each pair of $\hat{h_{i}^{*}}$
and $\hat{h_{j}^{*}}$ where $i\neq j$, the probability
$P(\hat{h_{i}^{*}}>\hat{h_{j}^{*}})$ is given by (see Remark E.2):

$P(\hat{h_{i}^{*}}>\hat{h_{j}^{*}})=1-\Phi\left[\
-(\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{j}^{*}}})/\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{j}^{*}}}^{2})}\
\right]$ (64)

Here, $\Phi$ represents the cumulative distribution function (CDF) of the
standard normal distribution. Since the ordering results in $Y$ is given by
the ordering of latent representation, we can use Equation (64) to calculate
matrix $P$ as follows.

Denote $\gamma_{ij}=P(\hat{h_{i}^{*}}<\hat{h_{j}^{*}})=\Phi\left[\
-(\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{j}^{*}}})/\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{j}^{*}}}^{2})}\
\right]$, we can calculate the value of each entry in $P$ (see Lemma D.2)
under the assumption of imperfect reconstruction as:

$\displaystyle p_{ij}$
$\displaystyle=\sum_{\mathcal{S}_{i},\mathcal{T}_{i}}\prod_{k\in\mathcal{S}_{i}}\Big{[}P(\hat{h_{j}^{*}}<\hat{h_{k}^{*}})\Big{]}\prod_{k\in\mathcal{T}_{i}}\Big{[}P(\hat{h_{j}^{*}}>\hat{h_{k}^{*}})\Big{]}$
(65)
$\displaystyle=\sum_{\mathcal{S}_{i},\mathcal{T}_{i}}\prod_{k\in\mathcal{S}_{i}}\Big{(}\gamma_{jk}\Big{)}\prod_{k\in\mathcal{T}_{i}}\Big{(}1-\gamma_{jk}\Big{)}$
(66)

where the summation is over all possible partitions of the index set
$\\{1,2,...,M\\}\setminus{i}$ into two disjoint subsets $\mathcal{S}_{i}$ and
$\mathcal{T}_{i}$ with $|\mathcal{S}_{i}|=i-1$ and $|\mathcal{T}_{i}|=M-i$.

An Intuitive Interpretation of Equation (66): To have the $j$-th element
$\bm{x_{j}}^{*}$ from the target sequence $Y^{*}$ to be the $i$-th element
$\bm{x_{i}}$ in the resulting sequence $Y$, the encoded 1-D representation of
$\bm{x_{j}}^{*}$ has to be smaller than the 1-D latent of $M-i$ number of
points and it has to be larger than the 1-D latent of $i-1$ number of points.
We need to consider all partitions to sum up the probabilities from all
possible orderings that has the $j$-th element from $Y^{*}$ to be placed at
$i$-th location in $Y$.

Note that Equation (66) assumes no ordering ambiguity is presented. If we want
to get a more accurate estimation of the error for latent sort algorithm, we
can add ordering ambiguity in it. The composite $P$ matrix can be estimated by
simply calculating a weighted summation over the two $P$ matrices from
Equation (14) and Equation (66), then normalize the summed matrix to be a
probability matrix.

###  C.6 Simplifying The Errors From Imperfect Reconstruction

Equation (66) gives an estimation of the errors of latent sort algorithm using
the matrix $P$. However, without simplification, it is intractable to
calculate the probability values due to the factoral number of possible
partitions $\mathcal{S}_{i}$ and $\mathcal{T}_{i}$. Here we will first provide
a way to simplify Equation (66), then we will use it to reveal the shortest
path property of latent sort algorithm.

The simplification leverages the sparsity in the summation in Equation (66),
based on the property that if any term in the product is zero, the product
becomes zero.

Let
$\Gamma_{j}=\prod_{k\in\mathcal{S}_{i}}\left(\gamma_{jk}\right)\prod_{k\in\mathcal{T}_{i}}\left(1-\gamma_{jk}\right)$.
In a partitions that yields $\mathcal{S}_{i}$ and $\mathcal{T}_{i}$ among all
possible partitions, if $\exists k\in\mathcal{S}_{i},\gamma_{jk}\to 0$ or
$\exists k\in\mathcal{T}_{i},\gamma_{jk}\to 1$, the value of $\Gamma_{j}\to
0$.

From Theorem C.3, we can see that the minimization on reconstruction loss,
empirically measured by $\|\bm{\sigma_{\hat{\bm{x_{i}}}}^{2}}\|$, can lead to
the collapse of the distribution of $\hat{h_{i}^{*}}$ to a sinlge value
$h_{i}$ (see Remark E.6). Thus, when we train the auto-encoder with
reconstruction loss, we can assumed the
$\|\bm{\sigma_{\hat{\bm{x_{i}}}}^{2}}\|$ is properly minimized (meaning it is
bounded by a small constant), and the value of $\mu_{\hat{h_{i}^{*}}}$ will
not deviate from $h_{i}$ too much, such that we can assume
$\mu_{\hat{h_{i}^{*}}}$ follows ascending order as $h_{i}$ does, i.e.,
$\mu_{\hat{h_{i}^{*}}}\leq\mu_{\hat{h_{j}^{*}}}\iff i\leq j$.

In addition, from Theorem C.2 and C.1, we know that neighboring points in the
original space $\mathbb{R}^{N}$ tend to have similar 1-D latent values.
Moreover, from the property of $\Phi$ as the CDF of normal distribution, the
value of $\Phi(x)$ quickly becomes infinitesimal when $x$ starts to deviate
from zero towards negative infinite, similarly $1-\Phi(x)$ quickly becomes
infinitesimal when $x$ starts to deviate from zero towards positive infinite.
As $\|\bm{\sigma_{\hat{\bm{x_{i}}}}^{2}}\|$ getting smaller in training, the
value of $\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{j}^{*}}}^{2})}$
is also getting smaller, thus the value of $\gamma_{ij}$ will be more
sensitive to the difference of mean values
$\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{j}^{*}}}$, pushing more and more values to
zero and one.

Therefore, we can replace the value of $y_{ij}$ to be 0 or 1 when $i$ and $j$
are not neighbors, and we can approximate $\gamma_{ij}$ as follows:

$\gamma_{ij}\approx\begin{cases}1,&\text{if}\ i<j-1\\\ \Phi\left[\
-(\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{j}^{*}}})/\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{j}^{*}}}^{2})}\
\right],&\text{if}\ j-1\leq i\leq j+1\\\ 0,&\text{if}\ i>j+1\end{cases}$ (67)

Thus using Equation (67) and the definition of $\mathcal{S}_{i}$ and
$\mathcal{T}_{i}$, we can simplify Equation (66) as follows:

$\displaystyle p_{ij}$
$\displaystyle=\sum_{\mathcal{S}_{i},\mathcal{T}_{i}}\prod_{k\in\mathcal{S}_{i}}\Big{(}\gamma_{jk}\Big{)}\prod_{k\in\mathcal{T}_{i}}\Big{(}1-\gamma_{jk}\Big{)}$
(68)
$\displaystyle=\begin{cases}\big{[}\gamma_{j,(j-1)}\big{]}\big{[}1-\gamma_{j,(j+1)}\big{]}+\big{[}\gamma_{j,(j+1)}\big{]}\big{[}1-\gamma_{j,(j-1)}\big{]},&i=j\\\
\big{[}\gamma_{j,(j-1)}\big{]}\big{[}\gamma_{j,(j+1)}\big{]},&i=j-1\\\
\big{[}1-\gamma_{j,(j-1)}\big{]}\big{[}1-\gamma_{j,(j+1)}\big{]},&i=j+1\\\
0,&i\notin[j-1,j+1]\end{cases}$ (69)

Equation (69) gives a simplification of $P$ matrix in the case of imperfect
reconstruction, where $P$ becomes a tridiagonal matrix with non-zero values on
the main diagonal and the diagonals immediately above and below it, which
looks like:

$P=\begin{pmatrix}p_{11}&p_{12}&0&\cdots&0\\\ p_{21}&p_{22}&p_{23}&\cdots&0\\\
0&p_{32}&p_{33}&\cdots&0\\\ \vdots&\vdots&\vdots&\ddots&\vdots\\\
0&0&0&\cdots&p_{MM}\\\ \end{pmatrix}$ (70)

We can do a quick check to verify that $P$ is still a probability matrix in
the simplied form (see Remark E.7).

Now we can use Equation (69) and Equation (70) to estimate the sorting error
$\mathcal{E}$. Lemma D.1 states the equivalence in minimization between
$\mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}$ and $\|P-I_{M}\|_{F}^{2}$, thus
we can calculate the difference between $P$ and $I_{M}$ using the simplified
equation by summing the main diagonal difference and the off-diagonal
difference:

$\displaystyle\|P-I_{M}\|_{F}^{2}$
$\displaystyle=\sum_{i=1}^{M}(p_{ii}-1)^{2}+\sum_{i=1}^{M-1}\Big{[}p_{i(i+1)}^{2}+p_{(i+1)i}^{2}\Big{]}$
(71)
$\displaystyle=\sum_{i=1}^{M}(p_{ii}-1)^{2}+\sum_{i=2}^{M}p_{(i-1)i}^{2}+\sum_{i=1}^{M-1}p_{(i+1)i}^{2}$
(72)

Where the main diagonal:

$\displaystyle(p_{ii}-1)^{2}$
$\displaystyle=\left\\{\big{[}\gamma_{i,(i-1)}\big{]}\big{[}1-\gamma_{i,(i+1)}\big{]}+\big{[}\gamma_{i,(i+1)}\big{]}\big{[}1-\gamma_{i,(i-1)}\big{]}-1\right\\}^{2}$
(73)

And the diagonals immediately above and below the main diagonal:

$\displaystyle p_{(i-1)i}^{2}$
$\displaystyle=\big{[}\gamma_{i,(i-1)}\big{]}^{2}\big{[}\gamma_{i,(i+1)}\big{]}^{2}$
(74) $\displaystyle p_{(i+1)i}^{2}$
$\displaystyle=\big{[}1-\gamma_{i,(i-1)}\big{]}^{2}\big{[}1-\gamma_{i,(i+1)}\big{]}^{2}$
(75)

Based on the above simplfication, now we can state the following theorem.

######  Theorem C.4.

Define $P$ using Equation (69). The minima of the error term
$\|P-I_{M}\|_{F}^{2}$ occur at $\gamma_{i,(i-1)}=0$ and $\gamma_{i,(i+1)}=1$
for all $i\in\\{2,3,...,M-1\\}$.

###### Proof.

Let $x=\gamma_{i,(i-1)}$ and $y=\gamma_{i,(i+1)}$. To find the minima of the
error term $\|P-I_{M}\|_{F}^{2}$, we need to differentiate the error term with
respect to $x$ and $y$, and find the values that minimize the error term. Here
is the error term we are working with:

$\displaystyle\|P-I_{M}\|_{F}^{2}=\sum_{i=1}^{M}(p_{ii}-1)^{2}+\sum_{i=2}^{M}p_{(i-1)i}^{2}+\sum_{i=1}^{M-1}p_{(i+1)i}^{2}$
(76)

Using Equation (70) and exclude the cases that makes $p_{ij}=0$, we have:

$\displaystyle\begin{cases}p_{ii}=x(1-y)+y(1-x),\\\ p_{(i-1)i}=xy,\\\
p_{(i+1)i}=(1-x)(1-y).\end{cases}$ (77)

Now, let’s differentiate the error term with respect to $x$ and $y$. First, we
will differentiate the error term with respect to $x$:

$\displaystyle\frac{\partial\|P-I_{M}\|_{F}^{2}}{\partial x}$
$\displaystyle=\sum_{i=1}^{M}2(p_{ii}-1)\frac{\partial p_{ii}}{\partial
x}+\sum_{i=2}^{M}\Big{[}2p_{(i-1)i}\frac{\partial p_{(i-1)i}}{\partial
x}\Big{]}$ (78) $\displaystyle\ \ \ \
+\sum_{i=1}^{M-1}\Big{[}2p_{(i+1)i}\frac{\partial p_{(i+1)i}}{\partial
x}\Big{]}$ (79)

Now, we differentiate $p_{ii}$, $p_{i(i+1)}$, and $p_{(i+1)i}$ with respect to
$x$:

$\displaystyle\begin{cases}\frac{\partial p_{ii}}{\partial x}=(1-y)-y=1-2y,\\\
\frac{\partial p_{(i-1)i}}{\partial x}=y,\\\ \frac{\partial
p_{(i+1)i}}{\partial x}=y-1.\end{cases}$ (80)

Plugging these derivatives back into the expression for the partial derivative
with respect to $x$, we get:

$\displaystyle\frac{\partial\|P-I_{M}\|_{F}^{2}}{\partial x}$
$\displaystyle=\sum_{i=1}^{M}\Big{[}2(p_{ii}-1)(1-2y)\Big{]}+\sum_{i=2}^{M}\Big{[}2p_{(i-1)i}(y)\Big{]}$
(81) $\displaystyle\ \ \ \ +\sum_{i=1}^{M-1}\Big{[}2p_{(i+1)i}(y-1)\Big{]}$
(82) $\displaystyle=\sum_{i=1}^{M}\Big{[}2(x+y-2xy-1)(1-2y)\Big{]}$ (83)
$\displaystyle\ \ \ \
+\sum_{i=2}^{M}\Big{[}2xy^{2}\Big{]}+\sum_{i=1}^{M-1}\Big{[}2(1-x-y+xy)(y-1)\Big{]}$
(84)

Since $x$ and $y$ are between $[0,1]$, for the above partial derivative to be
zero, either $x=0,y=1$ or $x=1,y=0$.

Next, we will differentiate the error term with respect to $y$. Notice that
symmetry in Equation (77), thus similarly we have:

$\displaystyle\frac{\partial\|P-I_{M}\|_{F}^{2}}{\partial y}$
$\displaystyle=\sum_{i=1}^{M}\Big{[}2(x+y-2xy-1)(1-2x)\Big{]}$ (85)
$\displaystyle\ \ \ \
+\sum_{i=2}^{M}\Big{[}2yx^{2}\Big{]}+\sum_{i=1}^{M-1}\Big{[}2(1-x-y+xy)(x-1)\Big{]}$
(86)

Again, to make the above derivative zero, we have $x=0,y=1$ and $x=1,y=0$ are
two possible solutions. Thus, by setting the two partial derivatives to zero,
we can see that the two minima of $\|P-I_{M}\|_{F}^{2}$ are $x=0,y=1$ and
$x=1,y=0$.

However, recall the definition $x=\gamma_{i,(i-1)}$ and $y=\gamma_{i,(i+1)}$,
and the definition of $\gamma_{ij}=P(\hat{h_{i}^{*}}<\hat{h_{j}^{*}})$. Since
Equation (67)(69)(70) are defined over every row of $P$, which correspond to
every position in the sorted sequence, we can see that the two minima
represent two different orderings in the sorted sequence, where $x=0,y=1$
gives ascending order and $x=1,y=0$ gives descending order. Since in our
discussion we only consider ascending order in the target sequence, $x=0,y=1$
is the only minimum that satisfy our definitions. ∎

Theorem C.4 implies that in order to minimize $\mathcal{E}$, we can minimize
$\|P-I_{M}\|_{F}^{2}$ instead (see Lemma D.1) by minimizing an equivalent
objective $\gamma_{i,(i-1)}-\gamma_{i,(i+1)}$ for all $i=2,...,M-1$, where
both $\gamma_{i,(i-1)}$ and $\gamma_{i,(i+1)}$ are between 0 and 1. Use the
definition of $\gamma_{ij}$, rewrite this alternative minimization problem on
$\forall i\in\\{2,3,...,M-1\\}$ as:

$\displaystyle\min_{f}$ $\displaystyle\ {\gamma_{i,(i-1)}-\gamma_{i,(i+1)}}$
(87) $\displaystyle=\min_{f}$ $\displaystyle\ {\Phi\left[\
-\frac{\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{i-1}^{*}}}}{\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{i-1}^{*}}}^{2})}}\
\right]}-{\Phi\left[\
-\frac{\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{i+1}^{*}}}}{\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{i+1}^{*}}}^{2})}}\
\right]}$ (88) $\displaystyle=\min_{f}$ $\displaystyle\ {\Phi\left[\
\frac{\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{i+1}^{*}}}}{\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{i+1}^{*}}}^{2})}}\
\right]}-{\Phi\left[\
\frac{\mu_{\hat{h_{i}^{*}}}-\mu_{\hat{h_{i-1}^{*}}}}{\sqrt{(\sigma_{\hat{h_{i}^{*}}}^{2}+\sigma_{\hat{h_{i-1}^{*}}}^{2})}}\
\right]}$ (89)

Assuming the 1-D latents are normalized to [0, 1]. From the property of CDF
$\Phi$, we know that (89) reaches minima -1 only when the reconstruction is
perfect such that the variances
$\sigma_{\hat{h_{i-1}^{*}}},\sigma_{\hat{h_{i}^{*}}},\sigma_{\hat{h_{i+1}^{*}}}$
are zero.

###  C.7 Optimal Sorting From The Perspective of Shortest Path Problem

In previous theories, we have discussed the error of $Y$ to the target sorting
$Y^{*}$ using the probability matrix $P$. However, the target sorting is not
yet defined. The optimal sorting need to have properties to minimize the
errors induced by imperfect $P$ matrices. Let us first discuss the desired
properties of $Y^{*}$.

Recall that $\mathbb{E}[Y]=PY^{*}$. From Equation (1) and (70), we have

$\displaystyle\ \ \ \ \mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}$ (90)
$\displaystyle=\|PY^{*}-Y^{*}\|_{F}^{2}$ (91)
$\displaystyle=\sum_{i=1}^{M}\left\|\sum_{j=1}^{M}p_{ij}\bm{x_{j}}^{*}-\bm{x_{i}}^{*}\right\|_{F}^{2}$
(92)
$\displaystyle=\sum_{i=1}^{M}\left\|p_{i(i-1))}\bm{x_{i-1}}^{*}+(p_{ii}-1)\bm{x_{i}}^{*}+p_{i(i+1))}\bm{x_{i+1}}^{*}\right\|_{F}^{2}$
(93)
$\displaystyle=\sum_{i=1}^{M}\Big{\|}p_{i(i-1))}\bm{x_{i-1}}^{*}+\big{[}-p_{i(i-1))}-p_{i(i+1))}\big{]}\bm{x_{i}}^{*}+p_{i(i+1))}\bm{x_{i+1}}^{*}\Big{\|}_{F}^{2}$
(94)
$\displaystyle=\sum_{i=1}^{M}\left\|p_{i(i-1))}(\bm{x_{i-1}}^{*}-\bm{x_{i}}^{*})+p_{i(i+1))}(\bm{x_{i+1}}^{*}-\bm{x_{i}}^{*})\right\|_{F}^{2}$
(95)

Recall that $f$ is a local DR that maps each point in $\mathcal{X}$ to its 1-D
latent representation independent from the other elements in the sequence.
That being said, we expect (95) to be minimized for any three points
$\\{\bm{x_{i-1}}^{*},\bm{x_{i}}^{*},\bm{x_{i+1}}^{*}\\}\subset\mathcal{X}$. To
make $\mathcal{E}\Big{(}\mathbb{E}[Y],Y^{*}\Big{)}$ minimized for all possible
$\mathcal{X}\subseteq\mathbb{R}^{N}$, we need to minimize the upper bound of
(95), giving by:

$\displaystyle\ \ \ \
\sum_{i=1}^{M}\left\|p_{i(i-1))}(\bm{x_{i-1}}-\bm{x_{i}})+p_{i(i+1))}(\bm{x_{i+1}}-\bm{x_{i}})\right\|_{F}^{2}$
(96)
$\displaystyle\leq\sum_{i=1}^{M}\left\|p_{i(i-1))}(\bm{x_{i-1}}^{*}-\bm{x_{i}}^{*})\right\|_{F}^{2}+\left\|p_{i(i+1))}(\bm{x_{i+1}}^{*}-\bm{x_{i}}^{*})\right\|_{F}^{2}$
(97)