The Average Mixing Matrix Signature

Abstract

Laplacian-based descriptors, such as the Heat Kernel Signature and the Wave Kernel Signature, allow one to embed the vertices of a graph onto a vectorial space, and have been successfully used to find the optimal matching between a pair of input graphs. While the HKS uses a heat diffusion process to probe the local structure of a graph, the WKS attempts to do the same through wave propagation. In this paper, we propose an alternative structural descriptor that is based on continuous-time quantum walks. More specifically, we characterise the structure of a graph using its average mixing matrix. The average mixing matrix is a doubly-stochastic matrix that encodes the time-averaged behaviour of a continuous-time quantum walk on the graph. We propose to use the rows of the average mixing matrix for increasing stopping times to develop a novel signature, the Average Mixing Matrix Signature (AMMS). We perform an extensive range of experiments and we show that the proposed signature is robust under structural perturbations of the original graphs and it outperforms both the HKS and WKS when used as a node descriptor in a graph matching task.

1 Introduction

Graph-based representations have been used with considerable success in computer vision in the abstraction and recognition of object shape and scene structure [4, 8, 13]. A fundamental problem in graph-based pattern recognition is that of recovering the set of correspondences (matching) between the vertices of two graphs. In computer vision, graph matching has been applied to a wide range of problems, from object categorisation [1, 4] to action recognition [14, 15]. More formally, in the graph matching problem the goal is to find a mapping between the nodes of two graphs such that the edge structure is preserved.

While there exists a wide range of methods to solve this problem, many graph matching algorithms greatly benefit from the use of local structural descriptors to maximise their performance. To this end, a structural descriptor or signature is assigned to each node of the graphs, effectively embedding the graphs nodes onto a vectorial space. Given an underlying correspondence between the node sets of the graphs, the assumption underpinning these approaches is that corresponding nodes will be mapped to points that are close in the signature space.

Two structural signatures that have proven to be particularly successful are the Heat Kernel Signature (HKS) [12] and the Wave Kernel Signature (WKS) [2]. Both signatures belong to the family of Laplacian-based spectral descriptors and were introduced for the analysis of non-rigid three-dimensional shapes. More specifically, they are similarly based on the idea of using the spectrum of the graph Laplacian (and more in general the Laplace-Beltrami operator on the shape surface) to characterise the points on the shape surface. While the Heat Kernel Signature is based on a heat diffusion process, the Wave Kernel Signature is based on the propagation of a wavefunction.

In this paper, we introduce a novel structural signature based on a quantum mechanical process taking place over the graph. In particular, we propose to probe the graph structure using a continuous-time quantum walk [6]. We make use of the average mixing matrix [5] to define a quantum analogue of the heat kernel, and we propose to take the rows of the average mixing matrix at increasing times t as our vertex descriptor. The average mixing matrix was introduced by Godsil, and it describes the time-averaged probability of a quantum walks starting from node x to visit node y at time t. The motivation behind our work is based on the fact that quantum walks have proven to be very successful in characterising graph structures [3, 9–11]. Moreover, using the average mixing matrix allows us to overcome to lack of convergence of quantum walks. We show that the proposed signature is robust to structural noise and outperforms state-of-the-art signatures in a graph matching task.

The remainder of the paper is organised as follows: Sect. 2 reviews the necessary quantum mechanical background and Sect. 3 introduces the proposed signature. The experimental evaluation is presented in Sect. 4 and Sect. 5 concludes the paper.

2 Quantum Walks and Average Mixing Matrix

In this section, we introduce the necessary quantum mechanical background. We start by reviewing the define of continuous-time quantum walks, and then we show how to compute the average mixing matrix of a graph.

2.1 Continuous-Time Quantum Walks

The continuous-time quantum walk represents the quantum analogue of the continuous-time random walk [6]. Let $G = (V,E)$ denote an undirected graph over n nodes, where V is the vertex set and $E \subseteq V \times V$ is the edge set. In a continuous-time random walk, $\mathbf {p}(t) \in \mathbb {R}^n$ denotes the state of a walk at time t. The state vectors evolves according to the heat equation

$$\begin{aligned} \mathbf {p}(t) = e^{-Lt} \mathbf {p}(0), \end{aligned}$$

(1)

where the graph Laplacian L is the infinitesimal generator matrix of the underlying continuous-time Markov process.

