An Efficient Probabilistic Approach for Graph
Similarity Search [Extended Technical Report]
Zijian Li #1 , Xun Jian #2 , Xiang Lian ∗3 , Lei Chen #4
#
1

Computer Science and Engineering Department, HKUST, Hong Kong, China

zlicb@cse.ust.hk
∗

2

xjian@cse.ust.hk

leichen@cse.ust.hk

Computer Science Department, Kent State University, Kent, USA
3

arXiv:1706.05476v2 [cs.DB] 23 Jan 2018

4

xlian@kent.edu

Abstract—Graph similarity search is a common and fundamental operation in graph databases. One of the most popular graph
similarity measures is the Graph Edit Distance (GED) mainly because of its broad applicability and high interpretability. Despite
its prevalence, exact GED computation is proved to be NP-hard,
which could result in unsatisfactory computational efficiency
on large graphs. However, exactly accurate search results are
usually unnecessary for real-world applications especially when
the responsiveness is far more important than the accuracy.
Thus, in this paper, we propose a novel probabilistic approach
to efficiently estimate GED, which is further leveraged for the
graph similarity search. Specifically, we first take branches as
elementary structures in graphs, and introduce a novel graph
similarity measure by comparing branches between graphs, i.e.,
Graph Branch Distance (GBD), which can be efficiently calculated in polynomial time. Then, we formulate the relationship
between GED and GBD by considering branch variations as
the result ascribed to graph edit operations, and model this
process by probabilistic approaches. By applying our model,
the GED between any two graphs can be efficiently estimated
by their GBD, and these estimations are finally utilized in the
graph similarity search. Extensive experiments show that our
approach has better accuracy, efficiency and scalability than
other comparable methods in the graph similarity search over
real and synthetic data sets.

I. I NTRODUCTION
Graph similarity search is a common and fundamental operation in graph databases, which has widespread applications
in various fields including bio-informatics, sociology, and
chemical analysis, over the past few decades. For evaluating
the similarity between graphs, Graph Edit Distance (GED)
[1] is one of the most prevalent measures because of its
wide applicability, that is, GED is capable of dealing with
various kinds of graphs including directed and undirected
graphs, labeled and unlabeled graphs, as well as simple graphs
and multi-graphs (which could have multiple edges between
two vertices). Furthermore, GED has high interpretability,
since it corresponds to some sequences of concrete graph edit
operations (including insertion of vertices and edges, etc.) of
minimal lengths, rather than implicit graph embedding utilized
in spectral [2] or kernel [3] measures. Example 1 illustrates
the basic idea of GED.
Example 1: Assume that we have two graphs G1 and G2
as shown in Figure 1. The label sets of graph G1 ’s vertices
and edges are {A, B, C} and {y, y, z}, respectively, and the
label sets of graph G2 ’s vertices and edges are {A, B, C}

and {x, y, z}, respectively. The Graph Edit Distance (GED)
between G1 and G2 is the minimal number of graph edit
operations to transform G1 into G2 . It can be proved that the
GED between G1 and G2 is 3, which can be achieved by 1
deleting the edge between v1 and v3 in G1 , and 2 inserting
an isolated vertex v4 with label A, and 3 inserting an edge
between v3 and v4 with label x.
With GED as the graph similarity measure, the graph
similarity search problem is formally stated as follows.
Problem Statement: (Graph Similarity Search) Given a
graph database D, a query graph Q, and a similarity threshold
τ̂ , the graph similarity search is to find a set of graphs D0 ⊆
D, where the graph edit distance (GED) between Q and each
graph in D0 is less than or equal to τ̂ .
A straightforward solution to the problem above is to check
exact GEDs for all pairs of Q and graphs in database D.
However, despite its prevalence, GED is proved to be NP-hard
for exact calculations [4], which may lead to unsatisfactory
computational efficiency when we conduct a similarity search
over large graphs. The most widely-applied approach for
computing exact GED is the A∗ algorithm [5], which aims to
search out the optimal matching between the vertices of two
graphs in a heuristic manner. Specifically, given two graphs
with n and m vertices, respectively, the time complexity of
A∗ algorithm is O(nm ) in the worst case.
Due to the hardness of computing exact GED (NP-hard)
[4], most existing works involving exact GED computation
[4] [5] only conducted experiments on graphs with less than
10 vertices. In addition, a recent study [6] indicates that the A∗
algorithm is incapable of computing GED between graphs with
more than 12 vertices, which can hardly satisfy the demand for
searching real-world graphs. For instance, a common requirement in bio-informatics is to search and compare molecular
structures of similar proteins [7]. However, the structures of
human proteins usually contain hundreds of vertices (i.e.,
atoms) [8], which obviously makes similarity search beyond
the capability of the approaches mentioned above. Another
observation is that many real-world applications do not always
require an exact similarity value, and an approximate one
with some quality guarantee is also acceptable especially in
real-time applications where the responsiveness is far more
important than the accuracy. Taking the protein search as an

v1
y
v2

C

v1

A
y
z
G1

y
B

v3 v2

C

v1

A
z

×

A v4

v1

y
B

1 Deleting edge hv1 , v3 i

Fig. 1.

A

v3 v2

C

A

y
z

B

2 Inserting vertex v4

v3 v2

C

A v4 u1
x

z

B

3 Inserting edge hv3 , v4 i

x
v3 u3

A u2

B
z

A

y
C

u4

G2

Graph edit process to transform G1 into G2 (Example 1)

example again, it is certainly more desirable for users to obtain
an approximate solution within a second, rather than to wait
for a couple of days to get the exact answer.
To address the problems above, many approaches have been
proposed to achieve an approximate GED between the query
graph Q and each graph in database D in polynomial time
[9], which can be leveraged to accelerate the graph similarity
search by trading accuracy for efficiency. Assuming that there
are two graphs with n and m vertices, respectively, where n ≥
m, one well-studied method [10] [11] for estimating the GED
between these two graphs is to solve a corresponding linear
sum assignment problem, which requires at least O(n3 ) time
for obtaining the global optimal value or O(n2 log n2 ) time for
a local optimal value by applying the greedy algorithm [12].
An alternative method is spectral seriation [13], which first
generates the vector representations of graphs by extracting
their leading eigenvalues of the adjacency matrix (O(n2 ) time)
[14], and then exploits a probabilistic model based on these
vectors to estimate GED in O(nm2 ) time.
To further enhance the efficiency of GED estimation and
better satisfy the demands for graph similarity search on large
graphs, we propose a novel probabilistic approach which aims
at estimating the GED with less time cost O(nd+τ̂ 3 ), where n
is the number of vertices, d is the average degree of the graphs
involved, and τ̂ is the similarity threshold in the stated graph
similarity search problem. Note that the similarity threshold
τ̂ is often set as a small value (i.e., τ̂ ≤ 10) and does not
increase with the number of vertices n in previous studies [4]
[15], thus, we can assume that τ̂ is a constant with regard to n
when the graph is sufficiently large. Moreover, most real-world
graphs studied in related works [11] [12] are scale-free graphs
[16], and we prove that the average degree d = O(log n) for
scale-free graphs. Therefore, under the assumptions above, the
time complexity of our approach is O(nd + τ̂ 3 ) = O(n log n)
for comparing two scale-free graphs, and O(|D|n log n) for
searching similar graphs in the graph database D, where |D|
is the number of graphs in database D.
Our method is mainly inspired by probabilistic modeling
approaches which are broadly utilized in similarity searches
on various types of data such as text and images [17]. The
basic idea of these methods is to model the formation of an
object as a stochastic process, and to establish probabilistic
connections between objects. In this paper, we follow this idea
and model the formation of graph branch distances (GBDs) as
the results of random graph editing processes, where GBD is
defined in Section III. By doing so, we prove the probabilistic
relationship between GED and GBD, which is finally utilized
to estimate GED by graph branch distance (GBD).
To summarize, we make the following contributions.
• We adopt branches [15] as elementary structures in graphs,
and define a novel graph similarity measure by compar-

ing branches between graphs, i.e., Graph Branch Distance
(GBD), which can be efficiently calculated in O(nd) time,
where n is the number of vertices, and d is the average
degree of the graphs involved.
• We build a probabilistic model which reveals the relationship between GED and GBD by considering branch
variations as the result ascribed to graph edit operations. By
applying our model, the GED between any two graphs can
be estimated by their GBD in O(τ̂ 3 ) time, where τ̂ is the
similarity threshold in the graph similarity search problem.
• We conduct extensive experiments to show our approach has
better accuracy, efficiency and scalability compared with the
related approximation methods over real and synthetic data.
The paper is organized as follows. In Section II, we formally
define the symbols and concepts which are used in this paper.
In Section III, we give definitions of branches and Graph
Branch Distance (GBD). In Section IV, we introduce the
extended graphs, which are exploited to simplify our model.
In Section V, we derive the probabilistic relation between
GBD and GED, which is leveraged in Section VI to perform
the graph similarity search. In Section VII, we demonstrate
the efficiency and effectiveness of our proposed approaches
through extensive experiments. We discuss related works in
Section VIII, and we conclude the paper in Section IX.
II. P RELIMINARIES
The graphs discussed in this paper are restricted to simple
labeled undirected graphs. Specifically, the i-th graph in
database D is denoted by: Gi , {Vi , Ei , L}, where Vi ,
{vi,1 , vi,2 , ..., vi,|Vi | } is the set of vertices, Ei , {ei,1 , ei,2 , ...,
ei,|Ei | } is the set of edges, while L is a general labelling
function. For any vertex vi,j ∈ Vi , its label is given by L(vi,j ).
Similarly, for any edge ei,j ∈ Ei , its label is given by L(ei,j ).
In addition, LV and LE are defined as the sets of all possible
labels for vertices and edges, respectively. We also define ε as
a virtual label, which will be used later in our approach. When
the label of a vertex (or edge) is ε, the vertex (or edge) is said
to be virtual and does not actually exist. Particularly, we have
ε∈
/ LV and ε ∈
/ LE . Note that our method can also handle
directed and weighted graphs by considering edge directions
and weights as special labels.
In this paper, we take Graph Edit Distance (GED) [1] as
the graph similarity measure, which is defined as follows.
Definition 1 (Graph Edit Distance): The edit distance between graphs G1 and G2 , denoted by GED(G1 , G2 ), is the
minimal number of graph edit operations which are necessary
to transform G1 into G2 , where the graph edit operations
(GEO) are restricted to the following six types:
• AV: Add one isolated vertex with a non-virtual label;
• DV: Delete one isolated vertex;
• RV: Relabel one vertex;
• AE: Add one edge with a non-virtual label;

