Graph Edit Distance Learning via Modeling Optimum Matchings with Constraints

1 Introduction

Graphs are ubiquitous and have been used in a wide range of applications, such as bioinformatics [Peng et al., 2015], software engineering [Xu et al., 2017], computer vision [Yan et al., 2020] and so on. A fundamental problem in graph data analysis is to compute the similarity between graphs. For example, in bioinformatics, molecules with similar graphs often have similar functions. Computing the similarity of molecular graphs is widely used in virtual screening of molecules. In software engineering, the control-flow of a code fragment can be modeled as a graph. Detecting the similarity of controlflow graphs is useful for code plagiarism detection and code vulnerability search.

Many graph similarity measures are proposed to measure the similarity of graphs. Among others, graph edit distance (GED), i.e., the minimum edit cost of transforming one graph to another, is widely used and domain-agnostic [Bai et al., 2020; Li et al., 2019; Chang et al., 2020; Zhao et al., 2018]. However, computing GED is NP-hard [Zeng et al., 2009].

A recent work [Blumenthal and Gamper, 2020] shows that the exact GED between graphs of more than 16 nodes cannot be computed reliably by the state-of-the-art algorithms in a reasonable time. Therefore, many inexact GED computation methods have been proposed [Blumenthal et al., 2020]. However, it is challenging to handcraft optimal algorithms, even by experts, for different definitions of edit cost.

Meanwhile, with the success of deep learning, learningbased graph similarity/matching computation methods have been receiving increasing research attention [Ma et al., 2019; Yan et al., 2020]. They learn graph similarity/matching models from data by minimizing the loss between the predicted graph similarity/matching and the ground-truth. When used for inference, the learned models can effectively predict the similarity/matching of unseen graphs.

It is well known that GED has two inherent characteristics: multiple optimum node matchings, i.e., two graphs have many node matchings that are optimum to produce their GED, and one-to-one node matching constraints, i.e., one node of a graph can match only one node of another graph in an optimum node matching. However, these two characteristics have not been well considered in the loss functions of existing graph similarity/matching learning works. For instance, GraphSim [Bai et al., 2020], the state-of-the-art graph similarity learning method, employs the mean-square error (MSE) of the predicted graph similarity as its loss function. This generic loss function allows GraphSim to support different types of graph similarity measures. Nevertheless, GraphSim can only make a suboptimal prediction of GED, as the inherent characteristics of GED are not encoded. On the other hand, GLMNet [Jiang et al., 2019] has encoded the one-to-one node matching constraints in the loss function. Let a binary matrix X denote the node matching between two input graphs and ˆX denote the prediction. A one-to-one node matching requires that each row and column of X should have just a 1. GLMNet proposes the constraint regularizer Lcr that forces the product of each pair of elements in each row or column of ˆX to be 0. Lcr itself cannot constrain the one-to-one node matching, as it says nothing to the 1's in ˆX. But, the sum of Lcr and the cross entropy between ˆX and X can constrain the one-to-one node matching. However, GED may have multiple optimum node matchings. If we just select one optimum node matching as the ground-truth, all other optimum node matchings are overkilled and the model trained is suboptimal. Moreover, it is non-trivial to modify the loss function of GLMNet to work with multiple optimum node matchings. The main reason is that Lcr needs the cross entropy to constrain the one-to-one node matching, but the cross entropy does not work when there are multiple optimum node matchings.

In this paper, we propose a novel GED-specific loss function to address the deficiency of existing works. First, we propose a new partial node matching-based regularizer Lm to encode multiple optimum node matchings. Lm is inspired by the Branch-and-Bound (BnB) search for exact GED computation and the Monte Carlo tree search (MCTS). The main ideas are as follows. i) We encode all the optimum node matchings covered by a sub-tree of the BnB tree in the ground-truth. The optimum node matchings are leaves of the BnB tree and each internal node is a partial node matching. We use the partial node matching of an internal node as the ground-truth and just penalize the mispredicted 1's in the partial node matching. All the optimum node matchings covered by the sub-tree are not penalized by the regularizer. ii) We select the optimal internal node to balance the mispenalty and the restraint of the regularizer. The reason is that with larger coverage of the internal node, fewer optimum node matchings not covered are mispenalized, which benefits the learning of the model. However, the regularizer is less restrictive, which is adverse to the learning of the model. An internal node is given a score. Inspired by the MCTS, the score is computed by random rollouts. If the rollouts from an internal node reach more optimum node matchings, the internal node has a higher score. If two internal nodes have the same score, we select the one having larger coverage. Using only one optimum node matching as the ground-truth is a special case of our regularizer.

