Source code for pygmtools.multi_graph_solvers

import functools
import importlib
import pygmtools
from pygmtools.utils import NOT_IMPLEMENTED_MSG, _check_shape, _get_shape, _unsqueeze, _squeeze, _check_data_type
import math


[docs]def cao(K, x0=None, qap_solver=None, mode='accu', max_iter=6, lambda_init=0.3, lambda_step=1.1, lambda_max=1.0, iter_boost=2, backend=None): r""" Composition based Affinity Optimization (CAO) solver for multi-graph matching. This solver builds a supergraph for matching update to incorporate the two aspects by optimizing the affinity score, meanwhile gradually infusing the consistency. Each update step is described as follows: .. math:: \arg \max_{k} (1-\lambda) J(\mathbf{X}_{ik} \mathbf{X}_{kj}) + \lambda C_p(\mathbf{X}_{ik} \mathbf{X}_{kj}) where :math:`J(\mathbf{X}_{ik} \mathbf{X}_{kj})` is the objective score, and :math:`C_p(\mathbf{X}_{ik} \mathbf{X}_{kj})` measures a consistency score compared to other matchings. These two terms are balanced by :math:`\lambda`, and :math:`\lambda` starts from a smaller number and gradually grows. :param K: :math:`(m\times m \times n^2 \times n^2)` the input affinity matrix, where ``K[i,j]`` is the affinity matrix of graph ``i`` and graph ``j`` (:math:`m`: number of nodes) :param x0: (optional) :math:`(m\times m \times n \times n)` the initial two-graph matching result, where ``X[i,j]`` is the matching matrix result of graph ``i`` and graph ``j``. If this argument is not given, ``qap_solver`` will be used to compute the two-graph matching result. :param qap_solver: (default: pygm.rrwm) a function object that accepts a batched affinity matrix and returns the matching matrices. It is suggested to use ``functools.partial`` and the QAP solvers provided in the :mod:`~pygmtools.classic_solvers` module (see examples below). :param mode: (default: ``'accu'``) the operation mode of this algorithm. Options: ``'accu', 'c', 'fast', 'pc'``, where ``'accu'`` is equivalent to ``'c'`` (accurate version) and ``'fast'`` is equivalent to ``'pc'`` (fast version). :param max_iter: (default: 6) max number of iterations :param lambda_init: (default: 0.3) initial value of :math:`\lambda`, with :math:`\lambda\in[0,1]` :param lambda_step: (default: 1.1) the increase step size of :math:`\lambda`, updated by ``lambda = step * lambda`` :param lambda_max: (default: 1.0) the max value of lambda :param iter_boost: (default: 2) to boost the convergence of the CAO algorithm, :math:`\lambda` will be forced to update every ``iter_boost`` iterations. :param backend: (default: ``pygmtools.BACKEND`` variable) the backend for computation. :return: :math:`(m\times m \times n \times n)` the multi-graph matching result .. note:: The input graphs must have the same number of nodes for this algorithm to work correctly. .. note:: Multi-graph matching methods process all graphs at once and do not support the additional batch dimension. Please note that this behavior is different from two-graph matching solvers in :mod:`~pygmtools.classic_solvers`. Example for Pytorch backend:: >>> import torch >>> import pygmtools as pygm >>> pygm.BACKEND = 'pytorch' >>> _ = torch.manual_seed(1) # Generate 10 isomorphic graphs >>> graph_num = 10 >>> As, X_gt = pygm.utils.generate_isomorphic_graphs(node_num=4, graph_num=10) >>> As_1, As_2 = [], [] >>> for i in range(graph_num): ... for j in range(graph_num): ... As_1.append(As[i]) ... As_2.append(As[j]) >>> As_1 = torch.stack(As_1, dim=0) >>> As_2 = torch.stack(As_2, dim=0) # Build affinity matrix >>> conn1, edge1, ne1 = pygm.utils.dense_to_sparse(As_1) >>> conn2, edge2, ne2 = pygm.utils.dense_to_sparse(As_2) >>> import functools >>> gaussian_aff = functools.partial(pygm.utils.gaussian_aff_fn, sigma=1.) # set affinity function >>> K = pygm.utils.build_aff_mat(None, edge1, conn1, None, edge2, conn2, None, None, None, None, edge_aff_fn=gaussian_aff) >>> K = K.reshape(graph_num, graph_num, 4*4, 4*4) >>> K.shape torch.Size([10, 10, 16, 16]) # Solve the multi-matching problem >>> X = pygm.cao(K) >>> (X * X_gt).sum() / X_gt.sum() tensor(1.) # Use the IPFP solver for two-graph matching >>> ipfp_func = functools.partial(pygmtools.ipfp, n1max=4, n2max=4) >>> X = pygm.cao(K, qap_solver=ipfp_func) >>> (X * X_gt).sum() / X_gt.sum() tensor(1.) # Run the faster version of CAO algorithm >>> X = pygm.cao(K, mode='fast') >>> (X * X_gt).sum() / X_gt.sum() tensor(1.) .. note:: If you find this graph matching solver useful in your research, please cite: :: @article{cao, title={Multi-graph matching via affinity optimization with graduated consistency regularization}, author={Yan, Junchi and Cho, Minsu and Zha, Hongyuan and Yang, Xiaokang and Chu, Stephen M}, journal={IEEE transactions on pattern analysis and machine intelligence}, volume={38}, number={6}, pages={1228--1242}, year={2015}, publisher={IEEE} } """ if backend is None: backend = pygmtools.BACKEND # check the correctness of input _check_data_type(K, backend) K_shape = _get_shape(K, backend) if not (len(K_shape) == 4 and K_shape[0] == K_shape[1] and K_shape[2] == K_shape[3]): raise ValueError(f"Unsupported input data shape: got K {K_shape}") num_graph, aff_size = K_shape[0], K_shape[2] num_node = int(math.sqrt(aff_size)) if not num_node ** 2 == aff_size: raise ValueError("The input affinity matrix is not supported. Please note that this function " "does not support matching with outliers or partial matching.") if not 0 <= lambda_init <= 1: raise ValueError(f"lambda_init must be in [0, 1], got lambda_init={lambda_init}") if not 0 <= lambda_max <= 1: raise ValueError(f"lambda_max must be in [0, 1], got lambda_max={lambda_max}") if not lambda_step > 1: raise ValueError(f"lambda_step must be >1, got lambda_step={lambda_step}") if x0 is not None: _check_data_type(x0, backend) x0_shape = _get_shape(x0, backend) if not len(x0_shape) == 4 and num_graph == x0_shape[0] == x0_shape[1] and num_node == x0_shape[2] == x0_shape[3]: raise ValueError(f"Unsupported input data shape: got K {K_shape} x0 {x0_shape}") else: if qap_solver is None: qap_solver = functools.partial(pygmtools.rrwm, n1max=num_node, n2max=num_node, backend=backend) x0 = qap_solver(K.reshape(num_graph ** 2, aff_size, aff_size)) x0 = pygmtools.hungarian(x0, backend=backend) x0 = x0.reshape(num_graph, num_graph, num_node, num_node) args = (K, x0, num_graph, num_node, max_iter, lambda_init, lambda_step, lambda_max, iter_boost) try: mod = importlib.import_module(f'pygmtools.{backend}_backend') if mode in ['pc', 'fast']: fn = mod.cao_fast_solver elif mode in ['c', 'accu']: fn = mod.cao_solver else: raise ValueError("Unknown value of mode: supported values ['c', 'accu', 'pc', 'fast']") except ModuleNotFoundError and AttributeError: raise NotImplementedError( NOT_IMPLEMENTED_MSG.format(backend) ) return fn(*args)
[docs]def mgm_floyd(K, x0=None, qap_solver=None, mode='accu', param_lambda=0.2, backend=None): r""" Multi-Graph Matching based on Floyd shortest path algorithm. A supergraph is considered by regarding each input graph as a node, and the matching between graphs are regraded as edges in the supergraph. Floyd algorithm is used to discover a shortest path on this supergraph for matching update. The length of edges on the supergraph is described as follows: .. math:: \arg \max_{k} (1-\lambda) J(\mathbf{X}_{ik} \mathbf{X}_{kj}) + \lambda C_p(\mathbf{X}_{ik} \mathbf{X}_{kj}) where :math:`J(\mathbf{X}_{ik} \mathbf{X}_{kj})` is the objective score, and :math:`C_p(\mathbf{X}_{ik} \mathbf{X}_{kj})` measures a consistency score compared to other matchings. These two terms are balanced by :math:`\lambda`. :param K: :math:`(m\times m \times n^2 \times n^2)` the input affinity matrix, where ``K[i,j]`` is the affinity matrix of graph ``i`` and graph ``j`` (:math:`m`: number of nodes) :param x0: (optional) :math:`(m\times m \times n \times n)` the initial two-graph matching result, where ``X[i,j]`` is the matching matrix result of graph ``i`` and graph ``j``. If this argument is not given, ``qap_solver`` will be used to compute the two-graph matching result. :param qap_solver: (default: pygm.rrwm) a function object that accepts a batched affinity matrix and returns the matching matrices. It is suggested to use ``functools.partial`` and the QAP solvers provided in the :mod:`~pygmtools.classic_solvers` module (see examples below). :param mode: (default: ``'accu'``) the operation mode of this algorithm. Options: ``'accu', 'c', 'fast', 'pc'``, where ``'accu'`` is equivalent to ``'c'`` (accurate version) and ``'fast'`` is equivalent to ``'pc'`` (fast version). :param param_lambda: (default: 0.3) value of :math:`\lambda`, with :math:`\lambda\in[0,1]` :param backend: (default: ``pygmtools.BACKEND`` variable) the backend for computation. :return: :math:`(m\times m \times n \times n)` the multi-graph matching result Example for Pytorch backend:: >>> import torch >>> import pygmtools as pygm >>> pygm.BACKEND = 'pytorch' >>> _ = torch.manual_seed(1) # Generate 10 isomorphic graphs >>> graph_num = 10 >>> As, X_gt = pygm.utils.generate_isomorphic_graphs(node_num=4, graph_num=10) >>> As_1, As_2 = [], [] >>> for i in range(graph_num): ... for j in range(graph_num): ... As_1.append(As[i]) ... As_2.append(As[j]) >>> As_1 = torch.stack(As_1, dim=0) >>> As_2 = torch.stack(As_2, dim=0) # Build affinity matrix >>> conn1, edge1, ne1 = pygm.utils.dense_to_sparse(As_1) >>> conn2, edge2, ne2 = pygm.utils.dense_to_sparse(As_2) >>> import functools >>> gaussian_aff = functools.partial(pygm.utils.gaussian_aff_fn, sigma=1.) # set affinity function >>> K = pygm.utils.build_aff_mat(None, edge1, conn1, None, edge2, conn2, None, None, None, None, edge_aff_fn=gaussian_aff) >>> K = K.reshape(graph_num, graph_num, 4*4, 4*4) >>> K.shape torch.Size([10, 10, 16, 16]) # Solve the multi-matching problem >>> X = pygm.mgm_floyd(K) >>> (X * X_gt).sum() / X_gt.sum() tensor(1.) # Use the IPFP solver for two-graph matching >>> ipfp_func = functools.partial(pygmtools.ipfp, n1max=4, n2max=4) >>> X = pygm.mgm_floyd(K, qap_solver=ipfp_func) >>> (X * X_gt).sum() / X_gt.sum() tensor(1.) # Run the faster version of CAO algorithm >>> X = pygm.mgm_floyd(K, mode='fast') >>> (X * X_gt).sum() / X_gt.sum() tensor(1.) .. note:: If you find this graph matching solver useful in your research, please cite: :: @article{mgm_floyd, title={Unifying offline and online multi-graph matching via finding shortest paths on supergraph}, author={Jiang, Zetian and Wang, Tianzhe and Yan, Junchi}, journal={IEEE transactions on pattern analysis and machine intelligence}, volume={43}, number={10}, pages={3648--3663}, year={2020}, publisher={IEEE} } """ if backend is None: backend = pygmtools.BACKEND # check the correctness of input _check_data_type(K, backend) K_shape = _get_shape(K, backend) if not (len(K_shape) == 4 and K_shape[0] == K_shape[1] and K_shape[2] == K_shape[3]): raise ValueError(f"Unsupported input data shape: got K {K_shape}") num_graph, aff_size = K_shape[0], K_shape[2] num_node = int(math.sqrt(aff_size)) if not num_node ** 2 == aff_size: raise ValueError("The input affinity matrix is not supported. Please note that this function " "does not support matching with outliers or partial matching.") if not 0 <= param_lambda <= 1: raise ValueError(f"param_lambda must be in [0, 1], got param_lambda={param_lambda}") if x0 is not None: _check_data_type(x0, backend) x0_shape = _get_shape(x0, backend) if not len(x0_shape) == 4 and num_graph == x0_shape[0] == x0_shape[1] and num_node == x0_shape[2] == x0_shape[3]: raise ValueError(f"Unsupported input data shape: got K {K_shape} x0 {x0_shape}") else: if qap_solver is None: qap_solver = functools.partial(pygmtools.rrwm, n1max=num_node, n2max=num_node, backend=backend) x0 = qap_solver(K.reshape(num_graph ** 2, aff_size, aff_size)) x0 = pygmtools.hungarian(x0, backend=backend) x0 = x0.reshape(num_graph, num_graph, num_node, num_node) args = (K, x0, num_graph, num_node, param_lambda) try: mod = importlib.import_module(f'pygmtools.{backend}_backend') if mode in ['pc', 'fast']: fn = mod.mgm_floyd_fast_solver elif mode in ['c', 'accu']: fn = mod.mgm_floyd_solver else: raise ValueError("Unknown value of mode: supported values ['c', 'accu', 'pc', 'fast']") except ModuleNotFoundError and AttributeError: raise NotImplementedError( NOT_IMPLEMENTED_MSG.format(backend) ) return fn(*args)
[docs]def gamgm(A, W, ns=None, n_univ=None, U0=None, sk_init_tau=0.5, sk_min_tau=0.1, sk_gamma=0.8, sk_iter=20, max_iter=100, param_lambda=1., verbose=False, backend=None): r""" Graduated Assignment-based multi-graph matching solver. Graduated assignment is a classic approach for hard assignment problems like graph matching, based on graduated annealing of Sinkhorn's temperature :math:`\tau` to enforce the matching constraint. The objective score is described as .. math:: \max_{\mathbf{X}_{i,j}, i,j\in [m]} \ \sum_{i,j\in [m]} \left( \lambda \ \mathrm{tr}(\mathbf{X}_{ij}^\top \mathbf{A}_{i} \mathbf{X}_{ij} \mathbf{A}_{j}) + \mathrm{tr}(\mathbf{X}_{ij}^\top \mathbf{W}_{ij})\right) Once the algorithm converges at a fixed :math:`\tau` value, :math:`\tau` shrinks as: .. math:: \tau = \tau \times \gamma and the iteration continues. At last, Hungarian algorithm is applied to ensure the result is a permutation matrix. .. note:: This algorithm is based on the Koopmans-Beckmann's QAP formulation and you should input the adjacency matrices ``A`` and node-wise similarity matrices ``W`` instead of the affinity matrices. :param A: :math:`(m\times n \times n)` the adjacency matrix (:math:`m`: number of nodes). The graphs may have different number of nodes (specified by the ``ns`` argument). :param W: :math:`(m\times m \times n \times n)` the node-wise similarity matrix, where ``W[i,j]`` is the similarity matrix :param ns: (optional) :math:`(m)` the number of nodes. If not given, it will be inferred based on the size of ``A``. :param n_univ: (optional) the size of the universe node set. If not given, it will be the largest number of nodes. :param U0: (optional) the initial multi-graph matching result. If not given, it will be randomly initialized. :param sk_init_tau: (default: 0.05) initial value of :math:`\tau` for Sinkhorn algorithm :param sk_min_tau: (default: 1.0e-3) minimal value of :math:`\tau` for Sinkhorn algorithm :param sk_gamma: (default: 0.8) the shrinking parameter of :math:`\tau`: :math:`\tau = \tau \times \gamma` :param sk_iter: (default: 200) max number of iterations for Sinkhorn algorithm :param max_iter: (default: 1000) max number of iterations for graduated assignment :param param_lambda: (default: 1) the weight :math:`\lambda` of the quadratic term :param verbose: (default: False) print verbose information for parameter tuning :param backend: (default: ``pygmtools.BACKEND`` variable) the backend for computation. :return: the multi-graph matching result (a :mod:`~pygmtools.utils.MultiMatchingResult` object) .. note:: Set ``verbose=True`` may help you tune the parameters. Example for Pytorch backend:: >>> import torch >>> import pygmtools as pygm >>> import itertools >>> import time >>> pygm.BACKEND = 'pytorch' >>> _ = torch.manual_seed(1) # Generate 10 isomorphic graphs >>> graph_num = 10 >>> As, X_gt, Fs = pygm.utils.generate_isomorphic_graphs(node_num=4, graph_num=10, node_feat_dim=20) # Compute node-wise similarity by inner-product and Sinkhorn >>> W = torch.matmul(Fs.unsqueeze(1), Fs.transpose(1, 2).unsqueeze(0)) >>> W = pygm.sinkhorn(W.reshape(graph_num ** 2, 4, 4)).reshape(graph_num, graph_num, 4, 4) # Solve the multi-matching problem >>> X = pygm.gamgm(As, W) >>> matched = 0 >>> for i, j in itertools.product(range(graph_num), repeat=2): ... matched += (X[i,j] * X_gt[i,j]).sum() >>> matched / X_gt.sum() tensor(1.) # This function supports graphs with different nodes (also known as partial matching) # In the following we ignore the last node from the last 5 graphs >>> ns = torch.tensor([4, 4, 4, 4, 4, 3, 3, 3, 3, 3], dtype=torch.int) >>> for i in range(graph_num): ... As[i, ns[i]:, :] = 0 ... As[i, :, ns[i]:] = 0 >>> for i, j in itertools.product(range(graph_num), repeat=2): ... X_gt[i, j, ns[i]:, :] = 0 ... X_gt[i, j, :, ns[j]:] = 0 ... W[i, j, ns[i]:, :] = 0 ... W[i, j, :, ns[j]:] = 0 # Partial matching is challenging and the following parameters are carefully tuned >>> X = pygm.gamgm(As, W, ns, n_univ=4, sk_init_tau=.1, sk_min_tau=0.01, param_lambda=0.3) # Check the partial matching result >>> matched = 0 >>> for i, j in itertools.product(range(graph_num), repeat=2): ... matched += (X[i,j] * X_gt[i, j, :ns[i], :ns[j]]).sum() >>> matched / X_gt.sum() tensor(1.) .. note:: If you find this graph matching solver useful in your research, please cite: :: @article{gamgm1, title={Graduated assignment algorithm for multiple graph matching based on a common labeling}, author={Sol{\'e}-Ribalta, Albert and Serratosa, Francesc}, journal={International Journal of Pattern Recognition and Artificial Intelligence}, volume={27}, number={01}, pages={1350001}, year={2013}, publisher={World Scientific} } @article{gamgm2, title={Graduated assignment for joint multi-graph matching and clustering with application to unsupervised graph matching network learning}, author={Wang, Runzhong and Yan, Junchi and Yang, Xiaokang}, journal={Advances in Neural Information Processing Systems}, volume={33}, pages={19908--19919}, year={2020} } This algorithm is originally proposed by paper ``gamgm1``, and further improved by paper ``gamgm2`` to fit modern computing architectures like GPU. """ if backend is None: backend = pygmtools.BACKEND # check the correctness of input _check_data_type(A, backend) A_shape = _get_shape(A, backend) if not (len(A_shape) == 3 and A_shape[1] == A_shape[2]): raise ValueError(f"Unsupported input data shape: got A {A_shape}") num_graph, max_node = A_shape[0], A_shape[1] _check_data_type(W, backend) W_shape = _get_shape(W, backend) if not (len(W_shape) == 4 and W_shape[0] == W_shape[1] == num_graph and W_shape[2] == W_shape[3] == max_node): raise ValueError(f"Unsupported input data shape: got A {A_shape}, W {W_shape}") if ns is not None: _check_data_type(ns, backend) ns_shape = _get_shape(ns, backend) if not (len(ns_shape) == 1 and ns_shape[0] == num_graph): raise ValueError(f"The size of ns mismatches the sizes of A and W: got ns {ns_shape}, A {A_shape}, W {W_shape}") if n_univ is None: n_univ = max_node if U0 is not None: _check_data_type(U0, backend) if not sk_init_tau > 0: raise ValueError(f"sk_init_tau must be >0, got sk_init_tau={sk_init_tau}") if not sk_min_tau > 0: raise ValueError(f"sk_min_tau must be >0, got sk_min_tau={sk_min_tau}") if not 0 < sk_gamma < 1: raise ValueError(f"sk_gamma must be in (0, 1), got sk_gamma={sk_gamma}") args = (A, W, ns, n_univ, U0, sk_init_tau, sk_min_tau, sk_gamma, sk_iter, max_iter, param_lambda, verbose) try: mod = importlib.import_module(f'pygmtools.{backend}_backend') fn = mod.gamgm except ModuleNotFoundError and AttributeError: raise NotImplementedError( NOT_IMPLEMENTED_MSG.format(backend) ) return fn(*args)