DE: Delete one edge;
RE: Relabel one edge.
Assume that a particular graph edit operation sequence
which transforms graph G1 into G2 is seqi , where i is the
unique ID of this sequence. Then, according to Definition 1,
GED(G1 , G2 ) is the minimal length for all possible operation
sequences, that is, GED(G1 , G2 ) = mini {|seqi |}, where
|seqi | is the length of the sequence seqi . The set of all
operation sequences from G1 to G2 of the minimal length
is defined as SEQ = {seqi |∀i, GED(G1 , G2 ) = |seqi |}.
III. B RANCH D ISTANCE B ETWEEN G RAPHS
To reduce the high cost of exact GED computations (NPhard [4]) in the graph similarity search, one widely-applied
strategy for pruning search results [4] [18] [19] [15] is to
exploit the differences between graph sub-structures as the
bounds of exact GED values. In this paper, we consider the
branches [15] as elementary graph units, which are defined as:
Definition 2 (Branches): The branch rooted at vertex v is
defined as B(v) = {L(v), N (v)}, where L(v) is the label of
vertex v, and N (v) is the sorted multiset containing all labels
of edges adjacent to v. The sorted multiset of all branches in
Gi is denoted by BGi = {B(vi,j ), ∀vi,j ∈ Vi }.
In practice, each branch B(v) is stored as a list of strings
whose first element is L(v) and the following elements are
strings in the sorted multiset N (v). In addition, BGi for each
graph Gi is stored in a sorted multiset, whose elements are
essentially lists of strings (i.e., branches) and are always sorted
ascendingly by the ordering algorithm in [20]. For a fair
comparison of the computational efficiency, we assume that
all auxiliary data structures in different methods, such as the
cost matrix [11], adjacency matrix [13], and our branches, are
pre-computed and stored with graphs.
To define the equality between two branches, we introduce
the concept of branch isomorphism as follows.
Definition 3 (Branch Isomorphism): Two branches B(v) =
{L(v), N (v)} and B(u) = {L(u), N (u)} are isomorphic if
and only if L(v) = L(u) and N (v) = N (u), which is denoted
by B(v) ' B(u).
Suppose that we have two branches B(v) and B(u) where
the degrees of vertices v and u are dv and du , respectively.
From previous discussions, the branch B(v) and B(u) are
stored as lists of lengths dv + 1 and du + 1, respectively.
Therefore, checking whether B(v) and B(u) are isomorphic
is essentially judging whether two lists of lengths dv + 1 and
du + 1 are identical, which can be done in dv time when
dv = du , and otherwise in one unit time since the length of a
list can be obtained in one unit time.
Finally, we define the Graph Branch Distance (GBD).
Definition 4 (Graph Branch Distance): The branch distance between graphs G1 and G2 , denoted by GBD(G1 , G2 ),
is defined as:
GBD(G1 , G2 ) = max{|BG1 |, |BG2 |} − |BG1 ∩ BG2 |
•

•

= max{|V1 |, |V2 |} − |BG1 ∩ BG2 |
(1)
where BG1 and BG2 are the multisets of all branches in graphs
G1 and G2 , respectively.

D
Gi
Vi
Ei
Q
L
LV
LE
ε
τ̂
seq
SEQ
GED
B(v)
L(v)
N (v)
BGi
GBD
G{k}
G01 , G02

,
,
,
,
,
,
,
,
,
,
,
,
,
,
,
,
,
,
,
,

TABLE I
TABLE OF N OTATIONS
{G1 , G2 , ...G|D| }, the graph database
{Vi , Ei , L}, i-th graph in database
{vi,1 , vi,2 , ..., vi,|Vi | }, the vertices in Gi
{ei,1 , ei,2 , ..., ei,|Ei | }, the edges in Gi
{VQ , EQ , L}, the query graph
labelling function for vertices and edges
the set of all possible vertex labels
the set of all possible edge labels
the virtual label
the similarity threshold
a graph edit operation (GEO) sequence
set of all GEO sequences of minimal lengths
Graph Edit Distance
the branch rooted at vertex v
the label of vertex v
sorted multiset of labels of edges adjacent to v
sorted multiset of all branches in graph Gi
Graph Branch Distance
extended graph of G with extension factor k
{0}
{|V |−|V1 |}
, G2 when |V1 | ≤ |V2 |
G1 2

The intuition of introducing GBD is as follows. The state-ofthe-art method [15] for pruning graph similarity search results
assumes that the difference between branches of two graphs
has a close relation to their GED. Therefore, in this paper, we
aim to use GBD to closely estimate the GED of two graphs.
Example 2 below illustrates the process of computing GBD.
Example 2: Assume that we have two graphs G1 and G2
as shown in Figure 1. The branches rooted at the vertices in
graphs G1 and G2 are as follows.
B(v1 ) = {A; y, y}, B(v2 ) = {C; y, z}, B(v3 ) = {B; y, z};
B(u1 ) = {B; x, z}, B(u2 ) = {A; y};
B(u3 ) = {A; x}, B(u4 ) = {C; y, z}.

The sorted multisets of branches in G1 and G2 are:
BG1 = {B(v1 ), B(v3 ), B(v2 )};
BG2 = {B(u3 ), B(u2 ), B(u1 ), B(u4 )}.

Therefore, according to Definition 4, we can obtain the
graph branch distance (GBD) between graphs G1 and G2 by
applying the Equation (1), which is:
GBD(G1 , G2 ) = max{|BG1 |, |BG2 |} − |BG1 ∩ BG2 | = 3,
where |BG1 | = 3 and |BG2 | = 4 are the sizes of multisets BG1
and BG2 , respectively. In addition, according to Definition 3,
the only pair of isomorphic branches between multisets BG1
and BG2 is B(v2 ) ' B(u4 ). Therefore, the intersection set
of multisets BG1 and BG2 is {B(v2 )}, whose size is |BG1 ∩
BG2 | = |{B(v2 )}| = 1.
Note that the time cost of computing the size of intersection
of two sorted multisets is max{m1 , m2 } [4], where m1 and
m2 are the sizes of two multisets, respectively. Therefore, the
GBD between query graph Q and any graph G ∈ D can be
computed in time:
Pn
(2)
i di = O(nd),
where n = max{|VQ |, |VG |}, di is the degree of i-th compared
vertex in G, and d is the average degree of graph G.
The GBD defined in this section is utilized to model the
graph edit process and further leveraged for estimating the
graph edit distance (GED) in Section V.

X

v1

ε

A

v4

u1

A u2

B

GED

y
v2

y
C

x
z
{1}

G1

Fig. 2.

B

v3

u3

S

R

GBD

y

z

A

C

Y

u4
Fig. 3.

{0}

G2

Extended Graphs for Example 3

IV. E XTENDED G RAPHS
In this section, we reduce the number of graph edit operation
types that need to be considered by extending graphs with
virtual vertices and edges, which helps to simplify our probabilistic model in Section V. Moreover, we show that the GED
and the GBD between the extended graphs stay the same as
the GED and the GBD between the original ones, respectively.
The definition of extended graphs is as follows.
Definition 5 (Extended Graphs): For graph G, its extended
graph, denoted by G{k} , is generated by first inserting k isolated virtual vertices into G, and then inserting a virtual edge
between each pair of non-adjacent vertices in G (including
virtual vertices), where k is the extension factor.
Example 3: For graphs G1 and G2 in Figure 1, their
{1}
{0}
extended graphs G1 and G2 are shown in Figure 2. The
virtual vertices are labelled by ε, while virtual edges are
represented by dashed lines. Note that when the extension
factor is 0, no virtual vertex will be inserted.
In particular, for any graph pair (G1 , G2 ) where |V1 | ≤ |V2 |,
{|V |−|V1 |}
{0}
we define G01 = G1 2
and G02 = G2 by extending G1
and G2 with extension factor |V2 | − |V1 | and 0, respectively.
Previous studies [21] [22] have shown that, for any graph edit
operation sequence seq which transforms the extended graph
G01 into G02 and has the minimal length, every operation in seq
is equivalent to a relabelling operation on G01 . Therefore, we
only need to consider graph edit operations of types RV and
RE when modeling the graph edit process of transforming the
extended graph G01 into G02 .
Finally, given the graphs G1 and G2 (for |V1 | ≤ |V2 |),
and their extended graphs G01 and G02 , we have the following
Theorems 1 and 2, which are utilized in the Section V.
Theorem 1: GED(G1 , G2 ) = GED(G01 , G02 )
Proof: Please refer to Appendix A.
Theorem 2: GBD(G1 , G2 ) = GBD(G01 , G02 )
Proof: Please refer to Appendix B.
Note that the extended graph is only a conceptual model for
reducing the number of graph edit operation types that need
to be considered. According to Theorems 1 and 2, whenever
the values of GED(G01 , G02 ) and GBD(G01 , G02 ) are required,
we can instead calculate GED(G1 , G2 ) and GBD(G1 , G2 ).
Therefore, in practice, we do not actually convert graphs into
their extended versions, and the calculations of GEDs and
GBDs are still conducted on original graphs rather than the
extended ones, which means that there is no overhead for
creating and maintaining extended graphs.
V. O UR P ROBABILISTIC M ODEL
In this section, we aim to solve the stated graph similarity search problem by estimating GED from GBD in a
probabilistic manner. According to Theorems 1 and 2, the
relation between original graphs’ GED and GBD must be

S
X
Y
Z
R

,
,
,
,
,

Z

Bayesian Network of Random Variables

TABLE II
N OTATION OF R ANDOM VARIABLES
Choice of operation sequence
Number of relabeled vertices
Number of relabeled edges
Number of vertices covered by relabeled edges
Number of vertices either relabeled or covered by relabeled edges

the same with the relation between extended graphs’ GED
and GBD. Therefore, as discussed in Section IV, we only
need to consider graph edit operations of types RV and RE
when building our probabilistic model.
To be more specific, we consider two given graphs G1 and
G2 where |V1 | ≤ |V2 | and their extended graphs G01 and G02 .
For simplicity, we denote GED(G01 , G02 ) and GBD(G01 , G02 )
by GED and GBD, respectively. As mentioned in Section I,
we model the formation of graph branch distances (GBDs) as
the results of random graph editing processes, and thus we
establish a probabilistic connection between GED and GBD.
The detailed steps to construct our model is as follows.
Step 1: We consider GED as an observed random variable.
Step 2: We randomly choose one graph edit operation (GEO)
sequence from all possible GEO sequences whose lengths are
equal to GED. We define a random variable X(ω), where ω
is a particular choice of operation sequence from SEQ, and
S(ω) = s iff ω choose the sequence with ID s, that is, seqs .
Step 3: We model the numbers of RV and RE operations
in the sequence chosen in Step 2 as random variables X(s)
and Y (s), respectively, where X(s) = x iff the number of
operations with type RV in seqs is x, and Y (s) = y iff the
number of operations with type RE in seqs is y. Note that,
when given GED = τ and X = x, we always have Y = τ −x.
Step 4: We define random variables Z(y) as the random
variable where y is a particular value of Y , and Z(y) = m iff
Y (s) = y and the number of vertices covered by relabeled
edges is m when conducting seqs . In addition, we define
R(x, m) as the random variable where x and m are particular
values of X and Z, respectively. That is, R(x, m) = r iff
X = x, Z = m and the number of vertices either relabelled
or covered by relabelled edges is r.
The reason we conduct Step 4 is because we want to model
branch variations by the number of branches rooted at the
vertices either relabelled or covered by relabelled edges, i.e.,
the random variable R.
Step 5: We consider GBD as the random variable dependent
on R, where their relation is proved in Appendix G.
The random variables defined in the five-step model above
are listed in Table II, and their relations among are represented
by a Bayesian Network, as shown in Figure 3. We use Example
4 below to better illustrate the key idea of our model.
Example 4: Given two extended graphs G01 and G02 , as
shown in Figure 4, where the virtual edges are represented by
dashed lines. The minimal number of graph edit operations to
transform G01 into G02 is 2, and the set of all possible graph

0.21

y

x
v2

u1

A

B

y
C

G01

Fig. 4.

A

v3

u2

Sampled
Inferred

0.18
0.15

x

B

C
G02

Extended Graphs for Example 4

u3

Frequency

v1

0.12
0.09
0.06
0.03
0

0

Fig. 5.

5
10
15
Graph Branch Distances

Inferred prior distribution

Fig. 6.

Inferred prior distribution

