二次分配#
- scipy.optimize.quadratic_assignment(A, B, method='faq', options=None)[源代码][源代码]#
近似求解二次分配问题和图匹配问题。
二次分配解决以下形式的问题:
\[\min_P & \ {\ \text{trace}(A^T P B P^T)}\ \mbox{s.t. } & {P \ \epsilon \ \mathcal{P}}\]其中 \(\mathcal{P}\) 是所有置换矩阵的集合,而 \(A\) 和 \(B\) 是方阵。
图匹配试图 最大化 相同的目標函數。這個算法可以被認為是找到兩個圖的節點對齊,這最小化了誘導邊的不一致數量,或者在加權圖的情況下,最小化了邊權重差的平方和。
注意,二次分配问题是NP难的。这里给出的结果是近似值,不保证是最优的。
- 参数:
- A2-D 数组, 方阵
上述目标函数中的方阵 \(A\)。
- B2-D 数组, 方阵
目标函数中的方阵 \(B\)。
- 方法str in {‘faq’, ‘2opt’} (默认: ‘faq’)
用于解决问题的算法。可用的有 ‘faq’ <optimize.qap-faq>`(默认)和 :ref:’2opt’ <optimize.qap-2opt>`。
- 选项dict, 可选
求解器选项的字典。所有求解器都支持以下选项:
- 最大化bool (默认: False)
如果为
True
,则最大化目标函数。- 部分匹配整数的二维数组,可选(默认:无)
修复了部分匹配问题。也称为“种子” [2]。
partial_match 的每一行指定了一对匹配的节点:A 的节点
partial_match[i, 0]
与 B 的节点partial_match[i, 1]
匹配。该数组的形状为(m, 2)
,其中m
不大于节点数 \(n\)。- rng{None, int,}
如果 seed 是 None(或 np.random),则使用
numpy.random.RandomState
单例。如果 seed 是 int,则使用新的RandomState
实例,并以 seed 为种子。如果 seed 已经是Generator
或RandomState
实例,则使用该实例。
对于特定方法的选项,请参见
show_options('quadratic_assignment')
。
- 返回:
- res优化结果
OptimizeResult
包含以下字段。- col_ind一维数组
对应于找到的节点 B 的最佳排列的列索引。
- 乐趣浮动
解决方案的目标值。
- 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, 和 C.E. Priebe, “快速近似二次规划用于图匹配,” PLOS one, 第10卷, 第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,” 维基百科。 https://en.wikipedia.org/wiki/2-opt
示例
>>> import numpy as np >>> from scipy.optimize import quadratic_assignment >>> 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) >>> print(res) fun: 3260 col_ind: [0 3 2 1] nit: 9
要查看返回的
col_ind
和fun
之间的关系,使用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) >>> 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 = {'partial_guess': guess}) >>> print(res) fun: 176 col_ind: [1 2 3 0] nit: 17