Second, we propose a plane intersection-based regularizer Lc to encode the one-to-one constraints for the node matchings encoded by Lm. Specifically, we leverage the constrained binary quadratic programming (CBQP) formulation of GED. The one-to-one node maching constraints in the CBQP can be modeled as the quadratic infeasibility penalty term Lqip in the loss function [Kochenberger et al., 2014]. However, Lqip becomes not effective to constrain the one-toone node matching after the binary variables in the CBQP are relaxed to be continuous in [0, 1]. We observe that the intersection of the solution plane of minimizing Lqip and that of minimizing the constraint regularizer Lcr is effective to constrain the one-to-one node matching. Therefore, we combine Lqip and Lcr as our plane intersection-based regularizer.

To learn the cross-graph representation of two input graphs G1 and G2, we build an association graph GX from the CBQP of G1 and G2. The nodes in GX model the variables of the CBQP and the edges model the relationship between the variables. GX can also be regarded as the "product" of G1 and G2. The graph neural network on GX is used to learn the embeddings of the nodes in GX. The node embeddings are fed into fully-connected layers to predict the GED. Our experiments show that our method is 4.2x-103.8x more accurate than the existing graph similarity learning methods and heuristic GED algorithms. In summary, this paper's contributions are as follows.

• We propose a novel partial node matching-based regularizer to encode multiple optimum node matchings. An optimal partial node matching selection method is proposed.
• We propose a plane intersection-based regularizer to impose the one-to-one constraints for the node matchings encoded by the partial node matching-based regularizer.
• We conduct extensive experiments to verify that our proposed techniques significantly outperform the state-ofthe-art methods on several benchmark datasets.

The rest of this paper is organized as follows. Section 2 summarizes the related work. Section 3 presents the preliminaries. The GED-specific loss function and model architecture are proposed in Section 4. Section 5 presents the experimental results and Section 6 concludes this paper.

2 RelatedWork

This section discusses the recent works that are closely related to this paper.

2.1 Graph Similarity Learning

Most recent graph similarity learning works use the graph neural networks (GNNs) to learn the embeddings of the graphs. The embeddings are fed into multilayer perceptrons (MLPs) or convolutional neural networks (CNNs) to predict the graph similarity. For example, the work [Wang et al., 2019b] uses a Siamese GNN model to predict the graph similarity by regression. The work [Xu et al., 2017] uses a Siamese GNN model to predict if two graphs are similar by binary classification. Different from fusing the embeddings of the two graphs at the very end of learning in above methods, the state-of-the-art works fuse the embeddings of the two graphs early via cross-graph learning. The work [Li et al., 2019] proposes a cross-graph attention network GMN. The loss function of GMN is the classification error of predicting if two graphs are similar. The work GraphSim [Bai et al., 2020] combines GNN and CNN. The loss function of GraphSim is the MSE of the predicted graph similarity. However, as the two inherent characteristics of GED have not been considered, GMN and GraphSim produce a suboptimal prediction. Interested readers are referred to [Ma et al., 2019] for a recent survey.

2.2 Deep GraphMatching

Graph matching is widely studied in computer vision to match images. For example, the work [Wang et al., 2019a] conducts cross-graph learning on the association graph built from the graphs extracted from two images. The embeddings of the nodes of the association graph are fed into the Sinkhorn layer to predict the node matching. The cross entropy of the predicted node matching and the ground-truth is used as the loss. In [Yu et al., 2020], the Hungarian attention is proposed to pay more attention to the mismatched nodes to address the overfitting issue caused by the cross entropy loss. The work [Jiang et al., 2019] introduces a constraint regularizer to encode the one-to-one node matching constraints.

The graph matching methods cannot be directly adopted for GED learning as they do not consider multiple optimum node matchings. In addition, graph matching has an inherent characteristic of neighborhood consensus [Fey et al., 2020], whereas GED does not.