edit operation sequences with the minimal length 2 is:
of GBDs on Fingerprint data set
of GEDs on Fingerprint data set
SEQ = {{op1 , op2 }1 , {op2 , op1 }2 , {op3 , op4 }3 , {op4 , op3 }4 }
model, where the expanded formula is:
P
P
P
where the subscript of each sequence is its ID, and
Λ 1 = x Ω1 m Ω2 r Ω 3 · Ω4
(8)
op1 = Relabelling the edge hv1 , v2 i to label y;
where
P
op2 = Relabelling the edge hv1 , v3 i to label x;
(9)
Ω1 = N1 s P r[X = x, Y = τ − x | S = s]
op3 = Relabelling the vertex v2 to label C.
Ω2 = P r [Z = m | Y = τ − x]
(10)
op4 = Relabelling the vertex v3 to label B.
Ω3 = P r [GBD = ϕ | R = r]
(11)
Then, the values of random variables in our model could be:
Ω4 = P r [R = r | X = x, Z = m]
(12)
1) First, we consider GED as a random variable, which has
Please refer to Appendix C for the closed forms of Ω1 ,
the value 2 in this example. (GED = 2)
Ω2 , Ω3 , and Ω4 . Due to the space limitation, please refer to
2) Second, we randomly choose one sequence from SEQ, Appendix D for the proof of Equation (8).
which is seq2 = {op2 , op1 }. Therefore, in this case the B. Calculating Λ2
random variable S = 2.
In this subsection, we aim to calculate Λ2 , which is essen3) Third, the numbers of RV and RE operations in seq2 are 0 tially to infer the prior distributions of GBDs. In practice, we
and 2, respectively. Therefore, the random variables X = 0 pre-compute the prior distribution of GBDs without knowing
and Y = 2 in this example.
the query graph Q in the graph similarity search problem. This
4) Then, after conducting operations in seq2 , the number of is because in most real-world scenarios [4] [15], the query
vertices covered by relabelled edges is 3, and the number graph Q ofter comes from the same population as graphs in
of vertices either relabelled or covered by relabelled edges database D. As a counter example, it is unusual to use a query
is 3. Therefore, the random variables Z = 3 and R = 3.
graph of protein structure to search for similar graphs in social
5) Finally, GBD is considered to be the random variable, networks, and vice versa. Therefore, we assume that GBDs
where GBD = 2 in this example.
between Q and each graph G in database D, follow the same
In this paper, we aim to infer the probability distribution of prior distributions as those among all graph pairs in D.
To calculate the prior distribution of GBDs, we first ranGED when given GBD, which is essentially to calculate the
domly sample α% of graph pairs from the database D, and
following probability for given constants τ̂ and ϕ:
calculate the GBD between each pair of sampled graphs.
P r[GED ≤ τ̂ | GBD = ϕ]
Then, we utilize the Gaussian Mixture Model (GMM) [25]
Pτ̂
= τ =0 P r[GED = τ | GBD = ϕ]
(3) to approximate the distribution of GBDs between all pairs of
By applying Bayes’ Rule, we have:
sampled graphs, whose probability density function is:
Pτ̂
PK
Λ3
P r[GED ≤ τ̂ | GBD = ϕ] = τ =0 Λ1 · Λ
,
(4)
f (φ) = i=1 πi · N (φ; µi , σi )
(13)
2
where
where K is the number of components in the mixture model
Λ1 (G01 , G02 ; τ, ϕ) = P r[GBD = ϕ | GED = τ ]
(5) defined by user, and N is the probability density function of
Λ2 (G01 , G02 ; ϕ) = P r[GBD = ϕ]
(6) the normal distribution. Here, πi , µi and σi are parameters
of the i-th component in the mixture model, which can be
0
0
Λ3 (G1 , G2 ; τ ) = P r[GED = τ ]
(7)
inferred from the GBDs over sampled graph pairs. Please refer
Therefore, the problem to solve becomes calculating the to [25] for the process of inferring parameters in GMM.
values of Λ1 , Λ2 and Λ3 when given extended graphs G01 and
Finally, we can compute the prior probability P r[GBD =
G02 , and constants τ and ϕ. In this following three subsections ϕ] by integrating the probability density function f (φ) on the
V-A, V-B and V-C, we illustrate how to utilize our model to adjacent interval of ϕ, i.e., [ϕ − 0.5, ϕ + 0.5].
R ϕ+0.5 PK
calculate Λ1 , Λ2 and Λ3 , respectively.
P r[GBD = ϕ] = ϕ−0.5 i=1 πi · N (φ; µi , σi ) dφ (14)
Note that this integration technique is commonly used in
A. Calculating Λ1
In this subsection, we aim to calculate Λ1 = P r[GBD = the field called continuity correction [26], which aims to
ϕ | GED = τ ]. The key idea to calculate Λ1 is to expand approximate discrete distributions by continuous distributions.
its formula by applying Chain Rule [24] according to the In addition, the choice of integral interval [ϕ − 0.5, ϕ + 0.5]
dependency relationship between random variables in our is a common practice in the continuity correction field [27].

Example 5: We randomly sample 60,000 pairs of graphs
from Fingerprint dataset of IAM Graph Database [28], where
the distribution of GBDs between all pairs of sampled graphs
is represented by the blue histogram in Figure 5. We infer
the GBD prior distribution of sampled graph pairs, which is
represented by the red line in Figure 5. This way, we can
compute and store the probability P r[GBD = ϕ] for each
possible value of variable ϕ by Equation (14), where the range
of variable ϕ is analyzed in Section VI.
C. Calculating Λ3
In this subsection, we focus on calculating Λ3 , i.e., the prior
distribution of GEDs. Recall that the GED computation is NPhard [4]. Therefore, sampling some graph pairs and calculating
the GED between each pair is infeasible especially for large
graphs, which means that we cannot simply infer the prior
distribution of GEDs from sampled graph pairs.
To address this problem, we utilize the Jeffreys prior [29]
as the prior distribution of GEDs. The Jeffreys prior is wellknown for its non-informative property [30], which means that
it provides very little additional information to the probabilistic model. Therefore, Jeffreys prior is a common choice in
Bayesian methods when we know little about the actual prior
distribution. Then, according to the definition of Jeffreys prior
[29], the value of probability P r[GED = τ ] is calculated by:
P r[GED = τ ]
qP
2τ
0
0
2
= C1
ϕ=0 P r[GBD = ϕ|GED = τ ] · Z (G1 , G2 ; τ, ϕ) (15)
qP
2τ
0
0
0
0
2
(16)
= C1
ϕ=0 Λ1 (G1 , G2 ; τ, ϕ) · Z (G1 , G2 ; τ, ϕ),
where Equation (15) comes from the definition of the Jeffreys
prior [29], which is the expected value of function Z 2 with
respect to the conditional distribution of variable GBD when
given GED = τ . Here, C is a constant for normalization, and
the function Z is defined in the following Equation (17).
∂{log P r[GBD|GED]}
Z(G01 , G02 ; τ, ϕ) =
(17)
∂{GED}
GED=τ,GBD=ϕ
where the symbol ∂ represents the partial derivative of a
function, and the symbol F (x)|x=k means to substitute value
k for variable x in the function F (x). Please refer to Appendix
C for the closed forms of Equation (16).
According to the closed form of Equation (16), we find
that the value of probability P r[GED = τ ] only depends on
the constant τ and the size of the extended graph G01 , i.e.,
(|V10 |). Therefore, for each data set, we pre-compute the value
of function P r[GED = τ ] for each possible value of τ and
|V10 |, and store these pre-computed values in a matrix for quick
looking up when searching similar graphs. In Section VI, we
discuss the ranges of τ and |V10 | in detail. Note that, if there
are k1 possible values of τ and k2 possible values of |V10 |,
then the normalization constant C in Equations (15) and (16)
will be C = 1/(k1 · k2 ).
Example 6: We show part of the Jeffreys prior distribution
of GEDs on Fingerprint data set [28] in Figure 6 as an
example, where the x-axis represents values of |V10 |, and
the y-axis represents values of τ . The gray scale of each

Algorithm 1 Graph Similarity Search with Graph Branch
Distance Approximation (GBDA)
Input: a query graph Q, a graph database D,
a similarity threshold τ̂ , and a probability threshold γ
Output: the search result D0
for each graph G ∈ D do
Step 1*: Pre-compute Λ2 and Λ3 for all possible inputs
Step 2: Calculate GBD(Q0 , G0 ) by Definition 4.
Step 3: Given GBD(Q0 , G0 ) = ϕ, we calculate
Φ =P r[GED(Q, G) ≤ τ̂ | GBD(Q, G) = ϕ]
τ̂
X

Λ3 (Q0 , G0 ; τ )
Λ2 (Q0 , G0 ; ϕ)
τ =0
where Λ1 (Q0 , G0 ; τ, ϕ) is calculated by Equation (8).
Step 4: Insert G into D0 if Φ ≥ γ
end for
=

Λ1 (Q0 , G0 ; τ, ϕ) ·

1 × 1 block in Figure 6 represents the corresponding value
of P r[GED = τ ].
VI. G RAPH S IMILARITY S EARCH WITH THE
P ROBABILISTIC M ODEL
In this section, we first elaborate on our graph similarity
search algorithm (i.e., GBDA) based on the model derived in
the previous section, which consists of two stages: offline preprocessing and online querying. Then, we study the time and
space complexity of these two stages in detail.
A. Graph Similarity Search Algorithm
Given a query graph Q, a graph database D, a similarity
threshold τ̂ , and a probability threshold γ, the search results
D0 is achieved by Algorithm 1, where Q0 and G0 are the
extended graphs of Q and G, respectively. Step 1 tagged with
symbol * is pre-computed in the offline stage, as discussed in
Sections V-B and V-C. We give the Example 7 below to better
illustrate the process of Algorithm 1.
Example 7: Assume that the graph G1 in Figure 1 is the
query graph Q, and G2 in Figure 1 is a graph in database
D. Given the similarity threshold τ̂ = 3 and the probability
threshold γ = 0.8, the process of determining whether G2
should be in the search result D0 is as follows.
1) First, we pre-compute Λ2 (Q0 , G02 ; ϕ) and Λ3 (Q0 , G02 ; τ ) by
inferring the prior distributions of GBDs and GEDs on
database D, respectively. Since this is a simulated example
and there is no concrete database D, we assume that
Λ3 (Q0 , G02 ; τ )/Λ2 (Q0 , G02 ; ϕ) ≡ 0.8 for all possible values
of τ and ϕ.
2) Second, from Example 2, we know GBD(Q, G2 ) = 3.
3) Then, according to Equation (8), we can calculate
Pτ̂
Λ (Q0 ,G0 ;τ )
Φ = τ =0 Λ1 (Q0 , G02 ; τ, ϕ) · Λ23(Q0 ,G02;ϕ)
2

= (0 + 0 + 0.5113 + 0.5631) × 0.8 = 0.8595 > γ = 0.8
4) Therefore, G2 is inserted into the search result D0 .
B. Complexity Analysis of Online Stage
The online querying stage in our approach includes Steps 2,
3 and 4 in Algorithm 1, where Step 4 clearly costs O(1) time
for each graph G. In addition, from the discussions in Section
III, Step 2 costs O(nd) time, where n = max{|VG |, |VQ |},
and d is the average degree of graph G.

In Step 3, since Λ2 and Λ3 have already been pre-computed
in the offline stage (i.e., Step 1), their values can be obtained
in O(1) time for each τ ∈ [0, τ̂ ] and graph G ∈ D.
Now we focus on analyzing the time complexity of computing Λ1 in Step 3 of Algorithm 1. Let the time for computing
Ω1 , Ω2 , Ω3 and Ω4 in Equation (8) be C1 , C2 , C3 and
C4 , respectively. According to Equation (8), the time for
computing Λ1 for
each P
τ ∈ [0, τ̂ ] and
graph G ∈ D is:
P
P
Ω1

Ω3 ·Ω4