Similarly to its classical counterpart, the state space of the continuous-time quantum walks is the vertex set of the graph. The classical state vector is replaced by a vector of complex amplitudes over V whose squared norm sums to unity, and as such the state of the system is not constrained to lie in a probability space, thus allowing interference to take place. The general state of the walk at time t is a complex linear combination of the basis states $\vert {u}\rangle$, i.e.,

$$\begin{aligned} \vert {\psi (t)}\rangle = \sum _{u\in V} \alpha _u(t) \vert {u}\rangle , \end{aligned}$$

(2)

where the amplitude $\alpha _u(t) \in \mathbb {C}$ and $\vert {\psi (t)}\rangle \in \mathbb {C}^{|V|}$ are both complex. Moreover, we have that $\alpha _u(t) \alpha _u^*(t)$ gives the probability that at time t the walker is at the vertex u, and thus $\sum _{u \in V} \alpha _u(t) \alpha ^{*}_u(t) = 1$ and $\alpha _u(t) \alpha ^{*}_u(t) \in [0,1]$, for all $u \in V$, $t \in \mathbb {R}^{+}$.

The evolution of the walk is governed by the Schrödinger equation

$$\begin{aligned} \frac{\partial }{\partial t}\vert {\psi (t)}\rangle = -iH\vert {\psi (t)}\rangle , \end{aligned}$$

(3)

where we denote the time-independent Hamiltonian as H. Generally speaking, a continuous-time quantum walk is induced whenever the structure of the graphs is reflected by the (0,1) pattern of the Hamiltonian. For example, we could take the adjacency matrix or the Laplacian. In the following we assume $H=A$.

Given an initial state $\vert {\psi (0)}\rangle$, solving the Schrödinger equation gives the expression of the state vector at time t,

$$\begin{aligned} \vert {\psi (t)}\rangle = U(t)\vert {\psi (0)}\rangle , \end{aligned}$$

(4)

where $U(t) = e^{-iAt}$ is the unitary operator governing the temporal evolution of the quantum walk. Equation 4 can be conveniently expressed in terms of the spectral decomposition of the adjacency matrix $A = \varPhi \varLambda \varPhi ^\top$, i.e., $\vert {\psi (t)}\rangle = \varPhi ^\top e^{-i \varLambda t} \varPhi \vert {\psi (0)}\rangle$, where $\varPhi$ denotes the $n \times n$ matrix $\varPhi =(\phi _1 | \phi _2 | ... | \phi _j | ... | \phi _n)$ with the ordered eigenvectors $\phi _j$s of H as columns and $\varLambda = \text{ diag }(\lambda _1, \lambda _2, ..., \lambda _j, ..., \lambda _n)$ is the $n \times n$ diagonal matrix with the ordered eigenvalues $\lambda _j$ of A as elements, and we have made use of the fact that $\text{ exp }[-iAt]=\varPhi ^\top \text{ exp }[-i \varLambda t] \varPhi$.

2.2 Average Mixing Matrix

Given a graph $G=(V,E)$ with adjacency matrix A and its associated unitary operator U(t), the behaviour of a continuous-time quantum walk at time t is captured by the mixing matrix [5]

$$\begin{aligned} M_G(t) = e^{iAt} \circ e^{-iAt} = U(-t) \circ U(t), \end{aligned}$$

(5)

where $A \circ B$ denotes the Schur-Hadamard product of two matrices A and B.

Note that, while in the classical case the probability distribution induced by a random walk converges to a steady state, this does not happen in the quantum case. However, we will show that we can enforce convergence by taking a time-average even if U(t) is norm-preserving. Let us define the average mixing matrix [5] at time T as

$$\begin{aligned} \widehat{M}_{G;T} = \frac{1}{T}\int _0^T M_G(t) \, \mathrm {d} t. \end{aligned}$$

(6)

The entry uv can be interpreted as the average probability that a walker is found at vertex v when starting the walk from vertex u. Let us define $P_{\lambda } = \sum _{k=1}^{\mu ({\lambda })} \phi _{\lambda ,k} \phi _{\lambda ,k}^\top$ to be the projection operator on the subspace spanned by the $\mu (\lambda )$ eigenvectors $\phi _{\lambda ,k}$ associated with the eigenvalue $\lambda$ of A. Given this set of projectors, the unitary operator inducing the quantum walk can be rewritten as

