我们终于从 Held Karp 松弛问题中跳出来，来到了 Asadpour 非对称旅行商问题算法的第二步。参考 Asadpour 论文中的算法 1，我们现在终于来到了第二步。
算法 1 一个 $O(\log n / \log \log n)$-逼近算法用于 ATSP
输入： 一个包含 $n$ 个点的集合 $V$ 和一个满足三角不等式的成本函数 $c\ :\ V \times V \rightarrow \mathbb{R}^+$。
输出： 由 $V$ 和 $c$ 描述的非对称旅行商问题实例的 $O(\log n / \log \log n)$-逼近。
- 求解 ATSP 实例的 Held-Karp LP 松弛问题，得到最优极值点解 $x^*$。定义 $z^*$ 如 (5) 所示，使其成为 $x^*$ 的对称化和缩小版本。向量 $z^*$ 可以被视为 $x^*$ 支撑图上的无向图的生成树多面体中的一个点，该图在忽略弧的方向后获得（参见第 3 节）。
- 设 $E$ 是在忽略弧方向后 $z^*$ 的支撑图。找到权重 ${\tilde{\gamma}}_{e \in E}$，使得生成树上的指数分布，$\tilde{p}(T) \propto \exp(\sum_{e \in T} \tilde{\gamma}_e)$ （近似）保留由 $z^*$ 强加的边际，即对于任何边 $e \in E$， $$\sum_{T \in \mathcal{T} : T \ni e} \tilde{p}(T) \leq (1 + \epsilon) z^*_e$$ 对于足够小的 $\epsilon$ 值。（在本文中，我们证明 $\epsilon = 0.2$ 满足我们的目的。参见第 7 节和第 8 节，了解如何计算这种分布。）
- 从 $\tilde{p}(.)$ 中采样 $2\lceil \log n \rceil$ 个生成树 $T_1, \dots, T_{2\lceil \log n \rceil}$。对于这些树中的每一个，对所有边进行定向，以使其相对于我们（非对称）成本函数 $c$ 的成本最小。设 $T^*$ 是在所有采样树中成本最小的树。
- 找到包含定向树 $\vec{T}^*$ 的最小成本积分循环。将此循环缩短为一个回路并输出。（参见第 4 节）。
第 7 节和第 8 节提供了两种不同的方法来找到所需的概率分布，其中第 7 节使用组合方法，第 8 节使用椭球法。考虑到科学 Python 生态系统中没有椭球求解器，并且我的导师和我已经决定不在此项目中实现一个，我将使用第 7 节中的方法。
第 7 节给出的算法如下
- 设置 $\gamma = \vec{0}$。
- 当存在边 $e$ 使得 $q_e(\gamma) > (1 + \epsilon) z_e$ 时
- 计算 $\delta$，如果我们将 $\gamma’$ 定义为 $\gamma_e’ = \gamma_e - \delta$，并且对于所有 $f \in E\ \backslash {e}$，$\gamma_f’ = \gamma_f$，那么 $q_e(\gamma’) = (1 + \epsilon/2)z_e$。
- 设置 $\gamma \leftarrow \gamma’$。
- 输出 $\tilde{\gamma} := \gamma$。
这种结构相当简单，但我们需要知道 $q_e(\gamma)$ 是什么以及如何计算 $\delta$。
找到 $\delta$ 非常容易，公式在 Asadpour 论文中给出（虽然我在写我的 GSoC 提案时并没有意识到这一点，并重新推导了 delta 的方程。幸运的是，我的公式与论文中的公式相符。）
$$ \delta = \ln \frac{q_e(\gamma)(1 - (1 + \epsilon / 2)z_e)}{(1 - q_e(\gamma))(1 + \epsilon / 2) z_e} $$
注意，$\delta$ 的公式依赖于 $q_e(\gamma)$。论文将 $q_e(\gamma)$ 定义为
$$ q_e(\gamma) = \frac{\sum_{T \ni e} \exp(\gamma(T))}{\sum_{T \in \mathcal{T}} \exp(\gamma(T))} $$
其中 $\gamma(T) = \sum_{f \in T} \gamma_f$。
我首先注意到，在分母中，求和是对图中所有生成树进行的，对于我们将要处理的完全图，这将是指数级，因此在这里使用“暴力”方法是无用的。幸运的是，Asadpour 和他的团队意识到我们可以利用基尔霍夫矩阵树定理的优势。
关于基尔霍夫矩阵树定理的旁注，我在此项目之前不熟悉该定理，因此我不得不阅读一些关于它的内容。基本上，如果你有一个拉普拉斯矩阵（邻接矩阵减去度矩阵），则任何余子式的绝对值都是图中生成树的数目。这对我来说完全出乎意料，我认为这种类型的联系的存在非常酷。
使用基尔霍夫定理的细节在第 5.3 节中给出。我们将使用由以下公式定义的加权拉普拉斯矩阵 $L$
$$ L_{i, j} = \left\{ \begin{array}{l l} -\lambda_e & e = (i, j) \in E \\\ \sum_{e \in \delta({i})} \lambda_e & i = j \\\ 0 & \text{otherwise} \end{array} \right. $$
其中 $\lambda_e = \exp(\gamma_e)$。
现在，我们知道将 Krichhoff 定理应用于 $L$ 将返回
$$ \sum_{t \in \mathcal{T}} \prod_{e \in T} \lambda_e $$
但 $q_e(\gamma)$ 的哪一部分是它？
如果我们应用 $\lambda_e = \exp(\gamma_e)$，我们发现
$$ \begin{array}{r c l} \sum_{T \in \mathcal{T}} \prod_{e \in T} \lambda_e &=& \sum_{T \in \mathcal{T}} \prod_{e \in T} \exp(\gamma_e) \\\ && \sum_{T \in \mathcal{T}} \exp\left(\sum_{e \in T} \gamma_e\right) \\\ && \sum_{T \in \mathcal{T}} \exp(\gamma(T)) \\\ \end{array} $$
因此，从第一行到第二行是一个令人困惑的步骤，但本质上我们正在利用指数的性质。回想一下 $\exp(x) = e^x$，所以我们可以将其写成 $\prod_{e \in T} e^{\gamma_e}$，但这会引入歧义，因为我们将有多个 $e$ 的含义。现在，对于所有 $e$ 的值，$e_1, e_2, \dots, e_{n-1}$ 在生成树 $T$ 中，该乘积可以展开为
$$ \prod_{e \in T} e^{\gamma_e} = e^{\gamma_{e_1}} \times e^{\gamma_{e_2}} \times \dots \times e^{\gamma_{e_{n-1}}} $$
每个指数因子具有相同的底数，因此我们可以将其折叠成
$$ e^{\gamma_{e_1} + \gamma_{e_2} + \dots + \gamma_{e_{n-1}}} $$
这也等于
$$ e^{\sum_{e \in T} \gamma_e} $$
但我们知道 $\sum_{e \in T} \gamma_e$ 是 $\gamma(T)$，因此它变成
$$ e^{\gamma(T)} = \exp(\gamma(T)) $$
一旦我们将其放回求和中，我们就得到了 $q_e(\gamma)$ 的分母，$\sum_{T \in \mathcal{T}} \exp(\gamma(T))$。
接下来，我们需要找到 $q_e(\gamma)$ 的分子。与之前一样，“暴力”方法的复杂度将是指数级的，因此我们必须找到更好的方法。好吧，分子和分母之间唯一的区别在于外层求和的条件，即 $T \in \mathcal{T}$ 被更改为 $T \ni e$ 或包含边 $e$ 的所有树。
这里也有一种方法可以使用基尔霍夫矩阵树定理。如果我们有一个图，其中每个生成树可以一对一地映射到原始图中的每个包含所需边 $e$ 的生成树。为了使生成树包含边 $e$，我们知道 $e$ 的端点 $(u, v)$ 将直接相互连接。因此，我们感兴趣的是所有到达顶点 $u$ 然后从顶点 $v$ 离开的生成树。（与到达顶点 $u$ 然后从同一个顶点离开的生成树相反）。从某种意义上说，我们将顶点 $u$ 和 $v$ 视为同一个顶点。我们可以通过从图中收缩 $e$ 来实际应用这一点，从而创建 $G / {e}$。此图中的每个生成树都可以从 $G / {e}$ 中唯一映射到 $G$ 中包含边 $e$ 的生成树。
从这里开始，证明 $L$ 的余子式实际上是 $q_e(\gamma)$ 的分子的逻辑与分母的逻辑是平行的。
此时，我们拥有创建 Asadpour 方法中下一个函数 spanning_tree_distribution() 的伪代码所需的所有信息。这里我将使用内部函数 q() 来查找 $q_e$。
def spanning_tree_distribution
input: z, the symmetrized and scaled output of the Held Karp relaxation.
output: gamma, the maximum entropy exponential distribution for sampling spanning trees
from the graph.
def q
input: e, the edge of interest
# Create the laplacian matrices
write lambda = exp(gamma) into the edges of G
G_laplace = laplacian(G, lambda)
G_e = nx.contracted_edge(G, e)
G_e_laplace = laplacian(G, lambda)
# Delete a row and column from each matrix to made a cofactor matrix
G_laplace.delete((0, 0))
G_e_laplace.delete((0, 0))
# Calculate the determinant of the cofactor matrices
det_G_laplace = G_laplace.det
det_G_e_laplace = G_e_laplace.det
# return q_e
return det_G_e_laplace / det_G_laplace
# initialize the gamma vector
gamma = 0 vector of length G.size
while true
# We will iterate over the edges in z until we complete the
# for loop without changing a value in gamma. This will mean
# that there is not an edge with q_e > 1.2 * z_e
valid_count = 0
# Search for an edge with q_e > 1.2 * z_e
for e in z
q_e = q(e)
z_e = z[e]
if q_e > 1.2 * z_e
delta = ln(q_e * (1 - 1.1 * z_e) / (1 - q_e) * 1.1 * z_e)
gamma[e] -= delta
else
valid_count += 1
if valid_count == number of edges in z
break
return gamma
下一步#
下一步是使用上面的伪代码作为大纲来实现函数 spanning_tree_distribution。我将从编写 q 开始，并使用我用来测试 Held Karp 松弛的相同图来对其进行测试。一旦 q 完成，函数的其余部分似乎相当简单。
我担心的一件事是我测试 spanning_tree_distribution 的能力。Asadpour 研究论文中没有给出任何例子，也没有其他容易获得的资源可以让我用来找到一个预言机。
我现在能想到的唯一方法是完成此函数，然后完成 sample_spanning_tree。一旦这两个函数都完成，我就可以采样大量的生成树，找到每棵树的实验概率，然后运行统计检验（例如 h-检验）来查看每棵树的概率是否接近 $\exp(\gamma(T))$，这就是所需的分布。另一种测试方法是使用分布中的边际，并手动检查
$$ \sum_{T \in \mathcal{T} : T \ni e} p(T) \leq (1 + \epsilon) z^*_e,\ \forall\ e \in E $$
其中 $p(T)$ 是来自采样树的实验数据。
这两种方法似乎都非常计算密集，并且由于它们是从概率分布中进行采样，因此它们可能会由于样本不太可能而随机失败。
参考文献#
A. Asadpour、M. X. Goemans、A. Mardry、S. O. Ghran 和 A. Saberi，《非对称旅行商问题的 $O(\log n / \log \log n)$-逼近算法》，《运筹学》，65 (2017)，第 1043-1061 页。