Ω2

}|
{
z}|{ z }| { z
(18)
xC1 + xmC2 + xmr(C3 · C4 )
where x, m and r are the summation subscripts in Equation
(8), and the ranges of x, m and r are:
• x ∈ [0, τ ]. Since x is the number of RV operations, x must
not be larger than the number of graph edit operations τ .
• m ∈ [0, 2τ ]. Since m is the number of vertices covered
by relabelled edges given the relabelled edge number Y =
τ − x, and each edge can cover at most two vertices, we
have 0 ≤ m ≤ 2(τ − x) ≤ 2τ .
• r ∈ [0, 3τ ]. Note that r is the number of vertices either
relabelled or covered by relabelled edges when the relabelled vertex number is X = x and the number of vertices
covered by relabelled edges is Z = m. Therefore, we have
0 ≤ r ≤ x + m ≤ τ + 2τ = 3τ .
In addition, according to the closed form of Equation (8)
in Appendix C, it is clear that C1 = C3 = C4 = O(1) and
C2 = O(m) = O(τ ). According to Equation (18), the time of
computing Λ1Pfor eachPτ ∈ [0,P
τ̂ ] and graph G ∈ D is:
Ω2

Ω1

Ω3 ·Ω4

z }| { z }| { z }| {
O(τ ) + O(τ 3 ) + O(τ 3 ) = O(τ 3 )
(19)
Moreover, from the above discussions about the ranges of
summation subscripts, for any τ ∈ [0, τ̂ ), we have:
Pτ̂ P2τ̂
x=0
m=0 Ω2 (m, x, τ̂ )
Pτ̂ P2τ̂
= x=0 m=0 P r[Z = m | Y = τ̂ − x]
(20)
Pτ P2τ
= x=0 m=0 P r[Z = m | Y = τ − x]
Pτ̂ −τ −1 P2τ
+ x=0
m=0 P r[Z = m | Y = τ̂ − x]
Pτ P2τ̂
+ x=0 m=2τ +1 P r[Z = m | Y = τ̂ − x]
Pτ̂ −τ −1 P2τ̂
+ x=0
(21)
m=2τ +1 P r[Z = m | Y = τ̂ − x]
Pτ P2τ
(22)
=f (m, x, τ̂ ) + x=0 m=0 Ω2 (m, x, τ )
where f (m, x, τ̂ ) is sum of last three terms in Equation
(21). Note that Equation (21) is a sum on four disjoint twodimensional intervals whose combination is the sum interval
of Equation (20).
P
Equation (22) means, the value of x,m Ω2 (m, x, τ ) where
τ <P
τ̂ have already been calculated in the process of computing
we can reduce redundant
x,m Ω2 (m, x, τ̂ ). Therefore,P
computations byP
only computing x,m Ω2 (m, x, τ̂ ) once to
obtain values of x,m Ω2 (m, x,P
τ ) for all τ < τ̂ . Similar conclusions can also be derived for x,m,r Ω3 (r, ϕ)·Ω4 (x, r, m),
where the detailed proofs are omitted here.
Therefore,
P
P the time cost ofPStep 3 in Algorithm 1 is:
Ω2

Ω3 ·Ω4

z }| {
O(τ̂ 3 )

Ω1

+

Λ2

z }| { z }| {
3
τ =0 {O(τ ) + O(1)} = O(τ̂ )

Pτ̂

(23)

for each graph G in database D.
Finally, we can obtain Theorem 3.
Theorem 3: The time of the whole online stage is:
Step 3
Step 4
Step 2
z }| { z }| { z }| {
3
O(nd) + O(τ̂ ) + O(1) = O(nd + τ̂ 3 )
(24)
for each graph G in database D, where n = max{|VG |, |VQ |},
d is the average degree of graph G, and τ̂ is the similarity
threshold in the graph similarity search problem.
Proof: Please refer to the discussions above.
Note that the similarity threshold τ̂ is often set as a small
value (i.e., τ̂ ≤ 10) and does not increase with the number of
vertices n in previous studies [4] [15], thus, we can assume that
τ̂ is a constant with regard to n when the graph is sufficiently
large. Moreover, most real-world graphs studied in related
works [11] [12] are scale-free graphs [16], whose average
degrees d = O(log n) as proved in Appendix K.
C. Complexity Analysis of Offline Stage
The offline pre-processing stage in our approach is Step 1
in Algorithm 1, which is essentially to pre-compute the prior
distributions of GEDs and GBDs respectively among all pairs
of graphs involved in the graph similarity search.
1) Complexity Analysis of Computing the Prior Distribution
of GBDs: As discussed in Section 5.1, the prior distributions
of GBDs can be pre-computed by the following four steps:
Step 1.1: Sample N graph pairs from the database D.
Step 1.2: Calculate GBD between each sampled graph pairs.
Step 1.3: Learn the Gaussian Mixture Model (GMM) of the
GBDs between sampled graph pairs.
Step 1.4: Calculate P r[GBD = ϕ] for each possible value
of ϕ by using Equation (14).
It is clear that Step 1.1 costs O(N ) time. From the discussions in Section III, Step 1.2 costs O(N · nd) time, where
n is the maximal number of vertices among the sampled
graphs, and d is the average degree of the sampled graphs.
The learning process of GMM in Step 1.3 costs O(N · K)
time [25], where K is the number of components in GMM,
and  is the maximal learning iterations for learning GMM.
As for Step 1.4, since ϕ is the value of GBD, from the
definition of GBD, the possible values of ϕ are essentially
{0, 1, 2..., n}, where n is the maximal number of vertices
among the sampled graphs. According to Equation (14),
computing P r[GBD = ϕ] for each ϕ costs O(K) time, where
K is the number of components in GMM derived in Step 1.3.
Thus, Step 1.4 costs O(nK) time.
Note that in the Gaussian Mixture Model, the component
number K and the maximal learning iterations  are fixed
constants. Therefore, the prior distributions of GBDs can be
calculated in time
Step 1.1
Step 1.2
Step 1.3
Step 1.4
z }| { z }| { z }| { z }| {
O(N ) + O(N nd) + O(N K) + O(nK) = O(N nd) (25)
where N is the number of graph pairs sampled in Step 1.1,
n is the maximal number of vertices among sampled graphs,
and d is the average degree of sampled graphs.
Based on the discussions above, the number of precomputed values of P r[GBD = ϕ] is at most n. Thus, the
space cost of storing the prior distribution of GBDs is O(n).

Data Set
AIDS
Finger
GREC
AASD
Syn-1
Syn-2

TABLE III
S TATISTICS OF DATA S ETS
|D|
|Q|
Vm
Em
1896
100
95
103
2159
114
26
26
1045
55
24
29
37995
100
93
99
3430
70
100K
1M
3430
70
100K
1M

d
2.1
1.7
2.1
2.1
9.6
9.4

Scale-free
Yes
Yes
Yes
Yes
Yes
No

Note: |D| is the number of graphs in database D. |Q| is the number
of query graphs. Vm and Em are the maximal numbers of vertices and
edges, respectively, while d is the average degree. K means thousand
and M means million.

2) Complexity Analysis of Computing the Prior Distribution
of GEDs: According to the discussions in Section V-C, computing the prior distribution of GEDs is essentially calculating
Equation (16) for each possible values of τ and |V10 |.
First, from the closed form of Equation (16) in Appendix
C, when τ and |V10 | are fixed values, it is clear that computing
d
dτ log Λ1 costs the same time as computing Λ1 , which is
O(τ̂ 3 ), where τ̂ is the user-defined similarity threshold. n and
d are the maximal number of vertices and the average degree
among all graphs, respectively.
Then, since one graph edit operation can at most change two
branches, when the GED between two graphs is τ , the possible
GBD values (i.e., ϕ in Equation (16)) between these two
graphs are {0, 1, 2..., 2τ }. Therefore, according to Equation
(16), the prior probability value of P r[GED = τ ] can be
calculated in time complexity O(2τ · τ 3 ) = O(τ 4 ) when τ
and |V10 | are fixed values.
Finally, recall that computing the GED prior distribution is
essentially calculating Equation (16) for all possible values of
τ and |V10 |, and it is clear that the possible values of τ are
{0, 1, 2..., τ̂ }, where τ̂ is the user-defined similarity threshold. In addition, the possible values of |V10 | are essentially
{1, 2, ...n}, where n is the maximal number of vertices among
all graphs involved the graph similarity search. Therefore, the
time of calculating the GED prior distribution is:
O(τ̂ · n · τ̂ 4 ) = O(nτ̂ 5 )
According to Section V-C, we need to store a matrix whose
rows represent possible values of τ , and columns represent
possible values of |V10 |. Therefore, the space cost of storing
the prior distribution of GEDs is O(τ̂ · n).
Finally, we have Theorem 4.
Theorem 4: The time complexity of the offline stage is:
GBD Prior

GED Prior

z }| { z }| {
O(N nd) + O(nτ̂ 5 ) = O(N nd + nτ̂ 5 )
and the space complexity of the offline stage is:
GBD Prior

GED Prior

z }| {
z }| {
O(n) + O(τ̂ · n) = O(n(1 + τ̂ ))
where N is the number of graph pairs sampled in Step 1.1,
n and d are the maximal number of vertices and the average
degree among all the graphs involved in the graph simialrity
search, respetively. τ̂ is the user-defined similarity threshold.
Proof: Please refer to the discussions above.
VII. E XPERIMENTS
A. Data Sets and Settings
We first present the experimental data sets for evaluating
our approaches. There are 4 real-world data sets (i.e., AIDS,

TABLE IV
C OSTS OF COMPUTING GBD PRIOR DISTRIBUTION
Data Set
AIDS
Finger
GREC
AASD
Syn-1
Time Costs
11.1s
7.5s
20.6s
232.4s
3.8h
Space Costs
0.06kb
0.04kb
0.10kb
1.21kb
13.3gb
TABLE V
C OSTS OF COMPUTING GED PRIOR DISTRIBUTION
Data Set
AIDS
Finger
GREC
AASD
Syn-1
Time Costs
70.32h
16.91h
15.40h
69.16h
6.31h
Space Costs
1.5kb
0.4kb
0.4kb
1.4kb
0.1kb

Syn-2
3.2h
0.3gb

Syn-2
6.31h
0.1kb

Note: h means hours and s means seconds. kb means KBytes and gb
means GBytes.

Fingerprint and GREC from the IAM Graph Database [28],
and AIDS Antiviral Screen Data (AASD) [31]), and 2 synthetic data sets (i.e., Syn-1 and Syn-2). The details about the
data sets are listed in Table III.
The 4 real-world data sets are widely-used for evaluating the
performance of GED estimation methods in previous works
[11] [12] [13]. In order to evaluate how well the GED is
approximated, we must know the exact value of GED, which
is NP-hard to compute [4]. Specifically, even the state-of-theart method [6] cannot compute an exact GED for graphs with
100 vertices within 48 hours on our machine (with 32 Intel
E5 2-core, 2.40 GHz CPUs and 128GB DDR3 RAMs).
However, we still manage to evaluate our proposed method
on large graphs. To address the problem above, we generate
2 sets of large random graphs (i.e., Syn-1 and Syn-2), where
the GED between each pair of graphs is known. Both data
sets Syn-1 and Syn-2 contain 7 subsets of graphs, where each
subset contains 500 graphs whose numbers of vertices are
1K, 2K, 5K, 10K, 20K, 50K, and 100K, respectively. The
difference between data sets Syn-1 and Syn-2 is that the graphs
in Syn-1 satisfy the scale-free property [16] while graphs in
Syn-2 are not. The algorithm of generating synthetic graphs
with known GEDs is described in Appendix I.
Note that, the scale-free property of real data sets is testified
by checking whether the degree distributions of the vertices
in real data sets follow the power-law distribution, while the
scale-free property of Syn-1 data set and the non-scale-free
property of Syn-2 data sets are ensured by our algorithm of
generating synthetic graphs in Appendix I.
For each real data set, we randomly select 5% graphs as
query graphs, while the remaining 95% graphs constitute the
graph database D. For each synthetic data set, we randomly
select 10 graphs from each of its subset as query graphs.
On real data sets, we evaluate our method with the similarity thresholds τ̂ = {1, 2, ..., 10}, which are commonly-used
values of the similarity thresholds in previous studies [4] [15].
On the synthetic data sets, we test our method with larger
similarity thresholds τ̂ = {10, 11, 12, ..., 30} to show that,
when the similarity threshold is larger than the commonlyused values (i.e., {1, 2, ..., 10}), our GBDA method is still
more efficient than the competitors on large graphs.
B. Evaluating Offline Stage
In this subsection, we evaluate the time and space costs of
the offline stage in our GBDA approach, which is essentially
to pre-compute the prior distributions of GEDs and GBDs,
on both real and synthetic data sets. Tables IV and V present
the time and space costs of estimating GBD and GED prior
distributions on different data sets, respectively, where the

number of graph pairs sampled to estimate the GBD prior
distribution is set to N = 100, 000.
The experimental results generally confirm the complexity
analysis in Section VI-C. Specifically, since the number of
sampled graph pairs N and the similarity threshold τ̂ are fixed
values in our experiments, the cost of inferring the GBD prior
distribution grows with n and d, while the cost of estimating
the GED prior distribution depends on n, where n and d are the
maximal number of vertices and the average degree among all
the graphs involved in the graph simialrity search, respetively.
Note that, the costs of computing the GED prior distribution
on synthetic graphs do not exactly follow the theoretical analysis. This is because the numbers of vertices in synthetic graphs
have only 7 possible values (i.e., 1K, 2K, 5K, 10K, 20K, 50K,
and 100K), instead of the worst-case range 1 ∼ 100K as
discussed in Section VI-C. In addition, although the synthetic
graphs are larger than the real ones, the number of possible
values of n on synthetic graphs are smaller than that on real
graphs, where n is the number of vertices in graphs. Therefore,
the costs of computing the GED prior distribution on synthetic
graphs are smaller than the costs on real data sets.
C. Evaluating Online Stage
In this subsection, we compare the efficiency and accuracy
of the online stage of our GBDA approach with three competitors (i.e., LSAP [11], Greedy-Sort-GED [12] and Graph
Seriation [13]), by conducting graph similarity search tasks
over both real and synthetic data sets. In addition, we analyze
how the efficiency and effectiveness (i.e., accuracy, recall and
F1-score) of our method are influenced by the parameters, such
as the similarity threshold τ̂ , the probability threshold γ, and
the number, n, of vertices.
1) The Efficiency Evaluation: In this subsection, we evaluate the query efficiency of our GBDA approach and three
competitors on real and synthetic data sets. Note that the time
cost of our GBDA methods depends on both n and τ̂ , while
the competitors’ time costs only depend on n, where n is the
number of vertices in graphs and τ̂ is the similarity threshold.
Therefore, the experiments in this subsection are conducted
under a fixed probability threshold γ = 0.9, since γ does not
affect the time costs of all the methods.
Specifically, given a specific method and its parameters
(e.g., τ̂ and γ), for each query graph Q, we utilize this method
to obtain a set of graphs similar to graph Q from each data
set, and we record the average query response time for each
data set, which are presented in Figures 7∼9. Particularly, for
a specific method under one specific parameter set (e.g., τ̂ and
γ), each query’s response time is recorded and counted only
once for each data set, and we present the average of response
times for all queries in the experimental figures.
The result in Figure 7 shows that our GBDA approach is
more efficient than the three competitors on all real data sets
where τ̂ is set to 1, 5 and 10, respectively. In addition, we
studied how the number, n, of vertices in graphs influences
the efficiency of our GBDA approach by comparing the query
response time on synthetic data sets with various similarity
thresholds τ̂ , where the results are shown in Figures 8 and 9.