$$\begin{aligned} U(t) = \sum _{\lambda } e^{-i\lambda t} P_{\lambda } \end{aligned}$$

(7)

Given Eq. 7, we can rewrite the equation for the mixing matrix as

$$\begin{aligned} M_G(t) = \sum _{\lambda _1 \in \varLambda } \sum _{\lambda _2 \in \varLambda } e^{-i(\lambda _1 - \lambda _2) t} P_{\lambda _1} \circ P_{\lambda _2}, \end{aligned}$$

(8)

and thus we can reformulate Eq. 6 as

$$\begin{aligned} \widehat{M}_{G;T} = \sum _{\lambda _1 \in \varLambda } \sum _{\lambda _2 \in \varLambda } P_{\lambda _1} \circ P_{\lambda _2} \frac{1}{T} \int _{0}^T \! e^{-i(\lambda _1 - \lambda _2) t} \, \mathrm {d} t, \end{aligned}$$

(9)

which has solution

$$\begin{aligned} \widehat{M}_{G;T} = \sum _{\lambda _1 \in \varLambda } \sum _{\lambda _2 \in \varLambda } P_{\lambda _1} \circ P_{\lambda _2} \frac{i(1-e^{iT(\lambda _2 - \lambda _1)})}{T(\lambda _2 - \lambda _1)}. \end{aligned}$$

(10)

In the limit $T \rightarrow \infty$, Eq. 10 becomes

$$\begin{aligned} \widehat{M}_{G;\infty } = \sum _{\lambda \in \tilde{\varLambda }} P_{\lambda } \circ P_{\lambda } \end{aligned}$$

(11)

where $\tilde{\varLambda }$ is the set of distinct eigenvalues of the adjacency matrix.

2.3 Properties of the Average Mixing Matrix

Given a graph $G = (V,E)$, let $\widehat{M}_{G;\infty }(i,j)$ denote the element of $\widehat{M}_{G;\infty }$ at the ith row and jth column. Recall that two nodes $u,v \in V$ are strongly cospectral if and only if $P_{\lambda }e_u = \pm P_{\lambda }e_v$ for all eigenvalues $\lambda$, where $e_u$ is an element of the standard basis of $\mathbb {R}^{|V|}$. It can be shown that the following properties hold for the matrix $\widehat{M}_{G;\infty }$:

• 1. $\widehat{M}_{G;\infty }$ is doubly stochastic and rational, and that $\widehat{M}_{G;\infty }$ is positive semidefinite whenever G is connected [5];

• 2. Let G and H be cospectral, regular, and switching equivalent graphs. Then $\widehat{M}_{G;\infty }(j,j)=\widehat{M}_{H;\infty }(j,j)$, for every $j=1,2,...,n$;

• 3. Let G and H be cospectral, walk-regular graphs. Then $\widehat{M}_{G;\infty }(j,j)=\widehat{M}_{H;\infty }(j,j)$, for every $j=1,2,...,n$;

• 4. two vertices i and j are strongly cospectral if and only if the two corresponding rows of the average mixing matrix are the same, i.e., $\widehat{M}_{G;\infty }(i,k)=\widehat{M}_{G;\infty }(j,k)$, for each k [5];

• 5. suppose G is a regular and connected graph and $\overline{G}$ is the complementary graph. Then $\widehat{M}_{G;\infty }=\widehat{M}_{\overline{G};\infty }$ if and only if $\overline{G}$ is connected.

We omit the proof of 2, 3 and 5 for lack of space, while the proofs of 1 and 4 can be found in  [5]. Note that 2, 3, 4 and 5 outline the conditions under which some diagonal element of the mixing matrix of can be duplicated, the rows of two mixing matrices are identical, and two mixing matrices have exactly the same elements, respectively. In other words, the above properties give an important indication of the limitations of the matrix $\widehat{M}_{G;\infty }$ to discriminate between different graphs. Although the same properties have not been proved for an arbitrary stopping time T, our experimental analysis in Sect. 4 suggests that some of these properties may still hold for $T < \infty$.

3 The Average Mixing Matrix Signature

Given a graph $G = (V,E)$ and a vertex $x \in V$, we define the Average Mixing Matrix Signature (AMMS) at x at time T as

$$\begin{aligned} AMMS(x,T) = \mathtt {sort}(\widehat{M}_{G;T}(x,-)), \end{aligned}$$

