1 Introduction↩︎

We study the Traveling Salesperson Problem (TSP), one of the key problems in combinatorial optimization: We are given a road network between \(n\) cities represented by an undirected edge-weighted graph \(G\). Our task is to find a closed walk in \(G\) of minimum total weight that visits every city (at least once).

Our main result is an algorithm that takes as input a TSP instance and additionally, for each edge, a number between \(0\) and \(1\) denoting how likely this edge is to belong to an optimal TSP tour. These additional numbers can be understood as predictions or a heatmap. The algorithm outputs a tour that is guaranteed to be at most \(2\eta\) longer than an optimal tour, where \(\eta\) denotes the L1 distance between the predictions vector and a groundtruth optimal tour.

Our contribution can be understood in two contexts: a theoretical context of the so-called algorithms with predictions (a.k.a. learning-augmented algorithms) [1], [2], and a practical context of the so-called solution search stage in the neural combinatorial optimization pipeline [3].

1.0.0.1 Algorithms with predictions.

This line of work tries to combine provable worst-case guarantees of classical algorithms with great beyond worst-case performance of ML models. The idea is to have a classical algorithm taking as additional input an uncertain hint – prediction of an ML model – and have its performance (e.g., running time or solution cost) bounded by a function of the error of the prediction. Despite almost 400 papers [2] developed over the last few years, approximation algorithms for TSP have not been studied in this paradigm. We address this gap in the literature.

1.0.0.2 Neural combinatorial optimization.

It is a big question if and how neural networks (NNs) can help solving combinatorial optimization problems with many hard constraints, which are usually hard to impose on NN’s output [4]. TSP has long served as a benchmark for combinatorial optimization techniques, and it has a long history of hand-optimized solvers [5], [6]. So it is interesting to understand if NNs can beat them. A common neural optimization approach to TSP – formalized by [3] into a five-stages pipeline – is to have a (graph) NN that assigns to each edge the probability that it belongs to an optimal solution (the so-called heatmap), followed by a “solution search” algorithm that turns these probabilities into an actual solution (a tour). Initially, researchers used very simple solution search algorithms, e.g., greedy search or beam search [7]–[10], which have no theoretical guarantees linking the NN performance to the final solution cost. Our algorithm is a drop-in replacement for the solution search stage, and tends to work better than greedy in practice (see Section 8).

More recently, researchers started using Monte Carlo Tree Search (MCTS) for solution search [11], making it the heaviest part of the pipeline, heavier than the NN itself. Actually, [12] argue that then NNs become useless, because one can as well replace the input to MCTS with a simple heuristic (the softmax of negative edge weight). Overall, MCTS does lead to better results, but this does not answer the original question of how NNs can help in combinatorial optimization.

To summarize, people started using MCTS in the pipeline because greedy search was not good enough, but with MCTS the NN becomes useless. Our algorithm is a better replacement for greedy search, but, unlike MCTS, it heavily uses the NN output while itself remains lightweight. Hence, our algorithm creates an opportunity for better NNs to become useful overall.

1.1 Our results↩︎

We propose algorithms for TSP that, given any vector of probabilities (a.k.a. heatmap) of each edge to belong to an optimal solution, computes a feasible TSP tour. This tour is provably near-optimal whenever the given heatmap is sufficiently accurate. We experimentally compare the solutions produced by our algorithms on Euclidean instances of TSP with existing heatmap-based techniques in the literature, considering heatmaps produced by NNs as well as by standard heuristics.

1.1.0.1 Definition of TSP.

Let \(G=(V,E)\) be a complete undirected graph (i.e., \(E=\binom{V}{2}\)) and \(w\colon E\to \mathbb{R}_+\) the weights of its edges satisfying the triangle inequality¹ (\(w(uv)\leq w(uz) + w(zv)\) for all \(u,v,z\in V\)). We want to find the optimal TSP tour, i.e., the Hamiltonian circuit on \(G\) with the smallest weight.

Consider a heatmap \(p\in [0,1]^E\). Ideally, we want \(p\) to be the characteristic vector \(\chi^{X^*}\) of some optimal TSP tour \(X^*\). We define the prediction error of \(p\) with respect to \(X^*\) as \[\eta = \eta(p,X^*)= \sum_{e\in E} \big|p_e - \chi^{X^*}(e)\big|\, w(e).\] Our algorithm satisfies the following:

In fact, our algorithm works with a prediction \(P\subseteq E\), where each \(e\in E\) is included in \(P\) independently with probability \(p\). Theorem 1 can be equivalently stated in terms of \(\eta(P,X^*) := \eta(\chi^P,X^*) = w(P\Delta X^*)\), where \(\mathbb{E}[\eta(P,X^*)] = \eta(p,X^*)\). In Section 7, we show that the linear dependence on \(\eta\) is necessary.

The statement of Theorem 1 holds for all TSP tours \(X^*\) simultaneously (optimal or not), which can be interpreted as follows: As far as \(p\) has a small error with respect to some TSP tour of a small weight, the algorithm will also find a tour of a small weight. An alternative way of stating this guarantee is \[\mathbb{E}[w(X)] \leq \min \{w(X^*) + 2\eta(p,X^*)\mid X^* \text{ is a tour in G}\}.\]

Our algorithm is based on the classical 1.5-approximation algorithm by [14]: It starts by finding a minimum spanning tree \(T\) biased towards edges suggested by \(p\). Then, it finds a minimum-weight perfect matching \(J\) of the odd-degree vertices in \(T\). While efficient implementations of algorithms for minimum-weight perfect matching exist [15], this step remains the major contribution to the running time of our algorithm. We can replace the exact algorithm for matching with a 2-approximation [16] to achieve a near-linear running time while losing a constant factor only in the dependence on \(\eta\). However, this requires, counterintuitively, biasing \(T\) towards the heaviest edges suggested by \(p\) and a more careful analysis. Note that the input size of a weighted complete graph on \(n\) vertices is \(\Theta(n^2)\). Even faster running times can be achieved on Euclidean and graphical instances which allow succinct representation, see below.

1.1.0.2 Euclidean TSP.

This is a special case of TSP which admits a succinct representation of the input complete graph \(G=(V,E)\). Each vertex \(v\in V\) is a point in the Euclidean plane and from such input of size \(\Theta(n)\) we can compute \(w(uv) = \lVert u-v\rVert_2\) for each pair \(u,v\in V\). Theorem 1 holds in this setting as stated and we can improve the running time of the algorithm in Theorem 2 to \(O((n+|P|)\mathop{\mathrm{polylog}}n)\) using the fast algorithm for Euclidean MST [17] and an approximation algorithm for Euclidean matching [18].

1.1.0.3 Graphical TSP.