The results show that our GBDA approach is more efficient
than the competitors on both scale-free and non-scale-free
graphs, where the similarity threshold τ̂ ≤ 20. Particularly,
when the similarity threshold τ̂ = 30, although our GBDA
approach costs more time than the other methods on graphs
with 1,000 vertices, our approach is faster than the competitors
on larger graphs with more than 2,000 vertices. Therefore,
when the similarity threshold is larger than the commonlyused values (i.e., {1, 2, ..., 10}), the time cost of our algorithm
is still smaller than the compared methods on large graphs
(i.e., graphs with more than 2,000 vertices).
Note that the competitors (i.e., LSAP, Greedy-Sort-GED and
Graph Seriation) can handle graphs with at most 20K vertices
on our machine. Specifically, when the graphs have more than
20K vertices, all the competitors consume more than 128 GB
memory on our machine, which exceeds the capacity of the
physical memory. However, our GBDA method can handle
graphs with 100K vertices efficiently, which confirms that our
method has better scalability (with respect to the number, n,
of vertices) than the competitors.
2) The Effectiveness Evaluation: We evaluate the effectiveness of our GBDA approach and three competitors by
comparing the precision, recall, and F1-score [32] of the
query results on each real data set with probability thresholds
γ = 0.7, 0.8 and 0.9, respectively.
The results in Figures 10∼13 show that our approach always
outperforms the other three competitors in precision on AIDS,
Fingerprint and GREC data sets, and achieves the highest
precisions among all methods where τ̂ = 1, 2, 3, 4, 5 and 7 on
AASD data set. The results in Figures 14∼17 show that our
method has the second highest recalls under most parameter
settings. Note that, LSAP method returns a lower bound of
GED [11], and therefore the recall of its search result is always
100%. However, if we evaluate the methods by the F1-score,
our method always outperforms the other three competitors
on AIDS, Fingerprint and GREC data sets, and achieves the
highest F1-scores among all methods when τ̂ = 1, 2, 3, 4, 5
and 7 on AASD data set. Therefore, the experimental results
confirm that the effectiveness of our method outperforms the
competitors’ under most parameter settings on the real data
sets in our experiments.
In addition, we study how the number, n, of vertices
in graphs influences the effectiveness (i.e., precision, recall
and F1-score) of our approach on the Syn-1 data set with
various probability thresholds γ and similarity thresholds τ̂ in
Appendix J.
The results in Appendix J show that, the precision of our
method outperforms the competitors on Syn-1 data set, where
the similarity threshold τ̂ = 15, 20, 25 and 30. Moreover,
there is no significant difference between the precision of
our method under various settings of the probability threshold
γ, which demonstrates the robustness of our method under
different settings of parameters. The recalls of our method are
slightly lower than the LSAP method, but much higher than
the Greedy-Sort-GED and Graph Seriation methods, where the
similarity threshold τ̂ = 20, 25 and 30. Finally, the F1-Scores

LSAP
greedysort
seriation

GBDA(τ̂ = 1)
GBDA(τ̂ = 5)
GBDA(τ̂ = 10)

GBDA(τ̂ = 10)
GBDA(τ̂ = 20)
GBDA(τ̂ = 30)

LSAP
greedysort
seriation

GBDA(τ̂ = 10)
GBDA(τ̂ = 20)
GBDA(τ̂ = 30)

LSAP
greedysort
seriation

102

102

101

101

10−7

10

−8

time(second)

time(second)

time(second)

10−6

100
10−1
10−2

10−9

10−2

−3

10−3

10
AIDS Fingerprint GREC

AASD

1K 10K 20K

AASD2

Fig. 8.

Time Costs on Real Graphs
GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

LSAP
greedysort
seriation

LSAP
greedysort
seriation

LSAP
greedysort
seriation

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

0.8

0.7

0.7

0.7

0.7

0.6

0.6

0.6

0.4

1.0
0.9

precision

precision

precision

0.8

0.4

0.5
0.4

0.6
0.5
0.4

0.3

0.3

0.3

0.3

0.2

0.2

0.2

0.2

0.1

0.1

0.1

1

2

3

4

5

6

7

8

9

10

τ̂

Fig. 10.

1

Precision vs. τ̂ on AIDS

3

4

5

6

7

8

9

10

0.1
1

2

3

4

5

τ̂

6

7

8

9

10

1

τ̂

Fig. 11. Precision vs. τ̂ on Fingerprint

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

LSAP
greedysort
seriation

2

Precision vs. τ̂ on GREC

Fig. 12.

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

LSAP
greedysort
seriation

0.9

0.8

0.8

0.7

0.7

0.7

0.7

0.6

0.6

0.6

0.6

0.4

0.4

recall

1.0

0.9

0.8

recall

1.0

0.9

recall

1.0

0.8

0.5
0.4
0.3

0.3

0.2

0.2

0.2

0.2

0.1

0.1

0.1

0.1

2

3

4

5

6

7

8

9

10

1

2

3

4

5

τ̂

6

7

8

9

10

1

2

3

4

5

τ̂

Recall vs. τ̂ on AIDS

LSAP
greedysort
seriation

LSAP
greedysort
seriation

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

Fig. 16.

6

7

8

9

10

1

0.9

0.8

0.8

0.8

0.8

0.7

0.7

0.7

0.7

0.6

0.6

0.6

0.6

0.5

0.5

0.5
0.4

0.4

0.3

0.3

0.3

0.3

0.2

0.2

0.2

0.2

0.1

0.1

0.1

0.1

4

5

6

7

8

9

10

τ̂

1

2

3

4

5

4

5

6

7

8

9

10

1

τ̂

Fig. 19. F1-Score vs. τ̂ on Fingerprint

of our method are mostly higher than the competitors. Therefore, the experimental results on synthetic graphs demonstrate
that our method is more effective than the competitors on large
graphs.
D. Comparing with Alternatives
In this subsection, we compare the effectiveness of our
GBDA approach and two variants of our method by comparing
the F1-score [32] of the query results on each real data set with
probability thresholds γ = 0.9. The variants of our GBDA
method, i.e., GBDA-V1 and GBDA-V2, are considered as
alternatives of GBDA method, which are illustrated as follows.
GBDA-V1: The method GBDA-V1 utilizes the average number of vertices among a sample of graphs from the graph
database as the parameter |V10 | when computing Λ1 and Λ3
in Algorithm 1, instead of using the number of vertices in the
extended query graph Q0 as the parameter |V10 |.
GBDA-V2: The method GBDA-V2 exploits the variant GBD
(VGBD) instead of the original GBD value when computing

6

7

8

9

10

Recall vs. τ̂ on AASD
GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

0.5

0.4

3

10

0.9

0.4

F1-Score vs. τ̂ on AIDS

9

1.0

f1

f1

1.0

0.9

f1

1.0

2

3

LSAP
greedysort
seriation

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

0.9

1

2

Fig. 17.

1.0

Fig. 18.

8

τ̂

Recall vs. τ̂ on GREC

LSAP
greedysort
seriation

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

7

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

τ̂

Recall vs. τ̂ on Fingerprint

Fig. 15.

6

0.4

0.3

1

5

0.5

0.3

Fig. 14.

4

Precision vs. τ̂ on AASD
LSAP
greedysort
seriation

0.9

0.5

3

τ̂

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

LSAP
greedysort
seriation

2

Fig. 13.

1.0

0.5

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

LSAP
greedysort
seriation

GBDA(γ = 0.70)
GBDA(γ = 0.80)
GBDA(γ = 0.90)

0.9

0.5

100K

Time Costs vs. n on Syn-2

0.8

0.5

50K

graph size

Fig. 9.

1.0

0.8

precision

1K 10K 20K

Time Costs vs. n on Syn-1

0.9

0.9

recall

100K

1.0

1.0

f1

50K

graph size

dataset

Fig. 7.

100
10−1

2

3

4

5

6

7

8

9

10

τ̂

Fig. 20.

F1-Score vs. τ̂ on GREC

1

2

3

4

5

6

7

8

9

10

τ̂

Fig. 21.

F1-Score vs. τ̂ on AASD

Λ1 and Λ2 in Algorithm 1, where VGBD is defined as:
V GBD(G1 , G2 ) = max{|V1 |, |V2 |} − w · |BG1 ∩ BG2 | (26)
where BG1 and BG2 are the multisets of all branches in graphs
G1 and G2 , respectively, |V1 | and |V2 | are numbers of vertices
in graphs G1 and G2 , respectively, and w is a user-defined
constant.
Specifically, we denote the number of sampled graphs in
method GBDA-V1 by α, and we test the method GBDA-V1
by setting α = 10, 50, and 100 on real data sets. In addition,
we evaluate the effectiveness of method GBDA-V2 on real
data sets by setting the parameter w = 0.1 and 0.5, where w
is defined in Equation (26).
The results in Figures 22∼25 show that our GBDA method
achieves higher F1-score than its variant GBDA-V1 for the
similarity threshold τ̂ ≤ 4, but generally has the same F1score as GBDA-V1 for τ̂ ≥ 5. In addition, from the results in
Figures 26∼29, we can find that the F1-score of our GBDA
method is higher than or almost the same as the GBDA-V2

