Fast Approximate Quadratic Programming for Graph Matching

Joshua T. Vogelstein; John M. Conroy; Vince Lyzinski; Louis J. Podrazik; Steven G. Kratzer; Eric T. Harley; Donniell E. Fishkind; R. Jacob Vogelstein; Carey E. Priebe

doi:10.1371/journal.pone.0121002

Abstract

Quadratic assignment problems arise in a wide variety of domains, spanning operations research, graph theory, computer vision, and neuroscience, to name a few. The graph matching problem is a special case of the quadratic assignment problem, and graph matching is increasingly important as graph-valued data is becoming more prominent. With the aim of efficiently and accurately matching the large graphs common in big data, we present our graph matching algorithm, the Fast Approximate Quadratic assignment algorithm. We empirically demonstrate that our algorithm is faster and achieves a lower objective value on over 80% of the QAPLIB benchmark library, compared with the previous state-of-the-art. Applying our algorithm to our motivating example, matching C. elegans connectomes (brain-graphs), we find that it efficiently achieves performance.

Citation: Vogelstein JT, Conroy JM, Lyzinski V, Podrazik LJ, Kratzer SG, Harley ET, et al. (2015) Fast Approximate Quadratic Programming for Graph Matching. PLoS ONE 10(4): e0121002. https://doi.org/10.1371/journal.pone.0121002

Academic Editor: Mark R. Muldoon, Manchester University, UNITED KINGDOM

Received: September 9, 2014; Accepted: February 9, 2015; Published: April 17, 2015

This is an open access article, free of all copyright, and may be freely reproduced, distributed, transmitted, modified, built upon, or otherwise used by anyone for any lawful purpose. The work is made available under the Creative Commons CC0 public domain dedication

Data Availability: The C. elegans connectome data can be found at https://github.com/jovo/FastApproximateQAP/tree/master/data/elegans. The QAPLIB can be found at http://anjos.mgi.polymtl.ca/qaplib/.