Here, \(G\) is a connected graph and \(w\colon E\to \mathbb{R}^+\) is arbitrary [19]. In this case, an optimal tour may need to traverse some edges twice and visit some vertices multiple times. Every instance of metric TSP can be seen as an instance of graphical TSP on the same complete graph, since repetitions of vertices can be eliminated without increasing the weight of the tour. On the other hand, any instance of a graphical TSP on a graph \(G\) can be modeled as an instance of metric TSP on the metric closure of \(G\) which is a complete graph whose edge weights represent the shortest path metric on \(G\). We refer to [20] for more details. The most interesting case is when \(G\) is sparse and can be represented in the space much smaller than its metric closure. Therefore, we want an algorithm with near-linear running time in the size of \(G\) instead of its metric closure.

Since tours in this setting may pass some edges several times, we assume that the predictor provides us information about the multiplicities. The algorithm receives a multiset \(P\) of edges which are likely to be used by an (optimal) tour \(X^*\) (that is also a multiset) which may be derived, e.g., from the output of some NN. Information about the multiplicities of the edges is necessary to achieve near-linear running time, see Section 6 for details. The error \(\eta\) of prediction \(P\) with respect to \(X^*\) is defined as \[\label{eq:grapherror} \eta := w(P \Delta X^*) = \sum_{e\in E} \lvert m_P(e) - m_{X^*}(e)\rvert\, w(e),\qquad{(1)}\] where \(m_S(e)\) denotes the multiplicity of \(e\) in the multiset \(S\).

Interestingly, if we do not require near-linear running time, predicting edge multiplicities is not necessary. Algorithm from Theorem 1 can be extended to this setting, preserving both its running time and approximation guarantees.

2 Related work↩︎

2.0.0.1 Algorithms with predictions.

Also known as learning-augmented algorithms (see the survey by [1] and updated list of papers by [2]), have been first studied for online problems [21], [22], with the premise that ML models may predict future data. Starting with the seminal work by [23], predictions are also used to improve offline (static) algorithms, which becomes useful when solving multiple similar instances and an ML model can discover a common structure. Among those, there are only few works that concern approximation algorithms. The closest to our work is the one by [24] about subset selection problems. Unfortunately, their algorithm does not preserve triangle inequality of the input instance and therefore it does not apply to TSP (see also Remark 4). [25] proposed an algorithm for permutation problems including TSP running in (large) polynomial time which finds an optimal solution in the setting with \(\epsilon\)-accurate predictions where errors occur independently at random.

2.0.0.2 Neural combinatorial optimization.

In this vast area of research [4], [26], TSP plays a central role. Papers most closely related to our work are discussed in Sections 1 and 8.

Beyond the heatmap-to-tour paradigm, several alternative approaches have emerged. Reinforcement learning (RL) has been explored extensively, with methods like POMO [27] and other policy-based approaches [28] showing promising results. Transformer-based architectures have also gained attention, including MatNet [29] and BQ-NCO [30]. Additional architectural innovations are explored in Sym-NCO [31] and DIMES [32].

More recently, Preference Optimization (a technique that transforms quantitative reward signals into qualitative preference comparisons) has gained attention as an alternative to RL-based approaches [33]–[35]. Concurrently, other works have explored hybrid architectures that blend global predictions with local solution construction [36], [37], while further methods divide the search space into sub-problems solved after by neural solvers in sequence [38].

Another emerging paradigm aims to train general-purpose neural solvers capable of addressing diverse combinatorial problems with a single model [39], [40]. For a comprehensive analysis of these various architectural choices and their roles in neural combinatorial optimization, particularly for TSP, we refer the reader to [41].

2.0.0.3 TSP.

The problem is known to be NP-hard, and even hard to approximate if distances do not satisfy the triangle inequality and nodes cannot be visited twice [42]. In the metric case, there are simple 2-approximation algorithms [43], and [14] and [44] proposed a 1.5-approximation algorithm, which was only improved recently to \(1.5-10^{-36}\) [45]. We refer to the textbooks by [46] and [20].

2.0.0.4 Fast algorithms for Euclidean graphs.

Euclidean setting often allows algorithms with improved performance. In particular, [47] proposed an approximation scheme for Euclidean TSP, Euclidean minimum spanning tree can be computed in time \(O(n\log n)\) using Delaunay triangulations [17], and [18] proposed an \((1+\epsilon)\)-algorithm for Euclidean min-weight perfect matching running in time \(O(\epsilon^{-3}n\log^6 n)\).

3 Notation and preliminaries↩︎

In some parts of the paper, we work with multisets, extending classical set operations as follows. Let \(F\) be a multiset from a universe \(U\), and \(m_F(e)\in \mathbb{N}\) denote the multiplicity of \(e\in U\) in \(F\). Then the following holds for any \(3\) multisets \(S\), \(A\), and \(B\) of \(U\) and any \(e\in U\). If \(S=A\cap B\), then \(m_S(e) = \min\{m_A(e),m_B(e)\}\). If \(S = A + B\), then \(m_S(e) = m_A(e)+m_B(e)\). If \(S = A - B\), then \(m_S(e) = \max\{m_A(e)-m_B(e), 0\}\). We also define the symmetric difference \(A\Delta B = (A-B)+(B-A)\).

Given an input graph \(G=(V,E)\), we denote by \(n=|V|\) the number of its vertices and by \(m=|E|\) the number of its edges. We say that edge weights \(w\colon E\to \mathbb{R}^+\) are metric if they satisfy triangle inequality, i.e., if we have \(w(uv)\leq w(uz) + w(zv)\) for each \(u,v,z\in V\). For a set \(F\subseteq E\), we denote by \(\mathop{\mathrm{odd}}(F)\) the set of vertices in \((V,F)\) with an odd degree. We say that \(G = (V,E)\) is Eulerian, if \(\mathop{\mathrm{odd}}(E) = \emptyset\). Given a set \(S\subseteq V\) of even size, we call \(J \subseteq E\) an \(S\)-join if \(\mathop{\mathrm{odd}}(J) = S\). We call \(J\subseteq E\) a perfect matching of \(S\) if, in \((V,J)\), the degree of each \(v\in S\) is equal to 1 and the degree of each \(v\in V\setminus S\) is equal to 0. Note that a perfect matching of \(S\) is a special case of an \(S\)-join. In graphs with metric weights, we can assume that all \(S\)-joins of minimum weight are perfect matchings.

A closed walk in \(G\) passing through each vertex at least once is called a tour, let \(X\) denote the multiset containing each edge \(e\in E\) as many times as it is traversed by the tour. It is well-known and easy to prove that \((V,X)\) is Eulerian. If \(G\) is a complete graph with metric edge weights and \(F\) is a multiset of edges such that \((V,F)\) is Eulerian, then there is a Hamiltonian cycle \(C\subseteq F\) of weight \(w(C)\leq w(F)\) which can be found using the standard short-cutting procedure [20].

4 Algorithms for metric TSP↩︎