V1(α = 50)
V1(α = 100)

V1(α = 50)
V1(α = 100)

GBDA
V1(α = 10)

0.9

0.9

0.9

0.8

0.8

0.8

0.8

0.7

0.7

0.7

0.7

0.6

0.6

0.6

0.6

0.5

0.5

0.9

0.5

0.5

0.4

0.4

0.4

0.4

0.3

0.3

0.3

0.3

0.2

0.2

0.2

0.2

0.1

0.1
2

3

4

5

6

7

8

9

10

0.1

0.1
1

2

3

4

5

τ̂

6

7

8

9

10

1

2

3

4

5

τ̂

Fig. 22. F1-Score vs. τ̂ on AIDS
(Compared with GBDA-V1)
GBDA
V2(w = 0.1)

GBDA
V2(w = 0.1)

V2(w = 0.5)

6

7

8

9

10

1

τ̂

Fig. 23. F1-Score vs. τ̂ on Fingerprint
(Compared with GBDA-V1)

Fig. 24. F1-Score vs. τ̂ on GREC
(Compared with GBDA-V1)
GBDA
V2(w = 0.1)

V2(w = 0.5)

0.9

0.8

0.8

0.8

0.7

0.7

0.7

0.7

0.6

0.6

0.6

0.6

0.5

0.5

0.5
0.4

0.4

0.3

0.3

0.3

0.3

0.2

0.2

0.2

0.2

0.1

0.1

0.1

0.1

6

7

8

9

τ̂

10

1

2

3

4

5

6

7

8

9

10

τ̂

Fig. 27. F1-Score vs. τ̂ on Fingerprint
(Compared with GBDA-V2)

method for the similarity threshold τ̂ ≤ 2 on all real data
sets, and slightly lower than GBDA-V2 method on Fingerprint
data set for τ̂ ≥ 3. In general, the experimental results show
that our GBDA method outperforms methods GBDA-V1 and
GBDA-V2 in most cases for similarity threshold τ̂ ≤ 5, and
performs better or almost the same as GBDA-V1 and GBDAV2 methods for the similarity threshold τ̂ ≥ 6.
VIII. R ELATED W ORKS
A. Exact GED Computation
The state-of-the-art method for exact GED computation is
the A∗ algorithm [5] and its variant [6], whose time costs are
exponential with respect to graph sizes [4]. To address this
NP-hardness problem, most graph similarity search method
are based on the filter-and-verification framework [4] [18] [19]
[15], which first filters out undesirable graphs from the graph
databases and then only verifies the remaining candidates. A
common filtering approach is to use the distance between substructures of two graphs as a lower bound of their GED, which
includes tree-based [18], path-based [33] , branch-based [15]
and partition-based [34] approaches. In this paper, we adopt
the branch structure [15] to build our model. However, we
re-define the distance between branches, since the original
definition [15] of branch distances requires O(n3 ) time for
computation while ours only requires O(nd) time. In addition,
a recent paper [35] propose a multi-layer indexing approach
to accelerate the filtering process based on their proposed
partition-based filtering method.
B. GED Estimation
In this paper, we focus on GED estimation approaches. One
well-studied method [10] [11] is to utilize the solution of a
linear sum assignment problem (LSAP) as an estimation of
GED. The LSAP is an optimization problem which can be
exactly solved by Hungarian method [36], or be approximately
solved by the greedy method [12] and the genetic algorithm
[37]. In our experiment, we compare our GBDA method with
the exact [11] and greedy [12] solutions of LSAP. Since the
exact solution of LSAP defined on two graphs is a lower
bound of their GED [11], the LSAP method can always

7

8

9

10

V2(w = 0.5)

0.5

0.4

5

6

0.9

0.4

4

5

1.0

f1

f1

0.9

0.8

f1

0.9

3

4

GBDA
V2(w = 0.1)

V2(w = 0.5)

1.0

2

3

τ̂

1.0

1

2

Fig. 25. F1-Score vs. τ̂ on AASD
(Compared with GBDA-V1)

1.0

Fig. 26. F1-Score vs. τ̂ on AIDS
(Compared with GBDA-V2)

V1(α = 50)
V1(α = 100)

1.0

f1

f1

1.0

1

f1

GBDA
V1(α = 10)

V1(α = 50)
V1(α = 100)

GBDA
V1(α = 10)

1.0

f1

f1

GBDA
V1(α = 10)
1.0

1

2

3

4

5

6

7

8

9

10

τ̂

Fig. 28. F1-Score vs. τ̂ on GREC
(Compared with GBDA-V2)

1

2

3

4

5

6

7

8

9

10

τ̂

Fig. 29. F1-Score vs. τ̂ on AASD
(Compared with GBDA-V2)

obtain all graphs whose GED to the query graph is no larger
than the similarity threshold and achieve 100% recall in
the similarity search tasks. On the other hand, the GreedySort-GED method [12] solves LSAP approximately and has
no bound to actual GED. However, the Greedy-Sort-GED
method generally achieves better estimations of GEDs [12]
and higher precisions in graph similarity search as shown in
our experiments. Another state-of-the-art approach compared
in our experiment is graph seriation [13], which first converts
graphs into one-dimensional vectors by extracting their leading
eigenvalues of the adjacency matrix, and then exploits a
probabilistic model based on these vectors to estimate the
GED. Although both graph seriation and our approach utilize
probabilistic models, the structure of our model is totally
different from the prior work [13]. In addition, their model
takes the leading eigenvalues of the adjacency matrix as the
inputs, while the inputs of our model are the GBDs. Moreover,
our GBDA method outperforms the competitors’ (i.e., the
LSAP [11], Greedy-Sort-GED [12] and Graph Seriation [13])
under most parameter settings on real data sets.
IX. C ONCLUSIONS
In this paper, we define the branch distance between two
graphs (GBD), and further prove that the GBD has a probabilistic relationship with the GED by considering branch variations as the result ascribed to graph edit operations and modeling this process by probabilistic approaches. Furthermore,
this relation between GED and GBD is leveraged to perform
graph similarity searches. Experimental results demonstrate
both the correctness and effectiveness of our approach, which
outperforms the comparable methods.
R EFERENCES
[1] A. Sanfeliu and K.-S. Fu, “A distance measure between attributed
relational graphs for pattern recognition,” IEEE transactions on systems,
man, and cybernetics, no. 3, pp. 353–362, 1983.
[2] T. Caelli and S. Kosinov, “An eigenspace projection clustering method
for inexact graph matching,” IEEE transactions on pattern analysis and
machine intelligence, vol. 26, no. 4, pp. 515–519, 2004.
[3] T. Gärtner, P. Flach, and S. Wrobel, “On graph kernels: Hardness results
and efficient alternatives,” in Learning Theory and Kernel Machines.
Springer, 2003, pp. 129–143.

[4] Z. Zeng, A. K. Tung, J. Wang, J. Feng, and L. Zhou, “Comparing
stars: on approximating graph edit distance,” Proceedings of the VLDB
Endowment, vol. 2, no. 1, pp. 25–36, 2009.
[5] P. E. Hart, N. J. Nilsson, and B. Raphael, “A formal basis for the heuristic
determination of minimum cost paths,” IEEE transactions on Systems
Science and Cybernetics, vol. 4, no. 2, pp. 100–107, 1968.
[6] K. Gouda and M. Hassaan, “Csi ged: An efficient approach for graph
edit similarity computation,” in 32nd IEEE International Conference on
Data Engineering, ICDE, 2016, pp. 265–276.
[7] R. Ibragimov, M. Malek, J. Guo, and J. Baumbach, “Gedevo: an evolutionary graph edit distance algorithm for biological network alignment,”
in OASIcs-OpenAccess Series in Informatics, vol. 34. Schloss DagstuhlLeibniz-Zentrum fuer Informatik, 2013.
[8] R. P. Ron Milo, Cell Biology by the Numbers. Garland Science, 2015.
[9] X. Gao, B. Xiao, D. Tao, and X. Li, “A survey of graph edit distance,”
Pattern Analysis and applications, vol. 13, no. 1, pp. 113–129, 2010.
[10] S. Bougleux, L. Brun, V. Carletti, P. Foggia, B. Gaüzère, and M. Vento,
“Graph edit distance as a quadratic assignment problem,” Pattern
Recognition Letters, 2016.
[11] K. Riesen and H. Bunke, “Approximate graph edit distance computation
by means of bipartite graph matching,” Image and Vision computing,
vol. 27, no. 7, pp. 950–959, 2009.
[12] K. Riesen, M. Ferrer, and H. Bunke, “Approximate graph edit distance
in quadratic time,” IEEE/ACM transactions on computational biology
and bioinformatics, 2015.
[13] A. Robles-Kelly and E. R. Hancock, “Graph edit distance from spectral
seriation,” IEEE transactions on pattern analysis and machine intelligence, vol. 27, no. 3, pp. 365–378, 2005.
[14] G. H. Golub and C. F. Van Loan, Matrix computations. JHU Press,
2012, vol. 3.
[15] W. Zheng, L. Zou, X. Lian, D. Wang, and D. Zhao, “Graph similarity
search with edit distance constraint in large graph databases,” in CIKM
’13: Proceedings of the 22nd ACM international conference on Information & Knowledge Management, 2013.
[16] Wikipedia, “Scale-free network,” 2017, https://en.wikipedia.org/w/index.
php?title=Scale-free network.
[17] D. J. Hu, “Latent dirichlet allocation for text, images, and music,”
University of California, San Diego. Retrieved April, vol. 26, 2009.
[18] G. Wang, B. Wang, X. Yang, and G. Yu, “Efficiently indexing large
sparse graphs for similarity search,” IEEE Transactions on Knowledge
and Data Engineering, vol. 24, no. 3, pp. 440–451, 2012.
[19] X. Zhao, C. Xiao, X. Lin, and W. Wang, “Efficient graph similarity joins
with edit distance constraints,” in IEEE 28th International Conference
on Data Engineering. IEEE, 2012, pp. 834–845.
[20] C. S. Library, “std::lexicographical compare,” 2017, http://www.
cplusplus.com/reference/algorithm/lexicographical compare/.
[21] D. Justice and A. Hero, “A binary linear programming formulation of
the graph edit distance,” IEEE Transactions on Pattern Analysis and
Machine Intelligence, vol. 28, no. 8, pp. 1200–1214, 2006.
[22] F. Serratosa, “Fast computation of bipartite graph matching,” Pattern
Recognition Letters, vol. 45, pp. 244–250, 2014.
[23] Z. Li, X. Jian, X. Lian, and L. Chen, “An efficient probabilistic approach
for graph similarity search (technical report),” 2017, http://www.cs.kent.
edu/∼xlian/TR/techreport.pdf.
[24] Wikipedia. (2017) Chain rule. [Online]. Available: https://en.wikipedia.
org/w/index.php?title=Chain rule (probability)&oldid=801944140
[25] N. E. Day, “Estimating the components of a mixture of normal distributions,” Biometrika, vol. 56, no. 3, pp. 463–474, 1969.
[26] Wikipedia, “Continuity correction - wikipedia,” 2017. [Online]. Available: https://en.wikipedia.org/w/index.php?title=Continuity correction
[27] M. H. DeGroot and M. J. Schervish, Probability and Statistics: Pearson
New International Edition. Pearson Higher Ed, 2013.
[28] K. Riesen and H. Bunke, “Iam graph database repository for graph based
pattern recognition and machine learning,” Structural, Syntactic, and
Statistical Pattern Recognition, pp. 287–297, 2008.
[29] H. Jeffreys, “An invariant form for the prior probability in estimation
problems,” in Proceedings of the Royal Society of London a: mathematical, physical and engineering sciences, vol. 186, no. 1007. The Royal
Society, 1946, pp. 453–461.
[30] Wikipedia, “Jeffreys prior,” 2017. [Online]. Available: https://en.
wikipedia.org/w/index.php?title=Jeffreys prior&oldid=777352107
[31] D. N. Zaharevitz. (2015) Aids antiviral screen data. [Online]. Available: https://wiki.nci.nih.gov/display/NCIDTPdata/AIDS+
Antiviral+Screen+Data