(12)

where $\widehat{M}_{G;T}(x,-)$ denotes the row of the average mixing matrix corresponding to x, and $\mathtt {sort}$ is a sorting function. Given a time interval $[T_{min},T_{max}]$, the AMMS is then created by concatenating the sorted rows for every $T \in [T_{min},T_{max}]$.

As Eq. 12 shows, for each stopping time T, we decide to take the whole row of the average mixing matrix rather than just the diagonal element $\widehat{M}_{G;T}(x,x)$. Recall from Sect. 2 that under some particular conditions two graphs can have the same diagonal entries of $\widehat{M}_{G;\infty }$. Although the same has not been proved for an arbitrary T, we take the entire row in an attempt to maximise the discriminatory power of the signature. Finally, to ensure the permutational invariance of the signature, for each T we sort the elements of the average mixing matrix rows.

Note that, while there is no specific reason to take all the elements of the sorted row rather than simply the first k, in the experimental evaluation we show that the best performance is usually achieved for large values of k, so we propose to take the entire row. Finally, from Eq. 10 we observe that the computational complexity of the proposed signature is $O(|V|^4)$.

4 Experimental Evaluation

To evaluate the proposed signature, we compare it on a graph matching problem against the HKS as well as the WKS. In the following, we refer to the percentage of correctly matched nodes as accuracy.

4.1 HKS and WKS

Heat Kernel Signature. The HKS is based on the analysis of a heat diffusion process on the shape governed by the heat equation

$$\begin{aligned} -\frac{\partial u(x,t)}{\partial t} = \varDelta u(x,t), \end{aligned}$$

(13)

where $\varDelta$ is the Laplace-Beltrami operator. Assuming at time $t=0$ all the heat is concentrated at a point x, it can be shown that the amount of heat that remains at x at time t is

$$\begin{aligned} HKS(x,t) = \sum _{i=0}^\infty e^{-\lambda _i t} \phi _i(x)^2 \end{aligned}$$

(14)

where $\lambda _i$ and $\phi _i$ denote the ith eigenvalue and the ith eigenfunction of the Laplace-Beltrami operator, respectively. The HKS at a point x is then constructed by evaluating Eq. 14 for different $t \in [t_{min},t_{max}]$. As suggested in [12], here we sample 100 points logarithmically over the time interval defined by $t_{min} = 4\ln (10)/\lambda _{max}$ and $t_{max} = 4\ln (10)/\lambda _{min}]$, where $\lambda _{min}$ and $\lambda _{max}$ denote the minimum and maximum eigenvalue of the graph Laplacian, respectively.

Wave Kernel Signature. The WKS is based on the analysis of wavefunction evolution on the shape surface governed by the Schrödinger equation

$$\begin{aligned} -\frac{\partial \psi (x,t)}{\partial t} = i\varDelta \psi (x,t). \end{aligned}$$

(15)

Although similar in the definition, note that the dynamics behind the WKS (oscillation) and the HKS (dissipation) are fundamentally different. The WKS evaluates the expected probability of measuring a quantum particle in the point x of the shape surface at any time, i.e.,

$$\begin{aligned} WKS(x,e) = \sum _{i=0}^\infty \phi _i(x)^2 f_e(\lambda _i)^2, \end{aligned}$$

(16)

where $\lambda _i$ and $\phi _i$ denote the ith eigenvalue and the ith eigenfunction of the Laplace-Beltrami operator, respectively, and $f_e$ is a log-normal energy distribution that depends on a parameters e. The WKS for a point x is then constructed by evaluating Eq. 16 for different values of $e \in [e_{min},e_{max}]$. Here we set all the parameters of the WKS as illustrated in [2].

4.2 Datasets

We use a dataset made of the images from the CMU house sequence, where each image is abstracted as a Delaunay graph over a set of corner feature points. All the resulting graphs are composed of 30 to 32 nodes. Figure 1 shows the ten images with the feature points and the graphs superimposed.

Fig. 1.

CMU house sequence with the feature points and the Delaunay graphs superimposed.

Fig. 2.

Evaluation of the sensitivity to the parameters values. Here k denotes the number of row elements used to construct the signature and $T_{max}$ is the maximum stopping time considered.