Our algorithms receive in input a prediction \(P\subseteq E\), which is ideally the set of edges of an optimal tour \(X^*\). Notice that in general we allow \(P\) to be any set of edges, in particular \(P\) might induce a disconnected graph, a graph with multiple cycles etc. In the settings where we are given a heatmap \(p\in [0,1]^E\) in input, we simply sample each edge \(e\in E\) independently with probability \(p_e\), and let \(P\) be the set of sampled edges. We remark that \(\mathbb{E}[\eta(P,X^*)]=\eta(p,X^*)\), therefore it is sufficient to upper bound the cost of the produced solution in terms of \(\eta=\eta(P,X^*)\).

4.1 Algorithm [alg:alg1] and proof of Theorem 1↩︎

We describe our algorithm proving Theorem 1. Our algorithm is based on the classical 1.5-approximation algorithm by [14]. First, it constructs a minimum spanning tree (MST) \(T\) on \(G\) which is biased towards edges in the prediction \(P\) in the sense that each edge in \(P\) is treated as if its weight was 0. This can be easily done in near-linear time. Having the spanning tree \(T\), it computes the minimum-weight \(\mathop{\mathrm{odd}}(T)\)-join \(J\) which, in the metric case, is the minimum weight perfect matching of \(\mathop{\mathrm{odd}}(T)\). We compute \(J\) using the classical algorithm [48], [49] with running time² \(O(n^3)\) treating predicted and not predicted edges equally, i.e., we always count with their original weights \(w\).

The multiset \(T+J\) where edges contained in both \(T\) and \(J\) are taken twice is connected (because \(T\) is a spanning tree) and Eulerian (each \(v\in V\) is contained in an even number of edges in \(T+J\)), since \(J\) is an \(\mathop{\mathrm{odd}}(T)\)-join, i.e., \(\mathop{\mathrm{odd}}(J) = \mathop{\mathrm{odd}}(T)\). Using the classical short-cutting procedure [20], we can construct a Hamiltonian circuit of weight no greater than the weight of \(T+J\) in time \(O(|T+J|)\). See Algorithm [alg:alg1] for the summary.

Theorem 1 is implied by the following statement and the fact that each step of Algorithm [alg:alg1] can be performed in near-linear time except for the construction of \(J\) which requires time \(O(n^3)\).

Figure 1 describes the steps of Algorithm [alg:alg1] on a simple example with a perfect prediction, showing that ideally we have \(T \subseteq X^*\), i.e., \(T = X^*\cap T\), and \(J = X^*\setminus T\). In our analysis, we show that prediction errors can make both \(T\) and \(J\) deviate from \(X^*\cap T\) and \(X^*\setminus T\) respectively, but the increase of their weight is proportional to the weight of the mispredicted edges as formulated in the lemmas below. We use the following notation. Given a prediction \(P\subseteq E\) and a tour \(X^*\), we denote \(H^+ = P \setminus X^*\) the set of false positives and \(H^- = X^*\setminus P\) the set of false negatives. We have \(\eta = w(H^-) + w(H^+)\) and \(P = (X^* \setminus H^-) \cup H^+\).

4.2 Near-linear time algorithm for metric TSP↩︎

Finding the min-weight perfect matching \(J\) is the most computationally intensive step in Algorithm [alg:alg1]. An obvious way to speed-up the algorithm is to replace this step with a 2-approximation algorithm by [16], which runs in near-linear time. Unfortunately, this does not lead to approximation ratio smaller than 1.5 even with very small \(\eta\). An example similar to Remark 5 can be found in Appendix 11.2.3. Here, we focus on the points in the analysis of the previous section which would fail, in order to motivate the steps performed by our next algorithm.

Consider the example in Figure 1. We have \(n\) points on the line and a perfect prediction \(P\) of the optimal tour \(X^*\) and a spanning tree \(T\) that omits the heaviest edge \(e\in X^*\) such that \(w(e) = w(X^*)/2\). Then, losing a factor of \(2\) in the bound of Lemma 2 due to the approximation ratio of [16], we would have \(w(T)+w(J)\leq w(X^*\cap T) + 2w(X^*\setminus T) = \frac{3}{2} w(X^*)\), i.e., an 1.5-approximation even with perfect predictions.

In order to avoid the previous example, we bias \(T\) towards the heaviest edges in \(P\). In particular, we look for the MST with the signs of the weights of the edges in \(P\) flipped. Intuitively, we are trying to make sure that the heaviest edges in \(P\) are contained in \(T\) and we pay for them only once, while the weight of \(X^*\setminus T\) where we suffer the multiplicative factor of 2 will be small.

Unfortunately, \(X^*\setminus T\) is always non-empty and contains a part of the optimal solution \(X^*\) whose weight we do not want to pay twice. Somewhat surprisingly, we can show that whenever \(|S|>2\) (\(\eta>0\) in such a case), \(w(J)\) can be bounded purely in terms of \(\eta\), see the lemma below. This means that the approximation ratio of [16] will affect only the constant in front of \(\eta\). Fortunately, this is already enough, since finding the exact min-weight perfect matching is trivial with \(|S|=2\) (it is just the edge connecting the two vertices in \(S\)). Our algorithm is summarized in Algorithm [alg:alg2]. All the steps can be performed in near-linear time.

Lemma 3 is the crucial part of our analysis. We bound the weight of \(J\) using an explicit \(\mathop{\mathrm{odd}}(T)\)-join \(J' = X^*\Delta T\). The proof is most difficult when \(\eta\) is relatively small and the most tricky part is in fact accounting for the weight of the lightest edge. We define \(e^*\) as the edge in \(X^*\) with the smallest weight. First, we ignore \(e^*\) and consider \(Y^*=X^*\setminus \{e^*\}\). We partition \(J'\) into parts as in Figure 2 (a), bounding the weight of each part separately, as shown in the figure.

The flipped sings of the weights of the edges in \(P\) are used in the bound \(w((Y^*\setminus T)\cap P)\leq w(T\setminus Y^*)\). However, this bound is not yet sufficient to prove Lemma 3. We need to bound \(w((X^*\setminus T)\cap P)\) instead, and there might be \(e^*\) missing on the left- or right-hand side. We need to consider several cases depending on whether \(e^*\) belongs to \(P\) or not and, if yes, whether it belongs to \(T\).

4.2.1 Analysis of Algorithm [alg:alg2]↩︎