There are also some works studying learning-based methods to solve the quadratic programming and other problems that can be modeled by the quadratic programming [Vesselinova et al., 2020]. However, the inherent characteristics of GED are also not considered.

3 Preliminaries

This paper studies undirected graphs, or simply called graphs. A graph is denoted as G = (V,E), where V and E are the node set and the edge set of G, respectively. The label of a node u and an edge (u, u′ ) is denoted by ℓ(u) and ℓ(u, u′ ), respectively. N(u) denotes the neighbors of u ∈ V . The size of G, denoted by |G|, is the number of nodes in G.

A graph G1 can be transformed to another graph G2 by a sequence of six types of edit operations: node/edge insertion, node/edge deletion, and node/edge relabeling. The sequence of edit operations is called as the edit path. The cost of the edit path is the total cost of the edit operations in it. Graph edit distance (GED) of G1 and G2, denoted by d(G1,G2), is the minimum cost of an edit path that transforms G1 to G2. In this paper, we assume that the cost of each edit operation is 1, yet our technique can be easily extended to other costs. In an optimum edit path, each node/edge of G1 is edited at most for one time, as otherwise the edit path's cost is not minimum.

We assume G1 and G2 have the same number of nodes, as otherwise we can add dummy nodes with a special label to the smaller graph as in [Chang et al., 2020]. An example of adding dummy nodes is shown in Figure 1(a)-(b).

An optimum edit path is equivalent to an optimum node matching between G1 and G2 [Blumenthal et al., 2020]. The intuition is as follows. If a node u or an edge (ui, uj) of G1 is matched to a node or an edge in G2 that has a different label, u or (ui, uj) is relabeled. If an edge (ui, uj) of G1 is matched to a non-existing edge of G2, (ui, uj) is deleted. If a non-existing edge of G1 is matched to an edge (vk, vl) of G2, (vk, vl) is inserted to G1.

Let a binary matrix X denote the node matching between G1 and G2. The optimum X can be computed by the constrained binary quadratic programming (CBQP) as follows.

where ci,k is the edit cost of matching ui in G1 and vk in G2. ci,k = 1 if ℓ(ui) ̸= ℓ(vk) and 0 otherwise. ci,j,k,l is the edit cost of matching the edge (ui, uj) in G1 and the edge (vk, vl) in G2. ci,j,k,l = 1 if ℓ(ui, uj) ̸= ℓ(vk, vl) and 0 otherwise.

Equations (2-4) constrain that the optimum node matching must be a one-to-one node matching. Otherwise, the node

(a) two graphs with GED = 2

Figure 2: A Brand-and-Bound (BnB) tree for computing the GED of the two graphs in Figure 1(b). Subtrees that do not cover any optimum node matching are pruned during the BnB search. The circled leaves are optimum node matchings.

matching cannot transform G1 to G2 or the edit cost is not the minimum. Two graphs can have multiple optimum node matchings, as shown in Figure 1(c), as there are many ways to transform a graph to another with the same minimum edit cost. A predicted node matching is shown in Figure 1(d).

The Brand-and-Bound (BnB) tree for computing the exact GED between G1 and G2 in Figure 1(b) is shown in Figure 2. Each internal node is a partial node matching. π is the order to processing the nodes in G1, which is decided by heuristics in the latest methods [Chang et al., 2020].

4 Graph Edit Distance Learning

Following the CBQP of GED, we formulate a GED learning problem. The key idea is to design a GED-specific loss function. Then, we use a graph neural network on the association graph built from the variables of the CBQP to minimize the GED-specific loss.

4.1 ProblemTransformation

We transform the CBQP to a graph learning problem, which can model the two inherent characteristics of GED. Specifically, the CBQP is equivalent to the unconstrained binary quadratic programming (UBQP) by moving the node matching constraints (Equations (2) and (3)) as the quadratic infeasibility penalty to the objective, as follows.

where Lqip is the quadratic infeasibility penalty and α is the weight.

We further relax the binary variables in the UBQP to [0, 1] to obtain an unconstrained quadratic programming (UQP) such that the variables are differentiable. Then, a GED learning problem can be obtained by minimizing the loss function corresponding to the objective of the UQP as follows.

where Ld is the mean-square error (MSE) of dist and the ground-truth GED.