[32] Wikipedia. (2017) F1 score. [Online]. Available: https://en.wikipedia.
org/w/index.php?title=F1 score&oldid=817504041
[33] X. Zhao, C. Xiao, X. Lin, W. Wang, and Y. Ishikawa, “Efficient
processing of graph similarity queries with edit distance constraints,”
The VLDB Journal, vol. 22, no. 6, pp. 727–752, 2013.
[34] X. Zhao, C. Xiao, X. Lin, Q. Liu, and W. Zhang, “A partition-based
approach to structure similarity search,” Proceedings of the VLDB
Endowment, vol. 7, no. 3, pp. 169–180, 2013.
[35] Y. Liang and P. Zhao, “Similarity search in graph databases: A multilayered indexing approach,” in Data Engineering (ICDE), 2017 IEEE
33rd International Conference on. IEEE, 2017, pp. 783–794.
[36] H. W. Kuhn, “The hungarian method for the assignment problem,” Naval
research logistics quarterly, vol. 2, no. 1-2, pp. 83–97, 1955.
[37] K. Riesen, A. Fischer, and H. Bunke, “Improving approximate graph edit
distance using genetic algorithms,” in Joint IAPR International Workshops on Statistical Techniques in Pattern Recognition and Structural
and Syntactic Pattern Recognition. Springer, 2014, pp. 63–72.
[38] Wikipedia, “Digamma function,” 2017. [Online]. Available: https://en.
wikipedia.org/w/index.php?title=Digamma function&oldid=767319315

A PPENDIX A
P ROOF OF T HEOREM 1
Proof: Since GED satisfies the triangle inequality, we have
GED(G1 , G2 ) ≤ GED(G1 , G01 ) + GED(G01 , G02 ) + GED(G02 , G2 )
GED(G01 , G02 ) ≤ GED(G01 , G1 ) + GED(G1 , G2 ) + GED(G2 , G02 )
According to Definition 1, adding a vertex or an edge with the virtual label
ε is not counted as a graph edit operation. Therefore, we have:
GED(G1 , G01 ) = GED(G01 , G1 ) = 0
GED(G2 , G02 ) = GED(G02 , G2 ) = 0
Therefore, GED(G1 , G2 ) ≤ GED(G01 , G02 ) ≤ GED(G1 , G2 ).
That is, GED(G1 , G2 ) = GED(G01 , G02 ).

A PPENDIX B
P ROOF OF T HEOREM 2
Proof: Let the sets of branches rooted at virtual vertices in G01 , G02 be
∆BG1 and ∆BG2 , respectively. We have:
|BG0 ∩ BG0 | = |BG1 ∩ BG2 | + |∆BG1 ∩ BG2 |
1

2

+ |BG1 ∩ ∆BG2 | + |∆BG1 ∩ ∆BG2 |
Since branches rooted at virtual vertices are not isomorphic to any other
branches, we also have:
|∆BG1 ∩ BG2 | = |BG1 ∩ ∆BG2 | = |∆BG1 ∩ ∆BG2 | = 0
Therefore, |BG0 ∩ BG0 | = |BG1 ∩ BG2 |.
1
2
From the definitions of branches and extended graphs, we have:
max{|V1 |, |V2 |} = max{|V10 |, |V20 |}
From the definition of GBD, we can obtain:
GBD(G1 , G2 ) = max{|V1 |, |V2 |} − |BG1 ∩ BG2 |
= max{|V10 |, |V20 |} − |BG0 ∩ BG0 | = GBD(G01 , G02 )
1

2

That is, GBD(G1 , G2 ) = GBD(G01 , G02 ).

A PPENDIX C
C LOSED F ORMS OF E QUATIONS
A. Closed Form of Equation (8)
Λ1 (G01 , G02 ; τ, ϕ) = P r[GBD = ϕ | GED = τ ]
X
X
X
=
Ω1 (x, τ )
Ω2 (m, x, τ )
Ω3 (r, ϕ)Ω4 (x, r, m)
x

m


Ω1 (x, τ ) = H x; |V10 | +

|V10 |
2



Ω3 (r, ϕ) =

|V10 |
, |V10 |, τ
2



−1

 Pm
t=0
τ −x
r  (D−1)ϕ
· Dr
r−ϕ

Ω2 (m, x, τ ) =

(27)

r

(−1)m−t

Ω4 (x, r, m) = H x + m − r; |V10 |, m, x
 M −K  M −1
H(x; M, K, N ) = K
x
N −x
N
|V10 |+|LE |−1
D = |LV | ·
|L |
E

(28)
 
t


|V10 | m
2
m
t
τ −x

(29)
(30)



(31)
(32)
(33)

B. Closed Formqof Equation (16)

1
P r[GED] = C
nP
P
Z = Λ1
x Ω1
m
1

P2τ

ϕ=0

∂Ω2
∂τ

P

Then:

Λ1 · Z 2

r Ω3 Ω4 +

(34)
∂Ω1
x ∂τ

P

P

m Ω2

P

r Ω3 Ω4

o

v(v−1)
v 1
2
· F1
x
τ −x

1 v(v−1)−1
P
v
2
· F2 · m
t
m
τ −x

d
dτ

Ω2 (m, x, τ ) =



Likewise, by marginalizing out Z from Θ1 (x, τ ) and applying Chain Rule:
X
Θ1 (x, τ ) =
{P r[GBD = ϕ | Z = m, X = x, Y = τ − x]
· P r[Z = m | X = x, Y = τ − x]}

(36)

τ

F3 · F4

(37)

where

From Bayesian Network in Figure 3, we have:
X
Θ1 (x, τ ) =
{P r [GBD = ϕ | X = x, Z = m]

F1 = H(τ ) − H( 12 v(v + 1) − 2τ ) − H(τ − x)
+ H(x − τ +

1
v(v
2

(45)

m

1 v(v+1)
2

Ω1 (x, τ ) =

{Ω1 (x, τ ) · Θ1 (x, τ )}

x

(35)

where Λ1 is defined in Equation (27), and Ω1 , Ω2 , Ω3 and Ω4 are defined
in Equations (28)∼(32), respectively. In addition, we have:
d
dτ

X

Ω=

− 1))

F2 = ψ(τ − x + 1) − ψ(x + 1 − τ + 21 v(v − 1))
 1 t(t−1)
2
F3 = (−1)m−t m
t
τ −x

m

·P r [Z = m | Y = τ − x]}

(38)
(39)

Define:
Ω2 (m, x, τ ) = P r [Z = m | Y = τ − x]
Θ2 (m, x, ϕ) = P r [GBD = ϕ | X = x, Z = m]

(40)

F4 = 1 + ψ(x + 1 − τ + 21 t(t − 1)) − ψ(τ − x + 1)
(41)
Here, v is short for |V10 |, τ is short for of GED, and x, m, t are summation
subscripts in Equation (27). H(n) is the n-th Harmonic Number, and ψ(·)
is the Digamma Function [38].

(46)
(47)

We have:
Θ1 (x, τ, ϕ) =

X

{Ω2 (m, x, τ ) · Θ2 (m, x, ϕ)}

(48)

m

Likewise, by marginalizing out R from Θ2 (m, x, ϕ) and then applying Chain
Rule:
X
Θ2 (m, x, ϕ) =
{P r [GBD = ϕ | R = r, X = x, Z = m]

A PPENDIX D
P ROOF OF E QUATION (8)
Define:

r

Ω = P r[GBD = ϕ | GED = τ ]

·P r [R = r | X = x, Z = m]}

(42)

From Bayesian Network in Figure 3, we have:
X
Θ2 (m, x, ϕ) =
{P r [GBD = ϕ | R = r]

By marginalizing out S from Ω:
X
Ω=
{P r[GBD = ϕ, S = s | GED = τ ]}

r

s

By applying Chain Rule on Ω:
X
Ω=
{P r[GBD = ϕ | S = s, GED = τ ] · P r[S = s | GED = τ ]}

· P r [R = r | X = x, Z = m]}
Define:

s

In our model, every GED operation sequence is randomly selected with same
probability, therefore,
P r[S = s | GED = τ ] = 1/|SEQ|

(49)
(50)

Finally, from Equations (42) ∼ (50) we can obtain:
X
Ω=
Ω1 (x, τ ) · Θ1 (x, τ )
x

Let N = |SEQ| and we have:
1 X
{P r [GBD = ϕ | S = s, GED = τ ]}
Ω=
N s

=

X

Ω1 (x, τ ) ·

x

=

By marginalizing out X and Y from Ω:
1 X
{P r[GBD = ϕ, X = x, Y = y | S = s, GED = τ ]}
Ω=
N s,x,y
Since Y = τ − x when given GED = τ and X = x:
1 X
{P r[GBD = ϕ, X = x, Y = τ − x | S = s, GED = τ ]}
Ω=
N s,x
By applying Chain Rule on Ω:
1 X
Ω=
{P r[GBD = ϕ | X = x, Y = τ − x, S = s, GED = τ ]
N s,x
· P r[X = x, Y = τ − x | S = s, GED = τ ]}
According to the Bayesian Network in Figure 3, we have:
1 XX
Ω=
{P r [GBD = ϕ | X = x, Y = τ − x]
N s x
·P r [X = x, Y = τ − x | S = s]}
Define:
Θ1 (x, τ ) = P r [GBD = ϕ | X = x, Y = τ − x]
1 X
Ω1 (x, τ ) =
P r [X = x, Y = τ − x | S = s]
N s

Ω3 (r, ϕ) = P r [GBD = ϕ | R = r]
Ω4 (x, r, m) = P r [R = r | X = x, Z = m]

(43)
(44)

X

X

Ω2 (m, x, τ ) · Θ2 (m, x, ϕ)

m

Ω1 (x, τ ) ·

x

X

Ω2 (m, x, τ ) ·

m

X

Ω3 (r, ϕ) · Ω4 (x, r, m)

r

Therefore, Equation (8) is proved, while the formulae for calculating Ω1 ,
Ω2 , Ω3 and Ω4 are given in Lemmas 1, 2, 3 and 4, respectively. Please refer
to Appendices E, F, G and H for Lemmas 1, 2, 3 and 4, respectively.

A PPENDIX E
L EMMA 1 AND ITS P ROOF
Lemma 1: Given GED = τ , for any integer x ∈ [0, τ ], we have:
P
1
Ω1 (x, τ ) = |SEQ|
s∈SEQ P r[X = x, Y = τ − x|S = s]


0 
= H x; |V10 | + |V21 | , |V10 |, τ

(51)

where function H(x; M, K, N ) is defined in Equation (??).
Proof: From the definitions, Ω1 (x) is the probability of a random graph
edit sequence seqs exactly relabelling x vertices and τ − x edges. Since
0 
the extended graph G01 is a complete graph, it has |E10 | = |V21 | edges.
0 
Therefore, the number of ways to choose x vertices for relabelling is |Vx1 | ,
 0 
|V1 | 
2
and the number of ways to choose τ − x edges for relabelling is
.
τ −x
Then we have:
 0 
|V1 | 
|V10 |
2


·
0 
x
τ −x
Ω1 (x, τ ) =
= H x; |V10 | + |V21 | , |V10 |, τ
 0 
|V1 | 
0
|V |+
1

τ

2

where function H(x; M, K, N ) is defined in Equation (??).