In addition, we consider a synthetic datasets consisting of 10 Delaunay graphs over 20 nodes. This dataset is generated as follows. We create a first graph by generating 20 computing the Delaunay triangulation of 20 uniformly randomly scattered two-dimensional points. Then, we generate 9 additional graphs by slightly perturbing the original points and recomputing the Delaunay triangulation.

4.3 Graph Matching Experimental Setup

For a pair of input graphs, we compute the AMMS for each graph node. Here we decide to sample 10 points logarithmically over the time interval [0.1, 1]. The choice of $T_{max} = 1$ is experimentally motivated in the next subsection. In general, we observe that on all the datasets considered in this paper the best performance is achieved for low values of $T_{max}$. Given the nodes signatures, we then compute the matrix of pairwise Euclidean distances between the signatures. Finally, we cast the graph matching problem as an assignment problem which can be solved in polynomial time using the well-known Hungarian method [7].

4.4 Sensitivity to Parameters

We first investigate the sensitivity of the proposed signatures to the number k of row elements that we consider for each stopping time $T \in [T_{min},T_{max}]$, as well the value of $T_{max}$. Staring from a seed graph G from the synthetic dataset, we generate 100 noisy copies $G'$ by removing or adding an edge with probability $p = 0.03$. Then, we compute the matching for each pair $(G,G')$ and each choice of the parameters k and $T_{max}$, and we plot the average accuracy in Fig. 2. Here we let $k = 1,2,\cdots ,n$, where n is the number of nodes of G, and we let $T_{max} = 1,2,\cdots ,100$. Figure 2 shows that the best accuracy is achieved for k larger than 10 and $T_{max} = 1$.

Fig. 3.

Average accuracy for increasing levels of Erdős-Rényi noise. Here p denotes the probability of adding or removing an edge.

4.5 Robustness to Erdős-Rényi Structural Noise

In order to evaluate the robustness of the AMMS to structural noise, we perform the following experiment. Given a graph G from the synthetic dataset, we create a series of noisy versions of G, for increasing amounts of Erdős-Rényi noise. To this end, we flip an edge of the graph with probability p and obtain the noisy graph $G'$. For each value of p, we repeat the noise addition process 100 times to obtain 100 noisy copies of G. Then, we compute the optimal matching using the proposed signature and we compare the results against those obtained using the HKS and the WKS. Figure 3(a) shows the average accuracy for increasing levels of structural noise. Note that the AMMS and the WKS are significantly less susceptible to noise than the HKS. It is also interesting to observe that as p grows over 0.8, the number of correctly matched nodes for the AMMS start to increase (see Fig. 3(b)). We posit that this is evidence that the average mixing matrix of a graph G is equivalent or very similar to that of its complement $\overline{G}$ also for arbitrary $T < \infty$. Note that for $T \rightarrow \infty$ this is true if G is regular and connected, and $\overline{G}$ is connected as well. Although this is not the case for the graphs used in these experiments, the geometric nature of the Delaunay triangulation yields indeed fairly regular graphs.

Table 1. Average accuracy (± standard error).

4.6 Graph Matching Accuracy

Table 1 shows the average accuracy on the CMU house and the synthetic datasets. We consider two alternative versions of the proposed signatures: (1) AMMS(diag) is the signature where only the diagonal elements $\widehat{M}_{G;T}(x,x)$ are used; (2) AMMS(row) is the signature where all the rows elements are used. In both the datasets considered, our signature performs significantly better than both the HKS and the WKS. Moreover, the AMMS(row) always outperforms the AMMS(diag), highlighting the importance of using all the elements of the row and suggesting once again that the properties listed in Sect. 2 hold for an arbitrary T. Finally, we note that the AMMS(diag) always outperforms the HKS but not the WKS.

5 Conclusion

In this paper, we have introduced a novel structural signature for graphs. We probed the structure of the graph using continuous-time quantum walks and we proposed to consider the rows of the average mixing matrix for increasing stopping times T as the node descriptor. The experimental results show that this signature can outperform state-of-the-art signatures like the Heat Kernel Signature and the Wave Kernel Signature in a graph matching task. Future work should focus on the optimisation of the signature parameters for general graphs, as well as reducing the computational complexity associated with the Average Mixing Matrix Signature. Recall that we observed that the best performance is usually achieved for low values of $T_{max}$. Therefore, one idea may be to reduce the computational complexity of the proposed signature by resorting to a linear approximation of its values.