scipy.optimize.

quadratic_assignment#

scipy.optimize.quadratic_assignment(A, B, method='faq', options=None)[源]#

近似求解二次分配问题和图匹配问题。

二次分配解决以下形式的问题

\[\begin{split}\min_P & \ {\ \text{trace}(A^T P B P^T)}\\ \mbox{s.t. } & {P \ \epsilon \ \mathcal{P}}\\\end{split}\]

其中 \(\mathcal{P}\) 是所有置换矩阵的集合,\(A\)\(B\) 是方阵。

图匹配尝试最大化相同的目标函数。该算法可以被视为寻找两个图的节点对齐方式,以最小化诱导边不一致的数量,或者,在加权图的情况下,最小化边权重平方差之和。

请注意,二次分配问题是 NP-hard 的。此处给出的结果是近似值,不保证是最佳的。

参数:
A2-D 数组,方阵

上述目标函数中的方阵 \(A\)

B2-D 数组,方阵

上述目标函数中的方阵 \(B\)

method字符串,取值范围 {‘faq’, ‘2opt’} (默认值: ‘faq’)

用于解决问题的算法。‘faq’(默认)和 ‘2opt’ 可用。

optionsdict,可选

求解器选项的字典。所有求解器都支持以下选项

maximize布尔值 (默认值: False)

如果为 True,则最大化目标函数。

partial_match整数的 2-D 数组,可选 (默认值: None)

固定匹配的一部分。也称为“种子”[2]

partial_match 的每一行指定一对匹配的节点:A 的节点 partial_match[i, 0]B 的节点 partial_match[i, 1] 匹配。该数组的形状为 (m, 2),其中 m 不大于节点数量 \(n\)

rngnumpy.random.Generator,可选

伪随机数生成器状态。当 rng 为 None 时,将使用操作系统的熵创建一个新的 numpy.random.Generator。非 numpy.random.Generator 类型的参数将传递给 numpy.random.default_rng 以实例化一个 Generator

版本 1.15.0 中有更改: 作为从使用 numpy.random.RandomStatenumpy.random.GeneratorSPEC-007 过渡的一部分。现在向此函数提供 np.random.RandomState 将发出 DeprecationWarning。在 SciPy 1.17 中,它的使用将引发异常。此外,依赖使用 np.random.seed 的全局状态将发出 FutureWarning。在 SciPy 1.17 中,将不再使用全局随机数生成器。使用整数类种子将引发 FutureWarning,在 SciPy 1.17 中,它将通过 np.random.default_rng 而非 np.random.RandomState 进行规范化。

有关特定于方法的选项,请参阅 show_options('quadratic_assignment')

返回:
resOptimizeResult

OptimizeResult 包含以下字段。

col_ind1-D 数组

与找到的 B 节点最佳置换相对应的列索引。

fun浮点数

解的目标值。

nit整数

优化期间执行的迭代次数。

说明

默认方法 ‘faq’ 使用快速近似 QAP 算法 [1];它通常提供速度和准确性的最佳组合。方法 ‘2opt’ 计算量可能很大,但可能是一个有用的替代方法,或者它可以用于改进由其他方法返回的解决方案。

参考文献

[1]

J.T. Vogelstein, J.M. Conroy, V. Lyzinski, L.J. Podrazik, S.G. Kratzer, E.T. Harley, D.E. Fishkind, R.J. Vogelstein, and C.E. Priebe, “Fast approximate quadratic programming for graph matching,” PLOS one, vol. 10, no. 4, p. e0121002, 2015, DOI:10.1371/journal.pone.0121002

[2]

D. Fishkind, S. Adali, H. Patsolic, L. Meng, D. Singh, V. Lyzinski, C. Priebe, “Seeded graph matching”, Pattern Recognit. 87 (2019): 203-215, DOI:10.1016/j.patcog.2018.09.014

[3]

“2-opt,” Wikipedia. https://en.wikipedia.org/wiki/2-opt

示例

>>> import numpy as np
>>> from scipy.optimize import quadratic_assignment
>>> rng = np.random.default_rng()
>>> A = np.array([[0, 80, 150, 170], [80, 0, 130, 100],
...               [150, 130, 0, 120], [170, 100, 120, 0]])
>>> B = np.array([[0, 5, 2, 7], [0, 0, 3, 8],
...               [0, 0, 0, 3], [0, 0, 0, 0]])
>>> res = quadratic_assignment(A, B, options={'rng': rng})
>>> print(res)
     fun: 3260
 col_ind: [0 3 2 1]
     nit: 9

要查看返回的 col_indfun 之间的关系,请使用 col_ind 形成找到的最佳置换矩阵,然后评估目标函数 \(f(P) = trace(A^T P B P^T )\)

>>> perm = res['col_ind']
>>> P = np.eye(len(A), dtype=int)[perm]
>>> fun = np.trace(A.T @ P @ B @ P.T)
>>> print(fun)
3260

或者,为了避免显式构建置换矩阵,可以直接置换距离矩阵的行和列。

>>> fun = np.trace(A.T @ B[perm][:, perm])
>>> print(fun)
3260

尽管通常不保证,但 quadratic_assignment 恰好找到了全局最优解。

>>> from itertools import permutations
>>> perm_opt, fun_opt = None, np.inf
>>> for perm in permutations([0, 1, 2, 3]):
...     perm = np.array(perm)
...     fun = np.trace(A.T @ B[perm][:, perm])
...     if fun < fun_opt:
...         fun_opt, perm_opt = fun, perm
>>> print(np.array_equal(perm_opt, res['col_ind']))
True

这是一个默认方法 ‘faq’ 未找到全局最优解的示例。

>>> A = np.array([[0, 5, 8, 6], [5, 0, 5, 1],
...               [8, 5, 0, 2], [6, 1, 2, 0]])
>>> B = np.array([[0, 1, 8, 4], [1, 0, 5, 2],
...               [8, 5, 0, 5], [4, 2, 5, 0]])
>>> res = quadratic_assignment(A, B, options={'rng': rng})
>>> print(res)
     fun: 178
 col_ind: [1 0 3 2]
     nit: 13

如果精度很重要,请考虑使用 ‘2opt’ 来改进解决方案。

>>> guess = np.array([np.arange(len(A)), res.col_ind]).T
>>> res = quadratic_assignment(A, B, method="2opt",
...     options = {'rng': rng, 'partial_guess': guess})
>>> print(res)
     fun: 176
 col_ind: [1 2 3 0]
     nit: 17