4.2 GED-specific Loss Function

Regularizer For Multiple Optimum Node Matchings

The main ideas of the regularizer are as follows. i) We use a partial node matching as the ground-truth, s.t. all the optimum node matchings covered by the partial node matching are not penalized. In this way, multiple optimum node matchings are encoded. ii) We select the optimal partial node matching to balance the mispenality and the restraint of the regularizer.

Specifically, let a binary matrix P denote the partial node matching of an internal node of the BnB tree. P has at most a 1 in each row and column. Our regularizer just penalizes the mispredicted 1's in P as follows.

where ˆX is the predicted node matching. Since the 0's in P are not penalized, all the optimum node matchings covered by the internal node are not penalized. Note that Lm is not the cross entropy and it is different from the Hungarian loss in [Yu et al., 2020], as only the mispredicted 1's are penalized.

We propose a random rollout-based method to select the optimal P . Let X denote an optimum node matching between G1 and G2. X is a leaf of the BnB tree. Let path denote the path from the root of the BnB tree to the leaf X. At the ith node in path, denoted by path[i], we conduct k random rollouts. In each rollout, we find a random path from path[i] to a leaf of the BnB tree. If the leaf reached is an optimum node matching, the rollout is called positive. The score of path[i] is the percentage of positive rollouts. If path[i] and path[j] with i < j has the same score, we select path[i] as the optimal internal node, because path[i] is closer to the root of the BnB tree and the subtree rooted at path[i] has a larger coverage.

Lm cannot constrain the node matchings to satisfy the oneto-one constraints. We next propose another regularizer to encode the one-to-one node matching constraints.

Regularizer To Constrain One-to-one Node Matching

When the variables in the quadratic programming are binary, Lqip can constrain the one-to-one node matching. However, after the variables are relaxed to be continuous in [0, 1], Lqip is not effective. We use an example in the 3D space to show the shortcoming.

Example 1. Suppose we have three binary variables x, y, z ∈ {0, 1} and we expect there is only one 1 in x, y, z. The solution points are sol = {(0, 0, 1), (0, 1, 0), (1, 0, 0)} in 3D space. If we use the constraint x+y+z = 1, the valid points are exactly sol. However, if we relax x, y, z to [0, 1] and minimize (1-(x+y+z))2 , all points on the plane x+y+z = 1

Figure 3: An example showing the effectiveness of Lc and Lm

are valid. sol are the intersection points of the plane with the X, Y and Z axes. It is not easy to find sol by only minimizing (1 - (x + y + z))2 with x, y, z ∈ [0, 1].

The work [Jiang et al., 2019] proposes the constraint regularizer Lcr to encode the one-to-one node matching constraints as follows.

where wi,j,k,l = 1 if i = j, k ̸= l or i ̸= j, k = l and 0 otherwise. Intuitively, Lcr is to force the product of each pair of elements in each row/column in X to be 0. However, Lcr is also not effective to constrain the one-to-one node matching after the continuous relaxation, as shown by the following example.

Example 2. Continuing Example 1, when x, y, z ∈ {0, 1}, sol can be obtained by minimizing xy +xz +yz. However, if we relax x, y, z to [0, 1], minimizing xy + xz + yz produces the X, Y and Z axes. It is also not easy to find sol. However, sol are exactly the intersection points of the solution of minimizing xy+xz+yz and that of minimizing (1-(x+y+z))2 .

From Examples 1 and 2, we can observe that after relaxing the variables to [0, 1], either Lqip or Lcr does not constrain the one-to-one node matching, but the intersection of their solution planes is a desired choice. Therefore, we propose a plane intersection-based regularizer for encoding the one-toone node matching constraints as follows.

Finally, our overall loss function is as follows.

where α and β are two hyperparameters.

Discussion

Lc and Lm are crucial for GED learning and should be used together to train an effective model. We illustrate the merits of our method by examples as shown in Figure 3.

• If we only use Lc and do not use Lm, we will obtain a suboptimal model. For example, given two graphs with GED = 1 as shown in Figure 3(a), the model converges at the node matching as shown in Figure 3(b) and the predicted GED is ∼2.3.