A PPENDIX F
L EMMA 2 AND ITS P ROOF
Lemma 2: Given GED = τ , for any integer x ∈ [0, τ ] and m ∈
[0, |V10 |], we have:
Ω2 (m, x, τ ) = P r [Z = m | Y = τ − x]
|V10 |
2



=

 −1



Pm

t=0

τ −x

 
t 
|V10 | m
2
t
m
τ −x

(−1)m−t

(52)

Proof: Let x0 = τ − x. Since the extended graph G01 is a complete
0 
graph, it has |E10 | = |V21 | edges. From definitions of Y and Z, Ω2 (m, x)
can be modelled by following problem:
• Randomly select x0 edges from a complete graph with |V10 | vertices and
|V10 |
edges, what is the probability of these edges exactly covering m
2
vertices.
Let VX be the vertices covered by an edge subset X ⊆ E10 . By the
inclusion-exclusion principle, for any vertex subset S ⊆ V10 where |S| = m,
the number of possible edge sets X which satisfies VX = S is:
 

kx0 ,m =

t
2

m
t

Pm

m−t
t=0 (−1)



x0

0 
Since the number of sets S ⊆ V10 of size m is |Vm1 | , the total number of
ways to pick an edge set X which satisfies |X| = x0 and |VX | = m is:

Kx0 ,m =

Pm

t=0

|V10 | m
m
t

(−1)m−t

 
t
2



x0

Also the number of ways to pick an edge set X ⊆ E10 is:
|V10 |
2
x0



Kx0 =

P

m

Kx0 ,m =

Kx0 ,m
Kx0

=

where the number of possible branch types D is the number of ways to
assign |LV | + 1 labels (including the virtual label) to the vertex in a branch,
multiplying the number of ways to assign |LE | + 1 labels to |V | − 1 edges
in the same branch. This is a variant of the classic object colouring problem
and we omit the detailed proof here.

A PPENDIX H
L EMMA 4 AND ITS P ROOF





Lemma 4: Given GBD = ϕ, for any integer r ∈ [0, |V10 |], we have:

Therefore, we have:
Ω2 (m, x, τ ) =

example. To model the probability of this situation, we define t = r − ϕ, so t
is essentially the size of B̃G1 . Here we assume that the occurrence probability
of each branch type is equal, and denote the number of branch types by D.
Therefore, Ω3 (r) can be abstracted as the ball-pair colouring problem as
follows.
• Given two lists of balls A1 , A2 of size r, where each ball has been
randomly coloured by one of the D colours. We define the ball pair as
two balls where one ball from A1 and another from A2 , where every ball
can only be in one pair. What is the probability that there are exactly t
ball pairs where the colour inside each pair is the same.
To solve this problem, we take the following steps to count the possible
ways to colour balls with the conditions above.
1) The number of ways to form ball pairs is r!;
2) When pairs are fixed, the number of ways to select t pairs and assign t
colours to them is rt Dt ;
3) For the remaining r −
t pairs, the number of ways to assign them different
colours inside each pair is [D(D − 1)]r−t .
Since there are totally r! · D2r ways to form ball pairs and assign colours
to all balls, we have:

 r  (D − 1)ϕ
r! rt Dt [D(D − 1)]r−t
·
=
Ω3 (r, ϕ) =
2r
r−ϕ
r! · D
Dr

Ω4 (x, r, m) = P r[R = r | X = x, Z = m]

= H x + m − r; |V10 |, m, x

Kτ −x,m
Kτ −x

So Lemma 2 is proved.

A PPENDIX G
L EMMA 3 AND ITS P ROOF
Lemma 3: Given GBD = ϕ, for any integer r ∈ [0, |V10 |], we have:
r  (D−1)ϕ
(53)
Ω3 (r, ϕ) = P r [GBD = ϕ | R = r] = r−ϕ
· Dr
where D is the number of all possible branch types, and
|V10 |+|LE |−1
D = |LV | ·
|L |

(54)

E

where H(x; M, K, N ) is defined in Equation (??).
Proof: Let t = x + m − r, so t is the number of vertices both relabelled
and covered by relabelled edges. Since the order of relabelling operations
does not affect the graph edit result, Ω4 (x, r, m) can be modelled by the
following problem:
• First randomly select m vertices from V10 and tag these vertices as special
ones. Then randomly select x vertices from V10 . What is the probability of
exactly selecting t special vertices in the second selection.
Since the number of ways to exactly select t special vertices in the second
 |V 0 |−m
0 
1
pick is |Vm1 | m
, and the number of all ways to pick m and x
t
x−t
0 
0 
vertices separately from V10 is |Vm1 | |Vx1 | , we have:

Proof:
Ω4 (x, r, m) =

v1

A
x

y
v3

B

C
Fig. 30.

u1

y
B

G1

v2

B

B
x

y
v4

u3

u2

y
A u4

C
G2

Graphs for the Proof of Lemma 3

From the definition, R = r means there are exactly r branches whose
vertices or edges are relabelled when transforming G1 into G2 during the
graph edit process. For simplicity, we define these branches as relabelled
branches. However, the value of R could be larger than the difference between
branches in G1 and G2 , i.e., GBD(G1 , G2 ), because it is possible that a
subset of relabelled branches are just a re-ordering of the original ones, which
are denoted by B̃G1 .
For instance, in the graph edit process of transforming G1 into G2 , as
shown in Figure 30, we need to relabel v1 by label B and v4 by label
A, which means that the number of relabelled branches R = 2. However,
GBD(G1 , G2 ) = 0 since this graph edit process just swaps the branches
rooted at v1 and the one rooted at v4 , so B̃G1 = {B(v1 ), B(v4 )} in this

(55)

=

|V10 | m |V10 |−m
m |V10 |−m
m
t
t
x−t
x−t
=
|V10 | |V10 |
|V10 |
m
x
x

H x + m − r; |V10 |, m, x

where function H(x; M, K, N ) is defined in Equation (??).

A PPENDIX I
G RAPH G ENERATING A LGORITHM
The algorithm for generating the synthetic graphs (i.e., Syn-1 and Syn-2
data sets) is as follows.
The algorithm aims to generate a set of graphs G, such that for ∀gi , gj ∈
G, the edit distance between gi and gj is known. In order to achieve this goal,
we first defined a valid modification center, then we designed the generation
algorithm which consists two phases: (1) Generate a random “qualified” graph
as a template; (2) Modify the template to generate the graph set G.
A modification center is a vertex vc in graph g such that ∀vi , vj ∈
{neighbors of vc }, i 6= j, the edit distance between g −e(i, c) and g −e(j, c)
is greater than 0, where e(i, c) is the edge between vertices vi and vc , and
e(j, c) is the edge between vertices vj and vc . For any modification center
vc , if we randomly modify its adjacent edges, the edit distance between the
original and modified graphs can be simply calculated by comparing the edges
adjacent to their modification centers in polynomial time.
However, identifying whether a vertex is a modification center is difficult.
Therefore, we propose a relatively efficient signature that can help us to filter

out some special cases that a vertex is certainly a modification center. For
vertex vc ’s neighbor vi , the signature is defined as {s0 , s1 , s2 , . . . , sn },
where s0 is a set contains vi ’s label, and other sets sk = {(vj .label
, e(i, j).label); ∀vj ∈ k-hop neighbors of vi }, k > 0. If a vertex is a modification center, then its neighbors’ signatures must be pair-wised different. That
is, if we find a vertex whose neighbors’ signatures are pair-wised different,
this vertex must be a modification center. If there is no such a vertex, we
re-generate the graph until success.
In our settings, the graph should be connected, and have at least one
modification center of degree at least d, to produce a set of graphs among
which the edit distance varies from 0 to d. To ensure the connectivity of the
graph, we force each vertex vi is connected to another vertex vj , where i > j.
After all vertices being connected, we then add remaining edges according to
the type of the graph. For random graphs, we randomly add edges between
in-adjacent vertices. For scale-free graphs, we add constant number of edges
to each vertex vi , where the neighbors of vi is picked from {vj ; ∀j < i}
with the probability proportional to their degrees.

A PPENDIX J
S UPPLEMENTAL F IGURES F OR E XPERIMENTS ON
S YNTHETIC DATA S ETS
Please refer to Figures 31∼42 on the last page.

A PPENDIX K
T HEOREM FOR AVERAGE D EGREE OF S CALE - FREE
G RAPHS
Theorem 5: The average degree of a scale free graph G is O(log n),
where n is the number of vertices in graph G.
Proof: The fraction of vertices with degree k in scale-free graphs is
C · k−δ [16], where 2 < δ < 3, and C is a constant. Therefore, the average
degree in scale-free graphs is:
n−1
n−1
X
X C
= C · H(n − 1) = O(log n)
(56)
d=
k · Ck−δ <
k
k=1
k=1
where H(n) is the n-th Harmonic Number.

LSAP
greedysort
seriation

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

1.0

0.9

0.9

0.8

0.8

0.8

0.7

0.7

0.7

0.6

0.6

0.6

0.5
0.4
0.3

precision

1.0

0.9

0.5
0.4
0.3

0.3
0.2

0.1

0.1

0.1

50K

100K

1K 10K 20K

LSAP
greedysort
seriation

50K

100K

1K 10K 20K

graph size

Fig. 31. Precision vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
15)

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

Precision vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
25)

0.9

0.7

0.8

0.8

0.6

0.7

0.7

0.5

0.6

0.6

recall

1.0

0.9

recall

1.0

0.8

0.5

0.5

0.3

0.4

0.2

0.3

0.3

0.1

0.2

0.2

0.1

0.1

1K 10K 20K

50K

0.4

100K

graph size

1K 10K 20K

Precision vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
30)
Fig. 34.

50K

100K

1K 10K 20K

graph size

Fig. 36. Recall vs. graph size on Syn-1

Data Set with various γ (τ̂ = 20)

LSAP
greedysort
seriation

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

LSAP
greedysort
seriation

0.8

0.8

0.8

0.7

0.7

0.7

0.6

0.6

0.6

0.5

f1

0.9

recall

0.9

0.9

0.5

0.4

0.4

0.3

0.3

0.3

0.2

0.2

0.2

0.1

0.1

0.1
1K 10K 20K
1K 10K 20K

50K

100K

1K 10K 20K

graph size

50K

Fig. 37. Recall vs. graph size on Syn-1

Fig. 38. Recall vs. graph size on Syn-1

Data Set with various γ (τ̂ = 25)

Data Set with various γ (τ̂ = 30)

LSAP
greedysort
seriation

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

F1-score vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
15)
LSAP
greedysort
seriation

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)
1.0

0.9

0.9

0.8

0.8

0.8

0.7

0.7

0.7

0.6

0.6

0.6

0.5

0.5

0.5

f1

1.0

0.9

0.4

0.4

0.4

0.3

0.3

0.3

0.2

0.2

0.2

0.1

0.1

0.1

1K 10K 20K

50K

100K

1K 10K 20K

graph size

F1-score vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
20)
Fig. 40.

50K

100K

1K 10K 20K

graph size

F1-score vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
25)

Fig. 41.

100K

Fig. 39.

1.0

f1

f1

LSAP
greedysort
seriation

50K

graph size

100K

graph size

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

1.0

1.0

0.4

100K

Data Set with various γ (τ̂ = 15)

1.0

0.5

50K

graph size

Fig. 35. Recall vs. graph size on Syn-1

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

LSAP
greedysort
seriation

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

LSAP
greedysort
seriation

0.9

0.4

100K

Fig. 33.

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

LSAP
greedysort
seriation

50K

graph size

Fig. 32. Precision vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
20)

1.0

precision

0.4

0.2

graph size

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

0.5

0.2

1K 10K 20K

recall

LSAP
greedysort
seriation

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

1.0

precision

precision

LSAP
greedysort
seriation

GBDA(γ = 0.60)
GBDA(γ = 0.70)
GBDA(γ = 0.80)

50K

100K

graph size

F1-score vs. graph size on
Syn-1 Data Set with various γ (τ̂ =
30)
Fig. 42.