Let \(X^* \subseteq E\) be the edges of an optimal tour of \(G\). Since \(G\) is a complete graph with metric weights, \(X^*\) is a Hamiltonian cycle. Given prediction \(P\subseteq E\), we denote \(P = (X^* \setminus H^-) \cup H^+\), where \(H^-\) and \(H^+\) are the sets of false negatives and false positives. We denote by \(J'= X^*\Delta T\) the reference \(S\)-join, with \(S=\mathop{\mathrm{odd}}(T)\), which we use to bound the weight of \(J\).

The following lemma holds for an arbitrary \(S\) but we will need it only for \(|S|=2\).

Consider next the case \(|S|>2\). In this case, \(J\) is only a 2-approximation of the optimal \(S\)-join and therefore we can claim only \(w(J)\leq 2w(J')\). Therefore, to prevent loosing factor of 2 on the weight of the edges in \(X^*\), we have to bound \(w(J')\) only in terms of \(H^+\) and \(H^-\). This is proved in the next lemma. If \(|S| > 2\), we have \[w(J') \leq 3w(H^+) + 5w(H^-).\] Note that, when \(|S|>2\), we have \(T \setminus X^* \neq \emptyset\) and at least one of \(H^+\) and \(H^-\) is non-empty. We use the following notation. Fix \(e^* \in X^*\) such that \(w(e^*)\) is minimal over the edges in \(X^*\) and denote \(Y^* = X^* \setminus \{e^*\}\). We will bound the weight of the \(S\)-join \(J'= X^*\Delta T\) piece by piece, as visualized in Figure 2. The weight of \(e^*\) will be bounded by case analysis.

Let \(G=(V, E)\) be a complete graph on \(n\) vertices with edge weights \(w\colon E\to \mathbb{R}_+\) satisfying triangle inequality. There is an algorithm running in time \(O(n^2\log n)\) which, receiving a prediction \(P \subseteq E\) with error \(\eta\) with respect to some TSP tour \(X^*\), outputs a TSP tour \(X\) of weight \[w(X) \leq w(X^*) + 7w(H^+) + 11w(H^-).\]

5 Euclidean TSP↩︎

In Euclidean TSP, the input is given as coordinates of the \(n\) vertices of the complete graph \(G\) in the Euclidean plane and the weight of the edge between two vertices \(u,v\in V\) is computed as a Euclidean distance between them. In particular, the size of the input is linear in \(O(n + |P|)\), which is \(O(n)\) for any prediction of reasonable quality. This is much more succinct representation than the vector of edge weights which has \(\Theta(n^2)\) coordinates.

First, we argue that we can find \(T\) in time \(O((n+|P|)\log n)\). To this end we recall the notion of Delaunay triangulation, which is a classical object in Computational Geometry [51]. It can be defined as any triangulation of the Delaunay Graph which has the following characterization.

We construct a Delaunay triangulation \(D\) of \(V\) (it has size \(O(n)\) and it can be found in time \(O(n\log n)\)), then we find \(T\) as an MST of \((V, P\cup D)\). It is well known that an Euclidean MST uses only the edges of \(D\) [17]. It is not difficult to show that if we modify the weights of edges in the set \(P\), the MST will use only edges in \(P\cup D\).

Finally, [18] show that a \((1+\epsilon)\)-approximation of the min-weight perfect Euclidean matching can be found in time \(O(\epsilon^{-3} n\log^6 n)\). We find \(J\) using their algorithm with \(\epsilon=2\). This way, each step of Algorithm [alg:alg2] can be implemented in time \(O((n+|P|) \mathop{\mathrm{polylog}}n)\) for Euclidean TSP.

6 Graphical TSP↩︎

In graphical TSP, the optimal TSP tour may traverse some edges twice unlike in metric TSP. This means that knowing the set of edges \(F\) used by the optimal tour, there is still a non-trivial step needed to identify which edges need to be traversed twice to construct the actual tour. Consider the example in Figure 5: the underlying set \(F\) of the optimal tour contains all edges drawn in the picture regardless of their color. In order to transform \(F\) into a tour, we need to perform a costly operation: compute the optimal T-join of the odd vertices in \((V,F)\) – edges drawn both black and orange in Figure 5 which will be traversed twice in the resulting tour. Turing \(F\) into an optimal solution requires computing an optimal \(\mathop{\mathrm{odd}}(F)\)-join which is computationally demanding.

In order to achieve near-linear running time, we need the help of the prediction when searching for the T-join. Therefore we require the prediction \(P\) to be a multiset containing the information on how many times is each edge traversed and we define the prediction error as in Equation ?? .

Our algorithm for graphical TSP is summarized in Algorithm [alg:alg3]. Similarly to Algorithm [alg:alg1], it finds an MST \(T\) biased towards edges in \(P\) (the weight of each \(e\in P\) is set to 0). Then, it finds a 2-approximate \(\mathop{\mathrm{odd}}(T)\)-join \(J\) using the algorithm by [16] which runs in time \(O(m\log n)\), where \(m=|E|\) is the number of edges in the input graph, biased towards the edges in \(P-T\), i.e., the edges in \(P\) not used by \(T\) or appearing in \(P\) more than once have weight 0.

Ideally, we would have a perfect prediction \(P=X^*= T+J\) as in Figure 5 and \(P-T = X^*-T\) would be used as a prediction for finding \(J\). In our analysis, we again consider a reference \(\mathop{\mathrm{odd}}(T)\)-join \(J'=X^*\Delta T\). Our crucial points are that the weight of \(P-T\) diverges from the ideal \(X^*-T\) by at most \(\eta\), while the weight of the unpredicted edges in \((X^*-T) - (P-T)\) is bounded by \(\eta\), see Figure 6 for a summary of our bounds. Note that here we consider error \(\eta\) w.r.t. to the tour \(X^*\) not w.r.t. the optimal T-join of \(T\). The key observation is that we loose factor of 2 only for the edges in \(J' - (P-T)\) whose total weight is \(O(\eta)\), while each edge in \(P-T\) is taken at most once.

Interestingly, Algorithm [alg:alg1] can be used for graphical TSP as stated, receiving only the prediction \(P\) about the ground set of the optimal tour \(X^*\). In this case, its running time will be \(O(m\log n) + O(|S|^3)\), where the number \(|S|\) of the odd vertices in \(T\) might be as large as \(n\). Details can be found in Section 10 in the Appendix.

6.1 Analysis of Algorithm [alg:alg3]↩︎

We denote \(\bar X^*\) the multiset of edges, where each edge is contained in \(\bar X^*\) as many times as it is traversed by the optimal tour. On the other hand, we denote \(X^*\) the ground set of \(\bar X^*\) containing all edges which are traversed by the optimal tour at least once.

We will also consider an auxiliary \(\mathop{\mathrm{odd}}(T)\)-join \(J' = (T\setminus \bar X^*) + (\bar X^* - T)\). In the following observations, we bound the weight of its parts as in Figure 7.

Let \(G=(V, E)\) be a graph on \(n\) vertices and \(m\) edges with weights \(w\colon E\to \mathbb{R}_+\). Let \(X^*\) be a TSP tour on \(G\) and \(P\) a prediction with error \(\eta\) with respect to \(X^*\). Then the cost of the tour \(X\) produced by Algorithm [alg:alg3] on \(G\) with prediction \(P\) is at most \(w(X^*) + 4\eta\).

7 Lower bound: linear dependence on \(\boldsymbol{\eta}\)↩︎

Algorithm [alg:alg1] has approximation ratio \(1+2\frac{\eta}{\mathop{\mathrm{OPT}}}\) where \(\mathop{\mathrm{OPT}}\) denotes the cost of an optimal solution. We next show that a linear dependence on \(\eta/\mathop{\mathrm{OPT}}\) is needed (though possibly with a factor smaller than the \(2\) that we get). The proof idea is that otherwise one could use such a learning-augmented algorithm for TSP to develop a too good approximation algorithm for TSP (without predictions).

We remark that by a similar construction as in the proof of Theorem 16, one can argue that a substantial improvement of the dependence on \(\eta/\mathop{\mathrm{OPT}}\) in the approximation factor is very difficult. For example, a learning-augmented algorithm for TSP with approximation factor \(1+0.49\frac{\eta}{\mathop{\mathrm{OPT}}}\) would imply a \(1.49+O(\varepsilon)\) approximation for TSP. The latter result would be a breakthrough in approximation algorithms given that only recently the \(1.5\) barrier was overcome and the best known algorithm has approximation ratio \(1.5-10^{-36}\) [45].

In Appendix 11, we provide further examples showing the tightness of our analysis and counterexamples motivating the design of our algorithms.

8 Empirical results↩︎

In this section, we present empirical results that support our theoretical approach. At a high level, we consider probability matrices (heatmaps) obtained from different sources (neural predictors, non-neural predictors, and synthetic constructions) and use them as input to our algorithm. We then evaluate whether this allows us to improve upon standard baselines in terms of optimality gap. Here, we provide a high-level overview of our experimental results; full details are deferred to Appendix 12.

8.0.0.1 The dataset.

We evaluate our model on three different datasets. Dataset \(\mathcal{U}\) is the test set of the well-known ML4CO benchmark [53]. The dataset consists of instances sampled uniformly from the square \([0,1]^2\) and contains instances of four different sizes: \(1280\) instances for each of \(n \in \{50,100\}\) and \(128\) instances for each of \(n \in \{500,1000\}\). Dataset \(\mathcal{T}_E\) contains 53 Euclidean 2D TSPLIB [54] instances with up to 1300 nodes. Finally, dataset \(\mathcal{T}_M\) contains 15 metric non-Euclidean instances from TSPLIB. For all instances, we obtain ground-truth optimal tour lengths using Concorde [5].

8.0.0.2 Data and code availability.

All instances used in our experiments, together with the corresponding predictions, are available on Zenodo.³ The code for reproducing our results and implementing the main algorithm is available on GitHub.⁴

8.0.0.3 Predictions.

We evaluate several types of predictors (Details are provided in Section 12.2). Among the Predictors from NNs we use the predictions from [3] (GNN4CO), which combines a graph neural network encoder with an autoregressive decoder that sequentially builds a tour, [55] (GNN-GLS), which predicts for each edge its regret, i.e., the cost of enforcing that edge relative to an optimal solution, and [56] (DIFUSCO), a diffusion-based graph model designed to generate high-quality feasible TSP solutions. We include SoftDist [12], a non-neural probabilistic baseline where, given a weighted graph, edge probabilities are defined as \(p_{ij} = \frac{e^{{-w(ij)}/{\tau}}}{\sum_{k \neq i} e^{{-w(ik)}/{\tau}}}\), with temperature parameter \(\tau\) set with guidance from the original paper. Finally, we evaluate a Synthetic Predictor (SP), constructed by applying perturbation to the optimal solution \(X^*\). The noise is governed by the parameter \(\varepsilon\).

8.0.0.4 Baselines.

Since our approach is largely inspired by Christofides’ method, we naturally include Christofides’ algorithm (frequently abbreviated as CHR) as a baseline. We evaluate the performance of our proposed algorithm against several baseline solution search techniques commonly employed to leverage neural heatmaps [26], [55], [56]. We implement two greedy variants: G1 builds a Hamiltonian cycle incrementally by starting from an arbitrary node and repeatedly adding the unvisited node with the maximum score \(q(e)\); G2 follows a greedy edge-selection strategy, adding edges in decreasing order of \(q(e)\) provided they do not create a subcycle or a node with degree \(3\), until the tour is complete. We also evaluate a beam search (BS), which explores multiple candidate tours in parallel by maintaining the top-\(k\) partial paths at each construction step

8.0.0.5 Algorithm [alg:alg1] vs \(\text{CHR}^+\).

To derive the edge subset \(P\) for Algorithm [alg:alg1], we evaluate two strategies: sampling \(k\) edges (Alg1, faithfully representing Algorithm [alg:alg1]), selecting the \(k\) edges with the highest probabilities (Alg1Top). Furthermore, we test a modified weighting approach, where we set \(w'(e) = w(e)(1 - q(e))\), which we denote as \(\text{CHR}^+\). We set \(k=n\), using 30 samples for Alg1 to account for stochasticity. Figure 14 (Appendix 12) shows that \(\text{CHR}^+\) consistently achieves superior performance; thus, we adopt it as our primary decoding strategy. More details can be found in Appendix 12.3.

8.0.0.6 Our algorithm behaves smoothly as the prediction error decreases.

To control the prediction error, we use SP as the predictor, with \(\varepsilon \in [0, 1]\). We first observe that this type of prediction, as expected, performs very poorly with greedy decoding, and hence, comparison is not made in this case (See Appendix 12.4). To study the behavior of our algorithm relative to Christofides, we run the full pipeline for \(n = 500\) on the entire dataset. Figure 8 demonstrates that our algorithm exhibits smooth behavior as the prediction error decreases, consistently improving upon CHR until the error reaches \(\sim \varepsilon = 0.25\). Details on this experiment can be found in Appendix 12.4.

8.0.0.7 Comparison of different solution search strategies with different predictors.

Table 1 summarizes performance on dataset \(\mathcal{U}\) via optimality gaps relative to Christofides. The discrepancy with the original DIFUSCO results (Figure 9) likely stems from inference-time configurations, such as graph sparsification, which were disabled in our setup. Notably, \(\text{CHR}^+\) is the only strategy that effectively exploits predicted heatmaps to consistently outperform the Christofides baseline. In contrast, G1, G2, and BS fail to yield systematic improvements. These trends remain consistent across TSPLIB instances, both Euclidean and non-Euclidean (Table 2 and Table 3). See Appendix 12.5 for further details.

8.0.0.8 2-opt improvement.

To further refine the solutions, we apply the 2-opt local search heuristic as a post-processing step, a standard practice for enhancing tour quality [55], [56]; All the details are provided in Appendix 12.6. Our results demonstrate that \(\text{CHR}^+\) remains the superior method while introducing negligible computational overhead compared to CHR. This efficiency is evidenced by the average runtimes reported in Table 1 and the comparable number of 2-opt swaps required to reach local optimality, as recorded in Table 4.

9 Conclusions↩︎

We introduce the first learning-augmented approximation algorithm for TSP. Our algorithm converts any edge heatmap, i.e., a vector of probabilities of each edge to belong to an optimal solution, to a feasible TSP tour of weight \(\leq \mathop{\mathrm{OPT}}+ 2\eta\), where \(\eta\) denotes the L1 distance between the heatmap and some optimal tour of weight \(\mathop{\mathrm{OPT}}\). As shown in our empirical results, our approach utilizes heatmaps more effectively than previous works. Moreover, its approximation performance improves smoothly with increasing accuracy of the provided heatmap, underscoring the potential for further advancements in high-quality heatmap generation models.

Acknowledgements↩︎

The work of E. Vercesi has been supported by the Swiss National Science Foundation (SNSF) project 200021-212929 / 1 “Computational methods for integrality gaps analysis”. The work of F. Grandoni has been partially supported by the SNSF projects 200021-20073 and 200021-236706.

10 Analysis of Algorithm [alg:alg1] in the setting of graphical TSP↩︎

We analyze Algorithm [alg:alg1] in the graphical TSP setting and show that it satisfies the same guarantees as in metric TSP setting, i.e., Theorem [alg:alg1]. Compared to Algorithm [alg:alg3], it does not require predicting multiplicities of edges in the optimal tour \(\bar X^*\). Algorithm [alg:alg1] receives only a prediction \(P\) of the support set \(X^*\) of edges used by \(\bar X^*\) at least once. On the other hand, the running time of Algorithm [alg:alg1] is \(O(n^3)\) compared to the near-linear running time of Algorithm [alg:alg3].

10.0.0.1 Notation:

For an arbitrary optimal tour \(\bar X^*\), we denote by \(X^* \subseteq E\) the support of \(\bar X^*\), i.e., the set of edges appearing in \(\bar X^*\) at least once. Notice that \(X^*\) is connected and \(\bar X^*\) is Eulerian. Given a prediction \(P\subseteq E\), we write \(P = (X^* \setminus H^-) \cup H^+\), where \(H^-\) and \(H^+\) are the sets of false negatives and false positives. Notice that \(\eta=w(H^-)+w(H^+)\).

11 Tightness of our analysis and further counterexamples↩︎

Here we describe construction showing the tightness of our results and our analysis.

11.1 Tightness of our analysis↩︎

Here we show that our analysis of Algorithm [alg:alg1] is tight, in the sense that the coefficient \(2\) in front of the prediction error \(\eta\) cannot be improved for the given algorithm. We provide two examples: one with false positives and one with false negatives. These two examples can be combined in a single input instance with the desired weight of false positives and false negatives. Note that Algorithm [alg:alg1] satisfies the same approximation guarantees in both metric TSP and graphical TSP setting. In what follows, we first describe an input instance of graphical TSP which, by taking the metric closure of the proposed graph, gives an input instance of metric TSP. We describe the final outputs of Algorithm [alg:alg1] in both setting, as they will differ slightly in the final short-cutting step.

11.1.0.1 Prediction with false positives.

Figure 10 describes a bad example for the dependence on false positives. Given prediction \(P=X^*\cup\{ab\}\), the set of thick orange edges is an MST with respect to the weights \(w'\) in Algorithm [alg:alg1]. In the following steps, Algorithm [alg:alg1] identifies \(S=\{a,b,c,d\}\) as the set of odd vertices and the edges \(J\) drawn with dark-blue dashed lines denote a min-weight \(S\)-join.

After short-cutting, the resulting tour in the graphical TSP setting will take the edge \(ab\), then the path from \(b\) to \(c\) of length \(M\), then taking edge \(cd\), then path from \(d\) to \(b\) of length \(M\), and then again edge \(ab\).

In the metric TSP setting, we take the edge \(ab\), then the path from \(b\) to \(c\) of length \(M\), then taking edge \(cd\), then path from \(d\) to \(v\) and finally the edge \(va\). This forms a Hamiltonian cycle of the same weight as the graphical tour, since the weight of the edge \(va\) in the metric closure is equal to the weight of the path \(v,b,a\).

Therefore, in both cases the weight of the tour produced by Algorithm [alg:alg1] will be \[w(X) = 1 + M + 2 + M + 1 = 2M + 2 + 2 = w(X^*) + 2w(H^+) = w(X^*) + 2\eta,\] showing that the approximation guarantee for Algorithm [alg:alg1] as stated in theorems 1 and 17 is tight.

11.1.0.2 Prediction with false negatives.

Figure 11 describes a tight example for the dependence on false negatives. It shows an input where Algorithm [alg:alg1] finds a tour with weight \(w(X) \geq w(X^*) + 2\frac{k}{k+1}\eta\), where \(2\frac{k}{k+1}\) can be arbitrarily close to 2 as \(k\) can be chosen as large as roughly \(n/10\). The graph \(G\) in Figure 11 contains \(k\) blocks, each containing 5 edges of length 1 and a single edge of length 2. The optimal tour \(X^*\) uses 4 edges of unit length in each block and its total length is \(M + 4k + 3\). The prediction error is \(\eta = k+1\): one edge of unit weight is missing in each block and another such edge on the right-hand side. On the other hand, Algorithm [alg:alg1] finds an MST \(T\) (fat orange) and an \(S\)-join for \(S=\mathop{\mathrm{odd}}(T)\) (fat dashed blue), using altogether 4 edges of unit length and one edge of length 2 in each block. In total, it produces a tour of length \[w(X) = M + k({\textcolor{orange}4} + {\textcolor{blue}2}) + {\textcolor{orange}2} + {\textcolor{blue}1} = w(X^*) + 2k = w(X^*) + 2\frac{k}{k+1}\eta.\] This shows that the approximation guarantee for Algorithm [alg:alg1] as stated in theorems 1 and 17 is tight up to the factor \(\frac{k}{k+1}\) which can be arbitrarily close to 1 as \(k\) increases.

11.2 Justification of the steps performed by our algorithms↩︎

11.2.1 Doubling \(\boldsymbol{T}\) instead of \(\boldsymbol{\mathop{\mathrm{odd}}(T)}\)-join↩︎

The example in Figure 12 shows that if, in Algorithm [alg:alg1], we take \(J=T\) instead of \(J= \mathop{\mathrm{odd}}(T)\)-join which is easy to find and inspired by the classical minimum spanning tree heuristic for TSP, see [20], the resulting algorithm would not have smooth dependence on \(\eta\). In particular, \(\eta = 2\) (there are two false negative edges of weight \(1\)) and the tour produced by the algorithm would have weight \(w(X) = 2M +2M -3 = 2w(X^*) -3\) can be chosen arbitrarily small.

11.2.2 \(\boldsymbol{J}\) has to be computed with respect to the original weights↩︎

If we computed both \(T\) and \(J\) with respect to the weight \(w'\), our algorithm would be equivalent to running Christofides algorithm with the weight of predicted edges set to zero, which is an approach proposed by [24] in their meta-algorithm. This, however does not work. Consider the example in Figure 12. There, the cheapest \(\mathop{\mathrm{odd}}(T)\)-join with respect to the weights \(w'\) would be \(T\) itself with \(w'(T)=0\). At the end, we would receive a tour \(X\) with weight \(w(X) = 2w(X^*)-3\) as argued in Section 11.2.1.

11.2.3 Optimal \(\boldsymbol{\mathop{\mathrm{odd}}(T)}\)-join in Algorihtm [alg:alg1] cannot be replaced by a 2-approximation↩︎

In Figure 13, we show an input graph which shows that Algorithm [alg:alg1], if using a 2-approximate \(\mathop{\mathrm{odd}}(T)\)-join instead of the optimal one, would achieve approximation ratio close to 1.5 even with arbitrary small error. In the picture, \(T\) is an MST and has five odd vertices in the left-hand side as well as in the right-hand side. Therefore, the optimal \(\mathop{\mathrm{odd}}(T)\)-join has to contain an edge (or a path) of length \(M\) connecting one vertex in the left-hand side to another vertex in the right-hand side. The other odd vertices are then matched within their copy of \(K_5\) in an arbitrary way. The total weight of the optimal \(\mathop{\mathrm{odd}}(T)\)-join is then \(M+4\). In Figure 13, there is an \(\mathop{\mathrm{odd}}(T)\)-join \(J\) with weight \(2M + 2\), i.e., it is a 2-approximate \(\mathop{\mathrm{odd}}(T)\)-join. A tour which starts at vertex \(v\) and continues over one of the dashed-blue edges in \(J\), even after short-cutting will have cost \(3M + O(1)\), while the optimal tour has cost \(2M + O(1)\).

11.3 Non-metric TSP↩︎

It is well known that general TSP without the triangle inequality is inapproximable: unless P=NP, no polynomial-time algorithm can achieve any finite approximation ratio for non-metric TSP [13]. Unfortunately, as we show below, this hardness persists even when predictions with arbitrarily small error are available.

Specifically, for every \(r>1\) and \(\epsilon >0\), we show that there is no polynomial-time \(r\)-approximation algorithm for non-metric TSP, even when given predictions of error at most \(\epsilon\), unless P=NP. The proof is by reduction from Hamiltonian Cycle.

Without loss of generality, assume \(\epsilon < 1/r\). Let \(G=(V,E)\) be an instance of Hamiltonian Cycle on \(n=|V|\) vertices. We construct an instance of non-metric TSP on the complete graph \(K = \left(V,\binom{V}{2}\right)\) with edge weights defined by \[w_e = \begin{cases} \epsilon/n & \text{if } e\in E,\\ 1 & \text{if } e\notin E. \end{cases}\] We also provide the prediction \(\hat{X} = \emptyset\).

If \(G\) contains a Hamiltonian cycle, then the corresponding tour in \(K\) uses only edges of weight \(\epsilon/n\), and therefore has total length \(\epsilon\). In this case, the prediction error is \(\eta=\epsilon\).

On the other hand, any tour that uses at least one edge outside \(E\) has cost strictly larger than \(1\). Since \(\epsilon < 1/r\), we have \(1 > r\epsilon\). Therefore, every non-optimal tour has length greater than \(k\) times the optimum. Consequently, any \(r\)-approximation algorithm must return the optimal tour itself, which necessarily corresponds to a Hamiltonian cycle in \(G\). Hence, such an algorithm would solve Hamiltonian Cycle in polynomial time, implying P=NP.

12 More details on the experiments described in Section 8↩︎

In this section, we include additional explanations and comments on our experiments.

12.1 Implementation details↩︎

The code was implemented in Python 3.13. For graph-related tasks, such as the computation of minimum spanning trees and minimum-weight \(S\)-joins, we used the networkx [57] package. We note that more efficient implementations exist for what is currently the main bottleneck of our procedure, namely, minimum-weight perfect matching, such as Blossom V [15]. However, for this paper, we prioritized simplicity over efficiency. All experiments were performed on a cluster featuring dual Intel Xeon E5-2650 v3 CPUs with 128–512 GB of RAM. For neural network inference, we used an NVIDIA GeForce GTX 1080 GPU.

12.2 On the score computation↩︎

In this section, we provide additional details on how the scores guiding the solution search methods are computed.

DIFUSCO:

The paper by [56] provides four checkpointed models, trained on instances with \(n \in \{50, 100, 500, 1000\}\). For inference, we follow a policy of selecting the smallest checkpoint dimension that is greater than or equal to the target size for each problem size under consideration. When the instances are larger than \(1000\) nodes, we used the \(1000\) checkpoint.

GNN-GLS:

The paper by [55] proposed three pretrained models, on data of different sizes: \(n_{\text{model}} \in \{20, 50, 100\}\). We adopt the following policy: \[n_{\text{model}} = \begin{cases} 20 & n \leq 20 \\ 50 & 20 < n \leq 50 \\ 100 & \text{otherwise.} \end{cases}\] The paper proposes to predict regrets, rather than probabilities, where a regret quantifies how undesirable it is to include a given edge in the candidate solution. We then convert regrets into probabilities through a nonincreasing transformation that incorporates a scaling component, rather than relying solely on normalization. The intuition is that regrets associated with “reasonable” edges (i.e., edges not trivially too long to be part of a good solution) should fall within the interval \(\big[0, \tfrac{1}{n} \big]\). To enforce this behavior, each regret value \(r_{ij}\) is mapped to a probability using \[\label{eq:from95regret95to95prob} p_{ij} = \max(1 - n \cdot r_{ij},\, 0),\tag{10}\] which guarantees that larger regrets lead to smaller (and eventually zero) probabilities.

In cases where all predicted regrets map to zero, we default to a uniform distribution weighted by the inverse cost, namely, \[p_{ij} = \frac{1}{w_{ij}}, \qquad \forall\, 1 \leq i < j \leq n.\]

GNN4CO:

The paper by [3] proposes two pre-trained models: one trained on instances of variable size and another trained on large cases with \(n = 200\) nodes. They also introduce two decoding strategies: an autoregressive (AR) decoder and a non-autoregressive (NAR) decoder, where the edge probabilities are predicted simultaneously. In our experiments, we made the following choices. We use only the architecture trained on variable-size instances, as it exhibits better generalization properties (see [3]). We adopt the AR decoding strategy, which yields better solution quality than the NAR decoder (see [3]). We use the same hyperparameter configuration as in the original paper and discard any predicted probabilities below 0.1, following their approach.

SoftDist:

To obtain meaningful probability values, we first normalize the edge costs. For each graph, for each edge \(e\), we define \[\bar{w}(e) := \frac{w(e)}{\max_{e'\in E} w(w)},\] so that all costs lie within the interval \([0, 1]\). To select the appropriate temperature, we choose the one corresponding to the closest instance size among those validated by [12]. For example, for an instance with \(612\) nodes, we use the temperature validated for \(n = 500\). In the case of tiled size, we select the temperature associated with the largest number of nodes, as we observed a correlation between instance size and the lowest optimal temperature (e.g, for \(n = 750\), we opt for the one associated to \(1000\)).

Synthetic predictors (SP)

Given \(\varepsilon \in (0,1)\) and letting \(OPT\) denote the optimal tour value, let \(X^*\) be the set of edges of an optimal tour. We sample a subset \(H^{-} \subseteq X^*\) such that \(\sum_{e \in H^{-}} w(e) \simeq \varepsilon \cdot OPT\) and define \(H^{+}\) analogously. We then set \(X' := (X^* \setminus H^{-}) \cup H^{+}\), and assign predictions equal to \(p = 0.9\) to the edges in \(X'\) and 0 otherwise.

For all predictors, we normalize the scores so that their maximum value is 1. An additional detail concerns normalization with respect to edge weights. For G1,G2, and BS we normalize the scores by the edge costs, as is standard in this line of work [56]. For \(\text{CHR}^+\), we do not apply this normalization. The reason is that our algorithm does not rely on heatmaps alone to construct a feasible solution, but explicitly incorporates the edge weights.

12.3 Algorithm [alg:alg1] vs \(\boldsymbol{\text{CHR}^+}\)↩︎

A careful analysis of Algorithm [alg:alg1] reveals that it requires a subset of edges \(P\) (of cardinality \(k\)) whose associated weights are set to zero. In practice, given a heatmap of edge probabilities, several strategies can be employed to derive \(P\). One approach is to sample \(k\) edges according to the edge probability distribution (Alg1), while another is to select the \(k\) edges with the highest probabilities (Alg1Top). Alternatively, rather than defining \(w'(e)\) as in Algorithm [alg:alg1], one can define the modified weight for every edge \(e\) as \(w'(e) = w(e)(1 - w(e))\) and execute Algorithm [alg:alg1]. Since this latter approach is essentially Christofides’ algorithm applied to a modified graph, we denote it as \(\text{CHR}^+\). While our theoretical analysis is developed for Alg1, we evaluate the three approaches. In all cases, we set \(k = n\). Under the Alg1 approach, we generate 30 distinct edge sets \(P\) for each instance to account for stochasticity; in contrast, the Alg1Top strategy defines \(P\) deterministically based on the provided heatmap. Using SoftDist as a representative predictor, we compute the resulting optimality gaps for each strategy to assess their relative performance. Figure 14 demonstrates that \(\text{CHR}^+\) consistently achieves superior performance. Consequently, we adopt \(\text{CHR}^+\) as the primary decoding strategy for all subsequent experiments. Figure 14 demonstrates the empirical superiority of \(\text{CHR}^+\) over Algorithm [alg:alg1] in effectively leveraging neural predictions.

12.4 Quality of the prediction error with respect to greedy↩︎

12.4.0.1 G2 decoding has worse error dependency then \(\text{CHR}^+\).

Figure 15 reports the average performance of Christofides, G2, and \(\text{CHR}^+\) on the first five instances of our dataset with \(n = 50\). The probabilities are generated using SP, and the sets \(H^+\) and \(H^-\) are sampled at random using 30 different seeds. Since the greedy approach performs significantly worse than the other methods, we exclude it from subsequent analyses. We report only G2 because our preliminary analysis indicates that G1 performs worse than G2.

12.4.0.2 \(\text{CHR}^+\) exhibits graceful degradation with respect to the prediction error.

Figure 8 shows how \(\text{CHR}^+\) performances are with respect to the prediction error \(\varepsilon\). To obtain this plot, we fix \(n = 500\) and sample \(H^+\) and \(H^-\) at random using 5 different seeds.

12.5 Comparison of different solution search methods↩︎

In this section, we report an extensive analysis of all the solution search methods on the dataset generated uniformly at random (\(\mathcal{U}\)) and the instances from TSPLIB (Euclidean, \(\mathcal{T}_E\), metric non-Euclidean, \(\mathcal{T}_M\)). See Table 1, Table 2, and Table 3. For our experiments, the beam search width is set to 50, providing an optimal trade-off between solution quality and computational efficiency. The results on dataset \(\mathcal{U}\) demonstrate that \(\text{CHR}^+\) is the most robust decoding strategy as the instance size \(n\) scales. While traditional solution search (G1, G2, BS) often sees its performance degrade significantly compared to the Christofides baseline as \(n\) increases, \(\text{CHR}^+\) consistently maintains or improves upon the baseline. Notably, for \(n=1000\), \(\text{CHR}^+\) paired with SoftDist achieves a gap of \(12.80\%\), outperforming plain Christofides (\(12.90\%\)), whereas G1 and G2 degrade to \(29.00\%\) and \(25.43\%\) respectively. When paired with the DIFUSCO predictor, \(\text{CHR}^+\) achieves state-of-the-art results on smaller instances, such as a \(0.10\%\) gap for \(n=50\). The empirical superiority of \(\text{CHR}^+\) generalizes to the Euclidean TSPLIB benchmark (dataset \(\mathcal{T}\)). Across nearly all instances, \(\text{CHR}^+\) produces the lowest integrality gaps. For example, on the pr76 instance, \(\text{CHR}^+\) with DIFUSCO yields a gap of \(1.91\%\), a significant improvement over the \(7.88\%\) gap provided by the standard Christofides algorithm. In the metric non-Euclidean setting, we can only evaluate SoftDist and GNNGLS, since the other two predictors are restricted to the Euclidean framework. In this regime, \(\text{CHR}^+\) is overall less effective. Although the performance gains are less pronounced than in the Euclidean case – primarily because the problem is more challenging – and in some instances closer to those of the plain Christofides algorithm, our CHR\(^+\) remains the most robust decoder across all instances. On average, it performs best when combined with GNNGLS, achieving an average gap of 8.49% compared to 8.70% for Christofides. This difference should be interpreted with some caution, as GNNGLS cannot compute heatmaps on large and difficult instances, which would likely increase the average gap.

12.6 2-opt improvement↩︎

We evaluate the impact of 2-opt in \(\text{CHR}^+\), G2, and BS. G1 is omitted due to its documented sub-optimality relative to G2; similarly, we restrict our predictors to SoftDist and DIFUSCO. The results for dataset \(\mathcal{U}\) are summarized in Table 4. \(\text{CHR}^+\) demonstrates superior performance across all dimensions except \(n = 50\), where BS achieves a lower gap—primarily due to a superior initialization tour. Crucially, the number of 2-opt swaps required to reach local optimality is comparable between \(\text{CHR}^+\) and the plain Christofides baseline. When considering the runtimes reported in Table 1, it is evident that \(\text{CHR}^+\) has, on average, the same running time as Christofides. Consequently, the \(\text{CHR}^+\) + 2-opt pipeline yields significantly higher solution quality than Christofides at an equivalent computational cost.

References↩︎

---

1. Without the triangle inequality, TSP is known to be inapproximable, i.e., no polynomial-time algorithm can achieve a finite approximation ratio for non-metric TSP, unless P=NP [13]. Unfortunately, this hardness persists even when predictions with arbitrarily small error are available, see Appendix 11.3.↩︎

2. One can also solve the problem in time \(O(n^{2.5}\log W)\) where \(W\) is the largest integer weight [50].↩︎

3. https://zenodo.org/records/20283463↩︎

4. https://github.com/eleonoravercesi/tour_from_neural_predictions↩︎