• If we use both Lc and Lm, but use the full node matching as the ground-truth in Lm, we will again obtain a suboptimal model. For example, if we use the samples in Figure 3(e) as the training data, the model converges at the node matching as shown in Figure 3(c) and the predicted GED is ∼1.5.
• We propose to use both Lc and Lm, and use the partial node matching as the ground-truth in Lm. We will obtain an effective model. For example, if we use samples (G1,G2,X5), (G1,G2,X5) and (G1,G2,X6) in the training data (X5 and X6 are partial node matchings shown in Figure 3(f)), the predicted GED is 0.99.

4.3 Graph Neural Network-basedModel

We follow the method of [Jiang et al., 2019] to build an association graph of G1 and G2 and then use graph convolution on the association graph to learn the representation crossing G1 and G2. The embeddings of the nodes are fed into multilayer perceptrons (MLPs) for GED prediction.

The association graph GX is built from the UBQP of G1 and G2. Specifically, for each variable in the UBQP, GX has a node. If ci,k = 1, the node Xi,k is of type 1 and 0 otherwise. Let t(·) denote the type of a node. If ci,k,j,l = 1, GX has an edge (Xi,k,Xj,l). Assume |VG1 | ≤ |VG2 | and let n = |VG2 |. GX has n2 nodes and |EG1 |n2 + |EG2 |n2 edges. If G1 and G2 are sparse, i.e., |EG1 | = O(|VG1 |) and |EG2 | = O(|VG2 |), GX is also sparse and GX has O(n3 ) edges.

The graph convolution on GX differentiates the types of the nodes as follows.

where u is a node in GX, ⃗hu is the embedding vector of u, N(u) is the neighbors of u in GX, L is the layer of graph convolution, W0 and W1 are trainable parameters for the nodes of type 0 and type 1, respectively, and σ is the activation function.

The time complexity of the graph convolution on GX is O(n3 ). Although it is larger than O(n2 ) of GMN and GraphSim, our model has a higher capacity and can produce more accurate GED prediction.

5 Experimental Evaluation

We evaluate our method against a number of state-of-the-art methods for GED computation, with the major goals of addressing the following questions:

• Q1: How effective is our method compared to the state-ofthe-art methods of GED computation, including both heuristic algorithms and neural network based models?
• Q2: How do the regularizer Lm for encoding multiple optimum node matchings and the regularizer Lc for node mapping constraints help with the GED learning?
• Q3: How our method can generalize to larger graphs that have different distributions of GED and are not seen during training?

Datasets. We use four real graph datasets AIDS, IMDB, LINUX and PTC that are from different domains in our experiments. The datasets are the same as those used in [Bai et al., 2020]. We use G≤30 and G>30 to denote the graphs of no more than and more than 30 nodes, respectively. We use G≤30 to answer Q1 and Q2, and use G>30 to answer Q3. Table 1 shows the statistics of the datasets. The GEDs and optimum node matchings are computed by the latest exact GED method [Chang et al., 2020].

Baseline methods. We follow GraphSim to compare with three well-known heuristic GED computation algorithms Beam [Neuhaus et al., 2006], Hung [Riesen and Bunke, 2009], and VJ [Fankhauser et al., 2011]. We also compare with two state-of-the-art GED learning methods GMN [Li et al., 2019] and GraphSim [Bai et al., 2020].

Parameter settings. We use two graph convolution layers and ReLU as the activation function. We use the one-hot encoding of node degree as the initial node embedding. The embedding dimensions are 32. The embedding of a node Xu,v in the association graph GX is fed into a fully-connected layer followed by a sigmoid function to predict Xu,v. The mean of the embeddings of all nodes in GX is fed into a fullyconnected layer to predict GED. The average mean-square error (MSE) of the predicted GED and the ground-truth on the test data is used as the performance metric. α and β are set by grid search. We set the batch size to 1 and use the Adam optimizer. The initial learning rate is 0.005 and reduced by 0.96 for each 5 epochs. We set the number of epochs to 600, and select the best model based on the lowest MSE of GED on validation data. The experiments are conducted using PyTorch on a server with Intel CPU Xeon Gold 6230R, 768G RAM, and a GPU card NVIDIA Tesla K80. The source code is available online. 1

5.1 Comparison with Baseline Methods