Funding: This work was partially supported by the Research Program in Applied Neuroscience (JV, RJV, CEP, http://www.jhuapl.edu/ourwork/red/an/), a National Security Science and Engineering Faculty Fellowship (CEP, http://www.acq.osd.mil/rd/basic_research/program_info/nsseff.html), Johns Hopkins University Human Language Technology Center of Excellence (DEF, JV, CEP, VL, http://hltcoe.jhu.edu/), and the XDATA program of the Defense Advanced Research Projects Agency (CEP, http://www.darpa.mil/Our_Work/I2O/Programs/XDATA.aspx) administered through Air Force Research Laboratory contract FA8750-12-2-0303. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Competing interests: The authors have declared that no competing interests exist.

1 Introduction

In its most general form, the graph matching problem (GMP)—finding an alignment of the vertices of two graphs which minimizes the number of induced edge disagreements—is equivalent to a quadratic assignment problem (QAP) [1]. QAPs were first devised by Koopmans and Beckmann in 1957 to solve a ubiquitous problem in distributed resource allocation [2], and many important problems in combinatorial optimization (for example, the “traveling salesman problem,” and the GMP) have been shown to be specialized QAPs. While QAPs are known to be NP-hard in general [3], they are widely applicable and there is a large literature devoted to their approximation and formulation; see [4] for a comprehensive literature survey. In casting the GMP as a QAP, we bring to bear a host of existing optimization theoretic tools and methodologies for addressing graph matching: Frank-Wolfe [5], gradient-descent [6], etc.

Graph matching has applications in a wide variety of disciplines, spanning computer vision, image analysis, pattern recognition, and neuroscience; see [7] for a comprehensive survey of the graph matching literature. We are motivated by applications in “connectomics,” an emerging discipline within neuroscience devoted to the study of brain-graphs, where vertices represent (collections of) neurons and edges represent connections between them [8, 9]. Analyzing connectomes is an important step for many neurobiological inference tasks. For example, it is becoming increasingly popular to diagnose neurological diseases via comparing brain images [10]. To date, however, these comparisons have largely rested on anatomical (e.g., shape) comparisons, not graph (e.g., structural) comparisons. This is despite the widely held doctrine that many psychiatric disorders are fundamentally “connectopathies,” i.e. disorders of the connections of the brain [11–14]. Thus, available tests for connectopic explanations of psychiatric disorders hinge on first choosing particular graph invariants to compare across populations, rather than comparing the graphs’ structure directly. Yet, recent results suggest that the graph invariant approach to classifying is both theoretically and practically inferior to comparing whole graphs via matching [15].

Part of the reason for the lack of publications that structurally compare brain-graphs is that existing algorithms for matching very large graphs are largely ineffective, often sacrificing accuracy for efficiency. Indeed, available human connectomes have 𝒪(10⁶) vertices and 𝒪(10⁸) edges, and building them leverages state-of-the-art advances in DT-MRI imagery, big data processing and computer vision [16]. Contrast this with the fact that the human brain consists of approximately 86 billion neurons [17]. At the other end of the spectrum, the small hermaphroditic Caenorhabditis elegans worm (C. elegans) has only 302 neurons, with a fully mapped connectome. Consequently, these graphs demand GM algorithms that both accurately match small–to–moderately sized graphs and also scale to match very large graphs.

When matching these large connectome graphs (and, more broadly, the large graphs common in big data [18]), we necessarily face an accuracy/efficiency trade-off when approximately solving these GMPs: slower algorithms could achieve better performance given more time (at the extreme, exhaustive search algorithms clearly have optimal performance). Our algorithm—the Fast Approximate QAP (FAQ) algorithm for approximate GM—achieves the best available trade-off between accuracy and efficiency, outperforming the existing state-of-the-art in both accuracy and efficiency on a large proportion of QAP benchmarks and biologically inspired network matching problems.

The remainder of this paper is organized as follows. Section 2.1 formally defines the QAP and a relaxation thereof that we will operate under. Section 2.3 defines graph matching, and explains how it can be posed as a QAP.

Section 3 describes the FAQ algorithm. Section 4 provides a number of theoretical and empirical results, and compares our algorithm to previous state-of-the-art algorithms. This section concludes with an analysis of FAQ on our motivating problem of matching worm brain connectomes. We conclude with a discussion of possible extensions of FAQ and related work in Section 5.

2 Preliminaries

In this section, we formally define the QAP and the GMP. We then show how the GMP can be recast as a special case of the QAP.

2.1 Quadratic Assignment Problems

We first define the general quadratic assignment problem. Let A = (a_uv), B = (b_uv) ∈ ℝ^n×n be two n × n real matrices. Let Π: = Π_n be the set of permutation functions (bijections) of the set [n] = {1,…, n}. We define the Koopmans-Beckmann (KB) version QAP via: (1)

Note that occasionally an additional linear function is added, though we drop it here for brevity.

Eq (1) can also be recast in matrix notation. Let 𝒫 be the set of n × n permutation matrices, $𝒫 = {P \in {0, 1}^{n \times n} : P^{⊤} 1 = P 1 = 1}$ , where 1 is the n-dimensional column vector consisting of all 1’s.

Thus, Eq (1) can be written more compactly in matrix notation via: (2) We hereafter refer to (2) as the QAP optimization function.

2.2 Relaxed Quadratic Assignment Problem

Eq (2) is a binary quadratic program with linear constraints. Because of the combinatorial nature of the feasible region, finding a global optimum of (2) is NP-hard. Rather than directly optimizing over the permutation matrices, we begin by relaxing the constraint set to the convex hull of 𝒫, the set of doubly stochastic matrices (i.e. the Birkhoff polytope), where ⪰ indicates an element-wise inequality. Relaxing 𝒫 to 𝒟 in (2) yields the relaxed quadratic assignment problem (rQAP): (3) Note that, although rQAP is a quadratic program with linear constraints, it is not necessarily convex. Indeed, the objective function, $f (P) = trace (A P B^{⊤} P^{⊤})$ , has a Hessian that is not necessarily positive definite: $\nabla^{2} f (P) = B \otimes A + B^{⊤} \otimes A^{⊤}$ , where ⊗ indicates the Kronecker product (note that if A and B are hollow matrices—as is common for graphs—then ∇² f(P) has trace equal to 0, and is indefinite).

While nonconvex quadratic optimization is, in general, NP-hard, relaxing the feasible region allows us to employ the tools of continuous optimization to search for a local optima of (3). These local optima can then be projected onto 𝒫, yielding an approximate solution of (2). We also note that when relaxed to 𝒟, the QAP optimization function is often multimodal, making initialization important when solving (3).

2.3 Graph Matching

A labeled graph G = (𝒱,ℰ) consists of a vertex set 𝒱 = [n], and an edge set $ℰ \subset (\begin{matrix} 𝒱 \\ 2 \end{matrix})$ in the undirected case, or ℰ ⊂ 𝒱 × 𝒱 in the directed case. For an n-vertex graph G, we define the associated adjacency matrix A = (a_uv) ∈ {0,1}^{n × n} to be the binary n × n matrix with a_uv = 1if {u, v} ∈ ℰ in the undirected setting, or (u, v) ∈ ℰ in the directed setting. Given a pair of n-vertex graphs G_A = (𝒱_A,ℰ_A) and G_B = (𝒱_B,ℰ_B), with respective adjacency matrices A and B, we consider the following two closely related problems:

Graph Isomorphism (GI): Does there exist a P ∈ 𝒫 such that $A = P B P^{⊤}$ .
Graph Matching: $\min_{P \in 𝒫} ‖ A - P B P^{⊤} ‖_{F}$ , where ‖·‖_F is the usual matrix Frobenius norm.

GI is one of few problems with unknown computational complexity in the NP-hierarchy [19]; indeed, if P ≠ NP, then GI might reside in an intermediate complexity class called GI-complete. Moreover, GI is, at worst, only moderately exponential, with complexity 𝒪(exp{n^1/2+o(1)}) [20]. On the other hand, the (harder) GMP—reducible to a QAP—is known to be NP-hard in general. Although polynomial time algorithms are available for GM (and GI) for large classes of problems (e.g., planar graphs, trees) [21], these algorithms often have lead constants which are very large. For example, if all vertices have degree less than k, there is a linear time algorithm for GI. However, the hidden constant in this algorithm is (512k³)! [22]. Because we are interested in solving GM for graphs with

\dot{\approx} 10^{6}

or more vertices, exact GM solutions will be computationally intractable. As such, we develop a fast approximate graph matching algorithm. Our approach is based on formulating GM as a quadratic assignment problem.

2.4 Graph Matching as a Quadratic Assignment Problem

Given a pair of n-vertex graphs G_A = (𝒱_A, ℰ_A) and G_B = (𝒱_B, ℰ_B), with respective adjacency matrices A and B, we can formally write the graph matching problem as an optimization problem: (4) Simple algebra yields that, (5) The GMP is then equivalent (i.e. same argmin) to (6) The objective function for GM is just the negative of the objective function for QAP. Thus, any descent algorithm for the former can be directly applied to the latter. Moreover, any QAP approximation algorithms also immediately yields an analogous GM approximation.

As is common in solving general QAPs, GM algorithms often begin by first relaxing (4) to a continuous problem (see, for example, [23]). The resulting problem is a convex quadratic program and can be efficiently exactly solved. The obtained solution is then projected back onto 𝒫 yielding an approximate solution to (4). Contrary to popular existing approaches, our FAQ algorithm first solves a relaxed version of (6) and subsequently projects the solution back onto 𝒫. This relaxation yields an indefinite quadratic program, and indefinite quadratic programs are are NP-hard to solve in general. However, recent theory indicates that the indefinite relaxation of (6), and not the convex relaxation of (4) is the provably correct approach [24]. Reflecting this theory, we find that FAQ obtains state-of-the-art performance in terms of both computational efficiency and objective function value for various QAPs (see Section 4).

3 Fast Approximate QAP Algorithm

Our algorithm, called FAQ, proceeds in three steps:

Choose a suitable initial position
Find a local solution to rQAP.
Project back onto the set of permutation matrices.

These steps are summarized in Algorithm 1. Below, we provide details for each step.

Algorithm 1 FAQ for finding a local optimum of rQAP

Require: Graphs (adjacency matrices) A and B as well as a stopping criteria

Ensure: $\hat{P}$ , an estimated permutation matrix

1: Choose an initialization, $P^{(0)} = 1 1^{⊤} / n$

2: while stopping criteria not met do

3: Compute the gradient of f at the current point via Eq (3)

4: Compute the direction Q⁽ⁱ⁾ by solving Eq (8) via the Hungarian algorithm

5: Compute the step size α⁽ⁱ⁾ by solving Eq (9)

6: Update P⁽ⁱ⁾ according to Eq (10)

7: end while

8: Obtain $\hat{P}$ by solving Eq (11) via the Hungarian algorithm.

Find a suitable initial position. While any doubly stochastic matrix would be a feasible initial point, we choose the noninformative “flat doubly stochastic matrix,” $J = 1 \cdot 1^{⊤} / n$ , i.e. the barycenter of the feasible region. Alternately, we could use multiple restarts, each initial point near J. Specifically, we could sample K, a random doubly stochastic matrix using 10 iterations of Sinkhorn balancing [25], and let P⁽⁰⁾ = (J+K)/2. Given this initial estimate (or estimates), we would then iterate the following five steps until convergence.
Find a local solution to rQAP. As mentioned above, rQAP is a quadratic problem with linear constraints. A number of off-the-shelf algorithms are readily available for finding local optima in such problems. We utilize the Frank-Wolfe algorithm (FW), a successive first-order optimization procedure originally devised to solve convex quadratic programs [5, 26]. Although FW is a relatively standard solver, especially as a subroutine for QAP algorithms [27], we provide a detailed view of applying FW to rQAP.
Given an initial position, P⁽⁰⁾, iterate the following four steps:
Step 1, Compute the gradient ∇f(P⁽ⁱ⁾): The gradient of $f (P) = - trace (A P B^{⊤} P^{⊤})$ with respect to P, evaluated at P⁽ⁱ⁾, is given by $\nabla f (P^{(i)}) = - A P^{(i)} B^{⊤} - A^{⊤} P^{(i)} B$ .
Step 2, Compute the search direction Q⁽ⁱ⁾: The search direction is given by the argument that minimizes a first-order Taylor series approximation to f(P) around the current estimate, P⁽ⁱ⁾: (7) Dropping terms independent of P, we obtain the following sub-problem: (8) which is equivalent to a Linear Assignment Problem (LAP), and can therefore be solved in O(n³) time via the “Hungarian Algorithm” of [28, 29]. Let Q⁽ⁱ⁾ indicate the argmin of Eq (8).
Step 3, Compute the step size α⁽ⁱ⁾: Given Q⁽ⁱ⁾, we minimize the original optimization problem, along the line segment from P⁽ⁱ⁾ to Q⁽ⁱ⁾: (9) This can be solved exactly, as f is a quadratic function of α. Let α⁽ⁱ⁾ indicate the argmin of Eq (9).
Step 4, Update P⁽ⁱ⁾: Finally, the next iterate is the doubly stochastic matrix (10)
Stopping criteria: Steps 1–4 are iterated until some stopping criterion is met. Often, a threshold ϵ > 0 or an iteration limit is given, and the algorithm iterates until the iteration limit is reached, ‖P⁽ⁱ⁾−P⁽ⁱ⁻¹⁾‖_F < ϵ, or ‖∇f(P⁽ⁱ⁾)‖_F < ϵ. In practice, the algorithm often converges with a modest number of FW iterates.
Project onto the set of permutation matrices. Let P^(final) be the doubly stochastic matrix resulting from the final iteration of FW. We project P^(final) onto the set of permutation matrices via (11) Note that Eq (11) is again equivalent to a LAP and can be solved in O(n³) time.

4 Results

Below we present a number of empirical and theoretical results demonstrating the state-of-the-art efficiency and accuracy of the FAQ algorithm.

4.1 Algorithm Complexity and leading constants

As mentioned above, GM is computationally difficult, and even in the special cases for which polynomial time algorithms are available, the leading constants are intractably large. Given a bounded number of FW iterates, the FAQ algorithm has complexity O(n³). However, a very large lead order constant could still render this algorithm practically infeasible. Fig 1 suggests that FAQ has O(n³) complexity, and also has very small leading constants (≈ 10⁻⁹). This suggests that this algorithm is feasible for matching even reasonably large graphs. Note that other state-of-the-art approximate graph matching algorithms also have cubic or worse time complexity in the number of vertices. We will describe these other algorithms and their time complexity in greater detail below.

Download:

Fig 1. Running time of FAQ as function of number of vertices.

Data was sampled from an Erdös-Rényi model with p = log(n)/n. Each dot represents a single simulation, with 100 simulations per n. The solid line is the best fit cubic function. Note the leading constant is ≈ 10⁻⁹. FAQ finds the optimal objective function value in every simulation.

https://doi.org/10.1371/journal.pone.0121002.g001

4.2 QAP Benchmark Accuracy

Having demonstrated both theoretically and empirically that FAQ has cubic time complexity, we evaluate the algorithms accuracy on a suite of standard benchmarks. More specifically, QAPLIB is a library of 137 quadratic assignment problems, ranging in size from 10 to 256 vertices [30]. Recent graph matching papers typically evaluate the performance of their algorithm on 16 of the benchmark QAPs that are known to be particularly difficult [23, 31]. We compare the results of FAQ to the results of four other state-of-the-art graph matching algorithms:

the PATH algorithm, which solves the relaxation of (4), and then solves a sequence of concave and convex problems to project the solution onto 𝒫 [23];
QCV which solves the relaxation of (4), and projects the obtained solution onto the closest permutation in 𝒫;
the RANK algorithm [32], a spectral graph matching procedure;
the Umeyama algorithm (denoted by U), another spectral graph matching procedure [1].

Note that the code for these four algorithms is freely available from the graphm package [23].

Performance is measured by minimizing the assignment cost $f (P) = trace (A P B^{⊤} P^{⊤})$ . We write ${\hat{f}}_{X}$ for the value of the local minimum of f obtained by the algorithm X ∈ {FAQ, PATH, QCV, RANK, U, all}, where “all” is just the best performer of all the non FAQ algorithms. Fig 2 plots the logarithm (base 10, here and elsewhere) of the relative accuracy, i.e. $\log_{10} ({\hat{f}}_{F AQ} / {\hat{f}}_{X})$ , for X ∈ {PATH, QCV, RANK, U, all}. FAQ performs significantly better than the other algorithms, outperforming all of them on ≈ 94% of the problems, often by nearly an order of magnitude in terms of relative error.

Download:

Fig 2. Relative accuracy—defined to be

\log_{10} ({\hat{f}}_{F A Q} / {\hat{f}}_{X})

—of all the four algorithms compared with FAQ.

Note that FAQ is better than all the other algorithms on ≈ 94% of the benchmarks. The abscissa is the log number of vertices. The gray dot indicates the mean improvement of FAQ over the other algorithms.

https://doi.org/10.1371/journal.pone.0121002.g002

4.3 QAP Benchmark Efficiency

The utility of an approximation algorithm depends not just on its accuracy, but also its efficiency. To empirically test these algorithms’ efficiency, we compare the wall time of each of the five algorithms on all 137 QAPS in QAPLIB in Fig 3. For each of the 5 algorithms, we fit an iteratively weighted least squares linear regression function (Matlab’s “robustfit”) to regress the logarithm of time (in seconds) onto the logarithm of the number of vertices being matched. The numbers beside the lines indicate the slopes of the regression functions.

Download:

Fig 3. Absolute wall time for running each of the five algorithms on all 137 QAPLIB benchmarks.

We fit a line on this log-log plot for each algorithm; the slope is displayed beside each line. The FAQ slope is much better than the PATH slope, and worse than the others. Note, however, the time for RANK and U appears to be superlinear on this log-log plot, suggesting that perhaps as the number of vertices increases, PATH might be faster.

https://doi.org/10.1371/journal.pone.0121002.g003

The figure demonstrates that the PATH algorithm has the largest slope, significantly larger than that of FAQ. QCV and FAQ have nearly identical slopes, which makes sense, given that the are solving very similar objective functions. Similarly, RANK and U have very similar slopes; they are both using spectral approaches. Note, however, that although the slope of RANK and U are smaller than that of FAQ, they both appear to be super linear in this log-log plot, suggesting that as the number of vertices increases, their compute time might exceed that of the other algorithms.

Combined with Fig 2, this figure suggests that FAQ achieves the state-of-the-art trade-off between efficiency and accuracy. Of note is that the FAQ algorithm has a relatively high variance in wall time for these problems. This is due to the number of Hungarian algorithms performed in the FW subroutine. We could fix the number of Hungarian algorithms, in which case the variance would decrease dramatically. However, in application, this variance is not problematic.

4.4 QAP Benchmark Accuracy/Efficiency Trade-off

In [23], the authors demonstrated that PATH outperformed both QCV and U on a variety of simulated and real examples. Fig 4 compares the performance of FAQ with PATH along both dimensions of performance—accuracy and efficiency—for all 137 benchmarks in the QAPLIB library. The right panel indicates that FAQ is both more accurate and more efficient on 80% of the problems (and is more accurate on 99% of the benchmarks). The middle plots the relative wall time of FAQ to PATH as a function of the number of vertices, also on a log-log scale. The gray line is the best fit slope on this plot, suggesting that FAQ is getting exponentially faster than PATH as the number of vertices gets larger.

Download:

Fig 4. Comparison of FAQ with PATH in terms of both accuracy and efficiency.

The left panel is the same as the left panel of Fig 2. The middle plots the relative wall time of FAQ to PATH as a function of the number of vertices, also on a log-log scale. The gray line is the best fit slope on this plot. Finally, the right panel plots log relative time versus log relative objective function value, demonstrating that FAQ outperforms PATH on both dimensions on 80% of the benchmarks.

https://doi.org/10.1371/journal.pone.0121002.g004

4.5 QAP Directed Benchmarks

Recently, Liu et al. [33] proposed a modification of the PATH algorithm that adjusted PATH to be more appropriate for directed graphs Note that our FAQ algorithm easily extend to directed or weighted graphs. Liu et al. compare the performance of their algorithm (EPATH) with U, QCV, and the GRAD algorithm of [34] on a set of 16 particularly difficult directed benchmarks from QAPLIB. The EPATH algorithm achieves at least as low objective value as the other algorithms on 15 of 16 benchmarks. Our algorithm, FAQ, performs at least as well as EPATH, U, QCV, and GRAD on all 16 benchmarks, and achieves the singular best performance on 8 of the benchmarks. Table 1 shows the numerical results comparing FAQ to EPATH and GRAD, the only algorithm considered in [33] to outperform EPATH. Note that some of the algorithms achieve the absolute minimum on 8 of the benchmark.

Download:

Table 1. Comparison of FAQ with optimal objective function value and previous state-of-the-art for directed graphs.

The best (lowest) value is in bold. Asterisks indicate achievement of the global minimum. The number of vertices for each problem is the number in its name (second column).

https://doi.org/10.1371/journal.pone.0121002.t001

4.6 Theoretical properties of FAQ

In addition to guarantees on computational time, we have a guarantee on performance:

Proposition 1 If A is the adjacency matrix of an asymmetric simple graph G, then Proof Let m denote the number of edges in G. As G is asymmetric, for any P ≠ I. By the Birkhoff-von Neuman Theorem, 𝒟 is the convex hull of 𝒫, i.e., for all D ∈ 𝒟, there exists constants {α_{D, P}}_{P ∈ 𝒫} such that D = ∑_{P ∈ 𝒫} α_{D, P} P and ∑_{P ∈ 𝒫} α_{D, P} = 1. Thus, if D is not the identity matrix, as D ≠ I implies that α_{D, P} > 0 for some P ≠ I. ■

Remark 1 Note that it trivially follows from Proposition 1 that if A and B are the adjacency matrices of asymmetric isomorphic simple graphs, then the minimum objective function value of rGMP is equal to the minimum objective function value of GMP.

Remark 2 For the convex quadratic GM relaxation, Eq (4), in general ${argmin}_{D \in 𝒟} ‖ A D - D A ‖_{F}^{2} \neq {I},$ even if G (the graph with adjacency matrix A) is asymmetric. Indeed, degree regular graphs provide a simple counterexample, as in this case $J \in {argmin}_{D \in 𝒟} ‖ A D - D A ‖_{F}^{2}$ . We will empirically explore the ramifications of this phenomena further in Section 4.8

This result says that, when solving the GI problem, nothing is lost by relaxing the indefinite GM problem as done by FAQ. Recent results also show this is almost surely the case (in a broad class of random graphs) when relaxing the indefinite GM problem, even in the non-isomorphic graph setting [24] (and is again almost surely not the case when relaxing the convex GM formulation). Combined, this serve to provide theoretical justification for our FAQ procedure.

4.7 Multiple Restarts

Due to the indefiniteness of the relaxation of (6), the feasible region may be rife with local minima. As a result, our FAQ procedure is sensitive to the chosen initial position. With this in mind, we propose a variant of the FAQ algorithm in which we run the FAQ procedure with multiple initializations. The algorithm outputs the best FAQ iterate over all the initializations. For each initialization, we sample K ∈ 𝒟, a random doubly stochastic matrix, using 10 iterations of Sinkhorn balancing [25], and let our initialization be P⁽⁰⁾ = (J+K)/2, where J is the doubly stochastic barycenter. Fixing the number of restarts, this variant of FAQ still has O(n³) complexity.

Table 2 shows the performance of running FAQ 3 and 100 times, reporting only the best result (indicated by FAQ₃ and FAQ₁₀₀, respectively), and comparing it to the best performing result of the five algorithms (including the original FAQ)—note that the best performing of the original five tested was always FAQ. Note that we only consider the 16 particularly difficult benchmarks for this evaluation. FAQ only required three restarts to outperform all other approximate algorithms on all 16 of 16 difficult benchmarks. Moreover, after 100 restarts, FAQ finds the absolute minimum on 3 of the 16 benchmarks. Note that no other algorithm ever achieved the absolute minimum on any of these benchmarks.

Download:

Table 2. Comparison of FAQ with optimal objective function value and the best result on the undirected benchmarks.

Note that FAQ restarted 100 times finds the optimal objective function value in 3 of 16 benchmarks, and that FAQ restarted 3 times finds a minimum better than the previous state-of-the-art on all 16 particularly difficult benchmarks.

https://doi.org/10.1371/journal.pone.0121002.t002

Remark 3 Another natural starting position for FAQ is the doubly stochastic solution of rGMP, the relaxed Eq (4). Promising results in this direction are pursued further in [24].

4.8 Brain-Graph Matching

The Caenorhabditis elegans (C. elegans) is a small worm (nematode) with 302 labeled neurons (in the hermaphroditic sex). The chemical connectome of C. elegans is a weighted directed graph on 279 vertices, with edge weights in {0,1,2,…} counting the number of directed chemical synapses between the neurons [35, 36]. We conduct the following synthetic experiments. For A = (A_uv) count the number of synapses (in {0,1,2,…}) from neuron u to neuron v in the C. elegans connectome, for all u, v. For k = 1,2,…,1000, we choose Q_(k) uniformly at random from 𝒫, and let $B_{(k)} = Q_{(k)} A {Q_{(k)}}^{⊤}$ , effectively shuffling the vertex labels of the connectome. Then, we match graphs A to B_(k). We define accuracy as the fraction of vertices correctly assigned (i.e. unshuffled).

Fig 5 displays the results of FAQ (initialized at J) along with U, QCV, and PATH. The left panel indicates that FAQ always perfectly unshuffles the chemical connectome, whereas none of the other algorithms ever perfectly unshuffles the graph. In light of Proposition 1, this is unsurprising. Indeed, there is a unique automorphism for this connectome, and the graph is asymmetric. For any choice of Q_(k), the indefinite problem (6) therefore has a unique solution, namely $Q_{(k)}^{⊤}$ . FAQ finds this global minima in all the cases. Contrast this with Eq (4)—the objective function of PATH and QCV—which could have multiple global minima in 𝒟. This could account for the high variance in the performance of QCV and PATH in Fig 5.

Download:

Fig 5. Performance of U, QCV, PATH, and FAQ on synthetic C. elegans connectome data, graph matching the true chemical connectome with permuted versions of itself.

Error is the fraction of vertices incorrectly matched. Circle indicates the median, thick black bars indicate the quartiles, thin black lines indicate extreme but non-outlier points, and plus signs are outliers. The left panel indicate error (fraction of misassigned vertices), and the right panel indicates wall time on a 2.2 GHz Apple MacBook. FAQ always obtained the optimal solution, whereas none of the other algorithms ever found the optimal. FAQ also ran very quickly, nearly as quickly as U and QCV, and much faster than PATH, even though the FAQ implementation is in Matlab, and the others are in C.

https://doi.org/10.1371/journal.pone.0121002.g005

The right panel compares the wall time of the various algorithms, running on an 2.2 GHz Apple MacBook. Note that we only have a Matlab implementation of FAQ, whereas the other algorithms are implemented in C. Unlike in the QAPLIB benchmarks, FAQ runs nearly as quickly as both U and QCV; and as expected, FAQ runs significantly faster than PATH. We also ran FAQ on a binarized symmeterized versions of the chemical connectome graph $\tilde{A}$ (i.e. ${\tilde{A}}_{u v} = 1$ if and only if A_uv ≥ 1 or A_vu ≥ 1). The resulting errors are nearly identical to those presented in Fig 5, although speed increased for FAQ by more than a factor of two. Note that the number of vertices in this brain-graph matching problem—279—is larger than the largest of the 137 benchmarks used above.

5 Discussion

This work presents the FAQ algorithm, a fast algorithm for approximately matching very large graphs. Our key insight was to continuously relax the indefinite formulation of the GMP in the FW implementation. We demonstrated that not only is FAQ cubic in time, but also its leading constant is quite small—10⁻⁹—suggesting that it can be used for graphs with hundreds or thousands of vertices (§4.1).

Moreover, FAQ achieves better accuracy than previous state-of-the-art approximate algorithms on on over 93% of the 137 QAPLIB benchmarks (§4.2), is faster than the state-of-the-art PATH algorithm (§4.3), and is both faster and achieves at least as low performance as PATH on over 80% of the tested benchmarks (§4.4), including both directed and undirected graph matching problems (§4.5). In addition to the theoretical guarantees of cubic run time, we provide theoretical justification for optimizing the indefinite GM formulation, Eq (6) as opposed to Eq (4) (§4.6). Indeed, the indefinite formulation (and not the convex formulation) has the property that when matching asymmetric isomorphic graphs, the unique global minimum of the indefinite relaxation is the isomorphism between the two graphs.

Because rQAP is indefinite, we also considered FAQ with multiple restarts, and achieve improved performance for the particularly difficult benchmarks using only three restarts (§4.7). Finally, we used FAQ to match permuted versions of the C. elegans connectomes (§4.8). Of the four state-of-the-art algorithms considered, FAQ achieved perfect performance 100% of the time, whereas none of the other three algorithms ever matched the connectomes perfectly. Moreover, FAQ ran comparably fast to U and QCV and significantly faster than PATH, even though FAQ is implemented in Matlab, and the others are implemented in C. Note that these connectomes have 279 vertices, more vertices than the largest QAP benchmarks.

5.1 Related Work

Our approach is quite similar to other approaches that have recently appeared in the literature. Perhaps its closest cousins include [23, 37–39], which are all of the “PATH” following or FW varieties. These algorithms begin by relaxing the convex objective function in (4), while FAQ begins by relaxing the indefinite objective function in (6). Although the convex relaxation is efficiently solvable, the obtained solution is almost surely incorrect (for a broad class of random graphs) and the correct solution is often not obtained even post projection [24]. The indefinite relaxation however, almost surely yields the correct solution when exactly solved (for a broad class of random graphs) [24]. With this in mind, it is unsurprising that FAQ outperforms PATH on nearly all benchmark problems. Others have considered similar relaxations to PATH, but usually in the context of finding lower bounds [40] or as subroutines for finding exact solutions [41]. Our work seems to be the first to utilize the precise algorithm described in Algorithm 1 to find fast approximate solutions to GMP and QAP.

We note that the work in [38] is a generalization of our FW approach. The authors of [38] apply FW with exact line search to find (12) where M is an n² × n² matrix whose M_{i, j;k,ℓ} entry measures the similarity between the vertex pairs (i, k) and (j,ℓ), and C is a constraint matrix enforcing the one–to–one (easily generalizing to many–to–one) matching of vertices. Our GM formulation (6), and our subsequent FW implementation, can be realized from (12) by setting M = −B^T ⊗ A^T, and having C enforce a one–to–one matching of the vertices. Our approach offers several distinct advantages over the general approach of [38], namely:

The ability to efficiently employ multiple random restarts for enhanced performance (see Section 4.7 for detail);
Our FW implementation exploits the Kronecker product structure of M which reduces the runtime of an O(n⁴) procedure to O(n³);
Our final step is projected to the nearest permutation matrix, which is not part of the algorithm outlined in [38];
In [38], they consider general M, and in so doing, they are considering all quadratic-objective optimizations (albeit without a linear and constant term). There are many quadratic formulations of GM, and our contribution is touting the efficiency and efficacy of solving our particular indefinite form.

5.2 Future Work

Even with the very small lead order constant for FAQ, as n increases, the computational burden gets quite high. For example, extrapolating the curve of Fig 1, this algorithm would take about 20 years to finish (on a standard laptop) when n = 100,000. We hope to be able to approximately solve rQAP on graphs much larger than that, given that the number of neurons even in a fly brain, for example, is ≈ 250,000. More efficient algorithms and/or implementations are required for such massive graph matching. Although a few other state-of-the-art algorithms were more efficient than FAQ, their accuracy was significantly worse. We are actively working on combining FAQ with dimensionality reduction procedures to achieve the desired scaling from FAQ [42].

We are also pursuing additional future work to generalize FAQ in a number of ways:

In addition to the theoretical results contained in Section 4.6, we have studied the properties of the convex and indefinite GMP relaxations in a very general random graph model [24]. Under some general assumptions on the random graph model, the indefinite relaxation of (6), and not the convex relaxation of (4) is the provably correct approach [24].
Many (brain-) graphs of interest will be errorfully observed [43], that is, vertices might be missing and putative edges might exhibit both false positives and negatives. Explicitly dealing with this error source is both theoretically and practically of interest [15].
For many brain-graph matching problems, the number of vertices will not be the same across the brains. Recent work from [23, 37] and [39] suggest that extensions in this direction would be both relatively straightforward and effective.
Often, a partial matching of the vertices is known a priori, and we can modify FAQ to leverage these seeded vertices to match the remaining unseeded vertices [44].
The most costly subroutine in FAQ is solving LAPs. Fortunately, the LAP is a linear optimization problem with linear constraints. As a result, a number of parallelized optimization strategies could be implemented on this problem [45].
Often, real data adjacency matrices have certain special properties, namely sparsity, which makes faster LAP subroutines [29] and more efficient algorithms (such as “active set” algorithms) readily available for further speed increases.
In many graph settings, we have some prior information that could easily be incorporated into the GM problem in the form of vertex attributes and features. For example, in brain graphs we know the position of the vertex in the brain, the vertex’s cell type, etc. These could be used to measure a “dissimilarity” between vertices and are easily incorporated into FAQ’s objective function to better match the graphs.
We are working to extend FAQ to match multiple graphs simultaneously.

5.3 Concluding Thoughts

In conclusion, this manuscript has presented the FAQ algorithm for approximately solving the quadratic assignment problem. FAQ is theoretically justified, fast, effective, and easily generalizable. Our algorithm achieves state-of-the-art matching performance and efficiency on a host of benchmark QAP problems and connectome data sets. Yet, the 𝒪(n³) complexity remains too slow to solve many problems of interest, and issues of scalability need be addressed. To facilitate further development and applications, all the code and data used in this manuscript is available from the first author’s website, http://jovo.me. We have further incorporated FAQ (as sgm.R) into the open-source R package, igraph, available for download at https://github.com/igraph/xdata-igraph/. MATLAB code is also available at https://github.com/jovo/FastApproximateQAP/tree/master/code/FAQ.

Acknowledgments

The authors would like to acknowledge Lav Varshney for providing the data.

Author Contributions

Conceived and designed the experiments: JV JC RJV SK LP CEP DEF EH. Performed the experiments: JC EH DEF JV SK LP. Analyzed the data: CEP JV JC DEF. Contributed reagents/materials/analysis tools: VL JV SK LP DEF EH. Wrote the paper: JV DEF VL.

References

1. Umeyama S (1988) An eigendecomposition approach to weighted graph matching problems. Pattern Analysis and Machine Intelligence, IEEE Transactions on 10: 695–703.
- View Article
- Google Scholar
2. Koopmans TC, Beckman M (1957) Assignment problems and the location of economic activities. The Econometric Society 25: 53–76.
- View Article
- Google Scholar
3. Papadimitriou CH, Steiglitz K (1998) Combinatorial Optimization: Algorithms and Complexity. Dover Publications, 496 pp.
4. Burkard RE (2013) The quadratic assignment problem. In: Handbook of Combinatorial Optimization, Springer. pp. 2741–2814.
5. Frank M, Wolfe P (1956) An algorithm for quadratic programming. Naval Research Logistics Quarterly 3: 95–110.
- View Article
- Google Scholar
6. Bazaraa MS, Sherali HD, Shetty CM (2013) Nonlinear programming: theory and algorithms. John Wiley & Sons.
7. Conte D, Foggia P, Sansone C, Vento M (2004) Thirty years of graph matching in pattern recognition. International Journal of Pattern Recognition and Artificial Intelligence 18: 265–298.
- View Article
- Google Scholar
8. Sporns O, Tononi G, Kotter R (2005) The human connectome: A structural description of the human brain. PLoS Computational Biology 1: e42. pmid:16201007
- View Article
- PubMed/NCBI
- Google Scholar
9. Hagmann P (2005) From diffusion MRI to brain connectomics. Ph.D. thesis, Institut de traitement des signaux.
10. Csernansky JG, Wang L, Joshi SC, Ratnanather JT, Miller MI (2004) Computational anatomy and neuropsychiatric disease: probabilistic assessment of variation and statistical inference of group difference, hemispheric asymmetry, and time-dependent change. NeuroImage 23 Suppl 1: S56–68. pmid:15501101
- View Article
- PubMed/NCBI
- Google Scholar
11. Kubicki M, McCarley R, Westin CF, Park HJ, Maier S, Kikinis R, et al. (2007) A review of diffusion tensor imaging studies in schizophrenia. Journal of Psychiatric Research 41: 15–30. pmid:16023676
- View Article
- PubMed/NCBI
- Google Scholar
12. Calhoun VD, Sui J, Kiehl KA, Turner JA, Allen E, Pearlson G (2011) Exploring the psychosis functional connectome: Aberrant intrinsic networks in schizophrenia and bipolar disorder. Frontiers in Psychiatry / Frontiers Research Foundation 2: 75.
- View Article
- Google Scholar
13. Fornito A, Bullmore ET (2012) Connectomic intermediate phenotypes for psychiatric disorders. Frontiers in Psychiatry / Frontiers Research Foundation 3: 32.
- View Article
- Google Scholar
14. Fornito A, Zalesky A, Pantelis C, Bullmore ET (2012) Schizophrenia, neuroimaging and connectomics. Neuroimage 62: 2296–2314. pmid:22387165
- View Article
- PubMed/NCBI
- Google Scholar
15. Vogelstein JT, Priebe CE (2014) Shuffled graph classification: Theory and connectome applications. Journal of Classification, In press.
- View Article
- Google Scholar
16. Gray Roncal W, Koterba ZH, Mhembere D, Kleissas DM, Vogelstein JT, Burns, R, et al. (2013) MIGRAINE: MRI graph reliability analysis and inference for connectomics. Global Conference on Signal and Information Processing.
17. Herculano-Houzel S (2012) The remarkable, yet not extraordinary, human brain as a scaled-up primate brain and its associated cost. Proceedings of the National Academy of Sciences of the United States of America 109 Suppl: 10661–8. pmid:22723358
- View Article
- PubMed/NCBI
- Google Scholar
18. Kolaczyk ED (2009) Statistical analysis of network data. Springer.
19. Fortin S (1996) The graph isomorphism problem. Technical Report, University of Alberta, Department of Computer Science.
20. Babai L (1981) Moderately exponential bound for graph isomorphism: 34–50.
21. Babali L, Erds P, Selkow SM (1980) Random graph isomorphism. SIAM Journal on Computing 9: 628–635.
- View Article
- Google Scholar
22. Chen J (1994) A linear-time algorithm for isomorphism of graphs of bounded average genus. SIAM Journal on Discrete Mathematics 7: 614–631.
- View Article
- Google Scholar
23. Zaslavskiy M, Bach FR, Vert JPP (2009) A path following algorithm for the graph matching problem. IEEE Transactions on Pattern Analysis and Machine Intelligence 31: 2227–2242. pmid:19834143
- View Article
- PubMed/NCBI
- Google Scholar
24. Lyzinski V, Fishkind DE, Fiori M, Vogelstein JT, Priebe CE, Sapiro, G (2014) Graph matching: Relax at your own risk. arXiv preprint, arXiv: 14053133.
25. Sinkhorn R (1964) A relationship between arbitrary positive matrices and doubly stochastic matrices. The Annals of Mathematical Statistics 35: 876–879.
- View Article
- Google Scholar
26. Bradley S, Hax A, Magnanti T (1977) Applied mathematical programming. Addison Wesley.
27. Anstreicher KM (2003) Recent advances in the solution of quadratic assignment. SIAM Journal on Optimization 97: 27–42.
- View Article
- Google Scholar
28. Kuhn HW (1955) The Hungarian method for the assignment problem. Naval Research Logistics Quarterly 2: 83–97.
- View Article
- Google Scholar
29. Jonker R, Volgenant A (1987) A shortest augmenting path algorithm for dense and sparse linear assignment problems. Computing 38: 325–340.
- View Article
- Google Scholar
30. Burkard RE, Karisch SE, Rendl F (1997) QAPLIB: A quadratic assignment problem library. Journal of Global Optimization 10: 391–403.
- View Article
- Google Scholar
31. Schellewald C, Roth S, Schnörr C (2001) Evaluation of convex optimization techniques for the weighted graph-matching problem in computer vision. In: Proceedings of the 23rd DAGMSymposium on Pattern Recognition. Springer-Verlag, pp. 361–368.
32. Singh R, Xu J, Berger B (2007) Pairwise global alignment of protein interaction networks by matching neighborhood topology: 16–31.
33. Liu ZY, Qiao H, Xu L (2012) An extended path following algorithm for graph matching problem. IEEE transactions on pattern analysis and machine intelligence 34: 1451–1456.
- View Article
- Google Scholar
34. Gold S, Rangarajan A (1996) A graduated assignment algorithm for graph matching. IEEE Transactions on Pattern Analysis and Machine Intelligence 18: 377–388.
- View Article
- Google Scholar
35. White JG, Southgate E, Thomson JN, Brenner S (1986) The structure of the nervous system of the nematode Caenorhabditis elegans. Philosophical Transactions of Royal Society London Series B, Biological Sciences 314: 1–340.
- View Article
- Google Scholar
36. Varshney LR, Chen BL, Paniagua E, Hall DH, Chklovskii DB (2011) Structural properties of the Caenorhabditis elegans neuronal network. PLoS Computational Biology 7: 1–41.
- View Article
- Google Scholar
37. Zaslavskiy M, Bach FR, Vert JPP (2010) Many-to-many graph matching: A continuous relaxation approach. Machine Learning and Knowledge Discovery in Databases 6323: 515–530.
- View Article
- Google Scholar
38. Leordeanu M, Hebert M, Sukthankar R (2009) An integer projected fixed point method for graph matching and map inference. Advances in Neural Information Processing Systems: 1114–1122,
- View Article
- Google Scholar
39. Escolano F, Hancock ER, Lozano M (2011) Graph matching through entropic manifold alignment. Computer Vision and Pattern Recognition.
- View Article
- Google Scholar
40. Anstreicher KM, Brixius NW (2001) A new bound for the quadratic assignment problem based on convex quadratic programming. Mathematical Programming 89: 341–357.
- View Article
- Google Scholar
41. Brixius NW (2000) Solving quadratic assignment problems using convex quadratic programming relaxations 1: Introduction. Management: 1–20.
- View Article
- Google Scholar
42. Lyzinski V, Sussman DL, Fishkind DE, Pao H, Priebe CE (2013) Seeded graph matching for large stochastic block model graphs. arXiv preprint: 1310.1297.
43. Priebe CE, Vogelstein J, Bock D (2013) Optimizing the quantity/quality trade-off in connectome inference. Communications in Statistics-Theory and Methods 42: 3455–3462.
- View Article
- Google Scholar
44. Lyzinski V, Fishkind D, Priebe C (2014) Seeded graph matching for correlated Erdos-Renyi graphs. Journal of Machine Learning Research To appear.
- View Article
- Google Scholar
45. Boyd S, Vandenberghe L (2004) Convex Optimization. Oxford University Press.

[ref1] 1. Umeyama S (1988) An eigendecomposition approach to weighted graph matching problems. Pattern Analysis and Machine Intelligence, IEEE Transactions on 10: 695–703.
View Article
Google Scholar

[2] View Article

[3] Google Scholar

[ref2] 2. Koopmans TC, Beckman M (1957) Assignment problems and the location of economic activities. The Econometric Society 25: 53–76.
View Article
Google Scholar

[5] View Article

[6] Google Scholar

[ref3] 3. Papadimitriou CH, Steiglitz K (1998) Combinatorial Optimization: Algorithms and Complexity. Dover Publications, 496 pp.

[ref4] 4. Burkard RE (2013) The quadratic assignment problem. In: Handbook of Combinatorial Optimization, Springer. pp. 2741–2814.

[ref5] 5. Frank M, Wolfe P (1956) An algorithm for quadratic programming. Naval Research Logistics Quarterly 3: 95–110.
View Article
Google Scholar

[10] View Article

[11] Google Scholar

[ref6] 6. Bazaraa MS, Sherali HD, Shetty CM (2013) Nonlinear programming: theory and algorithms. John Wiley & Sons.

[ref7] 7. Conte D, Foggia P, Sansone C, Vento M (2004) Thirty years of graph matching in pattern recognition. International Journal of Pattern Recognition and Artificial Intelligence 18: 265–298.
View Article
Google Scholar

[14] View Article

[15] Google Scholar

[ref8] 8. Sporns O, Tononi G, Kotter R (2005) The human connectome: A structural description of the human brain. PLoS Computational Biology 1: e42. pmid:16201007
View Article
PubMed/NCBI
Google Scholar

[17] View Article

[18] PubMed/NCBI

[19] Google Scholar

[ref9] 9. Hagmann P (2005) From diffusion MRI to brain connectomics. Ph.D. thesis, Institut de traitement des signaux.

[ref10] 10. Csernansky JG, Wang L, Joshi SC, Ratnanather JT, Miller MI (2004) Computational anatomy and neuropsychiatric disease: probabilistic assessment of variation and statistical inference of group difference, hemispheric asymmetry, and time-dependent change. NeuroImage 23 Suppl 1: S56–68. pmid:15501101
View Article
PubMed/NCBI
Google Scholar

[22] View Article

[23] PubMed/NCBI

[24] Google Scholar

[ref11] 11. Kubicki M, McCarley R, Westin CF, Park HJ, Maier S, Kikinis R, et al. (2007) A review of diffusion tensor imaging studies in schizophrenia. Journal of Psychiatric Research 41: 15–30. pmid:16023676
View Article
PubMed/NCBI
Google Scholar

[26] View Article

[27] PubMed/NCBI

[28] Google Scholar

[ref12] 12. Calhoun VD, Sui J, Kiehl KA, Turner JA, Allen E, Pearlson G (2011) Exploring the psychosis functional connectome: Aberrant intrinsic networks in schizophrenia and bipolar disorder. Frontiers in Psychiatry / Frontiers Research Foundation 2: 75.
View Article
Google Scholar

[30] View Article

[31] Google Scholar

[ref13] 13. Fornito A, Bullmore ET (2012) Connectomic intermediate phenotypes for psychiatric disorders. Frontiers in Psychiatry / Frontiers Research Foundation 3: 32.
View Article
Google Scholar

[33] View Article

[34] Google Scholar

[ref14] 14. Fornito A, Zalesky A, Pantelis C, Bullmore ET (2012) Schizophrenia, neuroimaging and connectomics. Neuroimage 62: 2296–2314. pmid:22387165
View Article
PubMed/NCBI
Google Scholar

[36] View Article

[37] PubMed/NCBI

[38] Google Scholar

[ref15] 15. Vogelstein JT, Priebe CE (2014) Shuffled graph classification: Theory and connectome applications. Journal of Classification, In press.
View Article
Google Scholar

[40] View Article

[41] Google Scholar

[ref16] 16. Gray Roncal W, Koterba ZH, Mhembere D, Kleissas DM, Vogelstein JT, Burns, R, et al. (2013) MIGRAINE: MRI graph reliability analysis and inference for connectomics. Global Conference on Signal and Information Processing.

[ref17] 17. Herculano-Houzel S (2012) The remarkable, yet not extraordinary, human brain as a scaled-up primate brain and its associated cost. Proceedings of the National Academy of Sciences of the United States of America 109 Suppl: 10661–8. pmid:22723358
View Article
PubMed/NCBI
Google Scholar

[44] View Article

[45] PubMed/NCBI

[46] Google Scholar

[ref18] 18. Kolaczyk ED (2009) Statistical analysis of network data. Springer.

[ref19] 19. Fortin S (1996) The graph isomorphism problem. Technical Report, University of Alberta, Department of Computer Science.

[ref20] 20. Babai L (1981) Moderately exponential bound for graph isomorphism: 34–50.

[ref21] 21. Babali L, Erds P, Selkow SM (1980) Random graph isomorphism. SIAM Journal on Computing 9: 628–635.
View Article
Google Scholar

[51] View Article

[52] Google Scholar

[ref22] 22. Chen J (1994) A linear-time algorithm for isomorphism of graphs of bounded average genus. SIAM Journal on Discrete Mathematics 7: 614–631.
View Article
Google Scholar

[54] View Article

[55] Google Scholar

[ref23] 23. Zaslavskiy M, Bach FR, Vert JPP (2009) A path following algorithm for the graph matching problem. IEEE Transactions on Pattern Analysis and Machine Intelligence 31: 2227–2242. pmid:19834143
View Article
PubMed/NCBI
Google Scholar

[57] View Article

[58] PubMed/NCBI

[59] Google Scholar

[ref24] 24. Lyzinski V, Fishkind DE, Fiori M, Vogelstein JT, Priebe CE, Sapiro, G (2014) Graph matching: Relax at your own risk. arXiv preprint, arXiv: 14053133.

[ref25] 25. Sinkhorn R (1964) A relationship between arbitrary positive matrices and doubly stochastic matrices. The Annals of Mathematical Statistics 35: 876–879.
View Article
Google Scholar

[62] View Article

[63] Google Scholar

[ref26] 26. Bradley S, Hax A, Magnanti T (1977) Applied mathematical programming. Addison Wesley.

[ref27] 27. Anstreicher KM (2003) Recent advances in the solution of quadratic assignment. SIAM Journal on Optimization 97: 27–42.
View Article
Google Scholar

[66] View Article

[67] Google Scholar

[ref28] 28. Kuhn HW (1955) The Hungarian method for the assignment problem. Naval Research Logistics Quarterly 2: 83–97.
View Article
Google Scholar

[69] View Article

[70] Google Scholar

[ref29] 29. Jonker R, Volgenant A (1987) A shortest augmenting path algorithm for dense and sparse linear assignment problems. Computing 38: 325–340.
View Article
Google Scholar

[72] View Article

[73] Google Scholar

[ref30] 30. Burkard RE, Karisch SE, Rendl F (1997) QAPLIB: A quadratic assignment problem library. Journal of Global Optimization 10: 391–403.
View Article
Google Scholar

[75] View Article

[76] Google Scholar

[ref31] 31. Schellewald C, Roth S, Schnörr C (2001) Evaluation of convex optimization techniques for the weighted graph-matching problem in computer vision. In: Proceedings of the 23rd DAGMSymposium on Pattern Recognition. Springer-Verlag, pp. 361–368.

[ref32] 32. Singh R, Xu J, Berger B (2007) Pairwise global alignment of protein interaction networks by matching neighborhood topology: 16–31.

[ref33] 33. Liu ZY, Qiao H, Xu L (2012) An extended path following algorithm for graph matching problem. IEEE transactions on pattern analysis and machine intelligence 34: 1451–1456.
View Article
Google Scholar

[80] View Article

[81] Google Scholar

[ref34] 34. Gold S, Rangarajan A (1996) A graduated assignment algorithm for graph matching. IEEE Transactions on Pattern Analysis and Machine Intelligence 18: 377–388.
View Article
Google Scholar

[83] View Article

[84] Google Scholar

[ref35] 35. White JG, Southgate E, Thomson JN, Brenner S (1986) The structure of the nervous system of the nematode Caenorhabditis elegans. Philosophical Transactions of Royal Society London Series B, Biological Sciences 314: 1–340.
View Article
Google Scholar

[86] View Article

[87] Google Scholar

[ref36] 36. Varshney LR, Chen BL, Paniagua E, Hall DH, Chklovskii DB (2011) Structural properties of the Caenorhabditis elegans neuronal network. PLoS Computational Biology 7: 1–41.
View Article
Google Scholar

[89] View Article

[90] Google Scholar

[ref37] 37. Zaslavskiy M, Bach FR, Vert JPP (2010) Many-to-many graph matching: A continuous relaxation approach. Machine Learning and Knowledge Discovery in Databases 6323: 515–530.
View Article
Google Scholar

[92] View Article

[93] Google Scholar

[ref38] 38. Leordeanu M, Hebert M, Sukthankar R (2009) An integer projected fixed point method for graph matching and map inference. Advances in Neural Information Processing Systems: 1114–1122,
View Article
Google Scholar

[95] View Article

[96] Google Scholar

[ref39] 39. Escolano F, Hancock ER, Lozano M (2011) Graph matching through entropic manifold alignment. Computer Vision and Pattern Recognition.
View Article
Google Scholar

[98] View Article

[99] Google Scholar

[ref40] 40. Anstreicher KM, Brixius NW (2001) A new bound for the quadratic assignment problem based on convex quadratic programming. Mathematical Programming 89: 341–357.
View Article
Google Scholar

[101] View Article

[102] Google Scholar

[ref41] 41. Brixius NW (2000) Solving quadratic assignment problems using convex quadratic programming relaxations 1: Introduction. Management: 1–20.
View Article
Google Scholar

[104] View Article

[105] Google Scholar

[ref42] 42. Lyzinski V, Sussman DL, Fishkind DE, Pao H, Priebe CE (2013) Seeded graph matching for large stochastic block model graphs. arXiv preprint: 1310.1297.

[ref43] 43. Priebe CE, Vogelstein J, Bock D (2013) Optimizing the quantity/quality trade-off in connectome inference. Communications in Statistics-Theory and Methods 42: 3455–3462.
View Article
Google Scholar

[108] View Article

[109] Google Scholar

[ref44] 44. Lyzinski V, Fishkind D, Priebe C (2014) Seeded graph matching for correlated Erdos-Renyi graphs. Journal of Machine Learning Research To appear.
View Article
Google Scholar

[111] View Article

[112] Google Scholar

[ref45] 45. Boyd S, Vandenberghe L (2004) Convex Optimization. Oxford University Press.

Figures

Abstract

1 Introduction

2 Preliminaries

2.1 Quadratic Assignment Problems

2.2 Relaxed Quadratic Assignment Problem

2.3 Graph Matching

2.4 Graph Matching as a Quadratic Assignment Problem

3 Fast Approximate QAP Algorithm

4 Results

4.1 Algorithm Complexity and leading constants

4.2 QAP Benchmark Accuracy

4.3 QAP Benchmark Efficiency

4.4 QAP Benchmark Accuracy/Efficiency Trade-off

4.5 QAP Directed Benchmarks

4.6 Theoretical properties of FAQ

4.7 Multiple Restarts

4.8 Brain-Graph Matching

5 Discussion

5.1 Related Work

5.2 Future Work

5.3 Concluding Thoughts

Acknowledgments

Author Contributions

References