To answer Q1 with this experiment, we sample 6k, 2k, 2k pairs of graphs in G≤30 as the training data, the validation data and the test data, respectively. The MSE of GED prediction on the test data is shown in Table 2. From Table 2, we can observe that our method significantly outperforms the state-of-the-art GED learning methods on all datasets. Specifically, the prediction error of our method is respectively 10.6x, 30.2x, 5.4x, and 103.8x smaller than GraphSim and even smaller than GMN on AIDS, LINUX, IMDB and PTC. Our method is 21.2x, 4.2x, and 4.9x better than the approximate algorithms on AIDS, LINUX and PTC, respectively, and is on par with them on IMDB. The IMDB dataset is harder as the graphs are much denser. The average density of IMDB is 0.83, whereas that of AIDS is only 0.09. We also compare with GLMNet. For example, the MSE of GLMNet on AIDS is 265.6. The average GED of the test data is 23.75, 22.92, 52.77, and 35.2 on AIDS, LINUX, IMDB, and PTC, respectively. Our MSEs 1.32 and 1.16 are very small comparing to the GEDs on AIDS and LINUX.

1 https://github.com/csypeng/graph edit distance learning

Besides the MSE, we examine the Spearman's rank correlation coefficient ρ, Kendall's rank correlation coefficient τ, and precision at 10 (p@10). For example, on AIDS, ρ, τ , and p@10 of our method is 0.89, 0.85, and 0.73, whereas those of GraphSim are 0.86, 0.81, and 0.61, respectively.

We also compare the efficiency and the result is shown in Table 3. From Table 3, we can observe that the approximate GED algorithms are faster than the learning-based methods. Our method is slower than the latest GED learning methods. The reason is that the time complexity of the graph convolution of our method is O(n3 ), where n is the number of nodes in the larger one of the two input graphs. The time complexities of GraphSim and GMN are O(n2 ). Our method is still practical as n is small in many real applications.

5.2 Ablation Study

To investigate Q2, in this experiment, we examine the effectiveness of our proposed regularizer Lm and Lc. The result is shown in Table 4. From Table 4, we can observe that Lm is important for GED learning. The best results on all datasets are achieved with using Lm in the loss function. We can observe that Lc is also important for GED learning. On LINUX and IMDB, the best results are achieved with incorporating Lc in the loss function.

We also examine the effectiveness of combining Lqip and Lcr. On AIDS, the MSE of using Ld + Lm + Lqip is 1.56 and the MSE of using Ld +Lm +Lcr is 1.61. It verifies that using either Lqip or Lcr is worse than combining them in the loss function.

5.3 Generalization Analysis

We answer Q3 in this experiment. We focus on AIDS and examine the graphs whose node number is more than 30 and no more than 50, denoted by G30 50, due to the time cost of computing the exact GEDs. G≤30 and G30 50 have very different distributions of GEDs: the mean of GEDs of G≤30 and G30 50 is 23.57 and 47.22; the standard deviation of GEDs of G≤30 and G30 50 is 7.83 and 14.17, respectively.

For G30 50, we use 6k pairs of graphs in G≤30 as the training data and 2k and 2k pairs of graphs in G30 50 as the validation and test data, respectively. The generalization performance is shown in Table 5. Comparing with Table 2, we can observe that MSE increases with the growth of graph size. It is reasonable as the distributions of GEDs are different. Our method outperforms existing methods with a small time cost.

6 Conclusion

In this paper, we propose a deep learning-based graph edit distance (GED) computation method. A novel GED-specific loss function is proposed to model the two inherent characteristics of GED. Specifically, we propose a partial node matching-based regularizer to encode multiple optimum node matchings. A partial node matching is used as the groundtruth and all the optimum node matchings covered by the partial node matching are encoded. We can also tune the number of optimum node matchings included to balance the overkill and restraint of the regularizer. We also propose a plane intersection-based regularizer to encode the one-to-one constraints for the node matchings included by the partial node matching-based regularizer. We build an association graph of the two input graphs and use the graph neural network on the association graph to learn the cross-graph representation. Our experiments show that our method is effective and significantly outperforms the state-of-the-art GED learning methods and the heuristic GED algorithms on real-world benchmark graphs.

Acknowledgements

This research was supported in part by the GRFs 12201518, 12201119, R2002-20F, 12200819 and CRF C6030-18GF.