<figure><img src="image_1.png" /></figure>

HAL
open science

# Graph sampling with applications to estimating the number of pattern embeddings and the parameters of a statistical relational model

Irma Ravkic, Martiň Znidaršič, Jan Ramon, Jesse Davis

► To cite this version:

Irma Ravkic, Martiň Znidaršič, Jan Ramon, Jesse Davis. Graph sampling with applications to estimating the number of pattern embeddings and the parameters of a statistical relational model. Data Mining and Knowledge Discovery, 2018, pp.36. 10.1007/s10618-018-0553-2. hal-01725971

HAL Id: hal-01725971

https://inria.hal.science/hal-01725971

Submitted on 7 Mar 2018

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# Graph sampling with applications to estimating the number of pattern embeddings and the parameters of a statistical relational model

Irma Ravkic · Martin Žnidaršič · Jan Ramon · Jesse Davis

Received: date / Accepted: date

**Abstract** Counting the number of times a pattern occurs in a database is a fundamental data mining problem. It is a subroutine in a diverse set of tasks ranging from pattern mining to supervised learning and probabilistic model learning. While a pattern and a database can take many forms, this paper focuses on the case where both the pattern and the database are graphs (networks). Unfortunately, in general, the problem of counting graph occurrences is \#P-complete. In contrast to earlier work, which focused on exact counting for simple (i.e., very short) patterns, we present a sampling approach for estimating the statistics of larger graph pattern occurrences. We perform an empirical evaluation on synthetic and real-world data that validates the proposed algorithm, illustrates its practical behavior and provides insight into the trade-off between its accuracy of estimation and computational efficiency.

**Keywords** graph sampling · graph pattern matching · parameter estimation · statistical relational learning

The first two authors contributed equally to this work and are ordered alphabetically.

**Irma Ravkic**

Department of Computer Science, KU Leuven, 3001 Heverlee, Leuven, Belgium

E-mail: irma.ravkic@cs.kuleuven.be

**Martin Žnidaršič**

Jožef Stefan Institute, Jamova cesta 39, SI-1000 Ljubljana, Slovenia

E-mail: martin.znidarsic@ijs.si

**Jan Ramon**

Department of Computer Science, KU Leuven, 3001 Heverlee, Leuven, Belgium

E-mail: jan.ramon@cs.kuleuven.be

**Jesse Davis**

Department of Computer Science, KU Leuven, 3001 Heverlee, Leuven, Belgium

E-mail: jesse.davis@cs.kuleuven.be
# 1 Introduction

Counting the number of times a pattern occurs in a database is an essential step in many applications. In frequent pattern mining, computing the support of an itemset, a sequence or a graph requires counting the number of times it occurs in the data. In supervised learning, many algorithms measure the correlations between features and target variables, which requires iterating over the examples that satisfy a given pattern. In probabilistic modeling, estimating the parameters of a model (e.g., constructing a conditional probability table for a Bayesian network) requires counting the number of times a set of variables takes on a specific combination of values in the data. Furthermore, mining and learning applications typically require addressing this counting problem for a large number of patterns. For example, pattern mining approaches compute the support for progressively more complex patterns and probabilistic model structure learning algorithms require estimating the parameters for a large number of candidate structures.

Unfortunately, in general, counting the number of occurrences of a graph pattern in a network is \#P-complete. Most earlier work addressed problems such as graph characterization (Shervashidze et al 2009; Bordino et al 2008) or frequent pattern discovery (Inokuchi et al 2003; Wernicke 2005), which consider all possible subgraphs. However, these approaches tend to focus only on small patterns.

In this paper, we propose a new algorithm for approximately counting the number of embeddings of a pattern in a graph based on a theoretical sampling idea suggested by Fürer and Kasiviswanathan (2014) for larger unlabeled undirected graphs. For practical use in a data mining context, the theoretical algorithm of Fürer and Kasiviswanathan (2014) has several weaknesses and open questions. First, it assumes that all graphs are unlabeled and undirected. Second, it assumes that the input pattern has an ordered-bipartite decomposition (OBD) whereas not all interesting patterns will have an OBD. Third, it assumes that the pattern's OBD is provided to the algorithm. From a practical perspective, it is necessary to automatically construct an OBD (if one exists) for an input pattern. Finally, only the count of the number of embeddings is estimated.

This paper attempts to address the aforementioned concerns and makes the following contributions. First, we generalize the approach of Fürer and Kasiviswanathan (2014) to labeled and directed graphs. Second, we propose a heuristic approach to automatically search for an OBD for a given pattern, whereas it was previously assumed that the OBD was given. Third, while Fürer and Kasiviswanathan (2014) provide only an error bound, in this paper we also estimate the accuracy of computed statistics from the sample. Fourth, we implement and empirically evaluate the approach on two tasks whereas the work of Fürer and Kasiviswanathan (2014) was purely theoretical. The first task investigates the accuracy of the estimated counts. The second task evaluates the suitability of using the approach for parameter learning within statistical relational learning (SRL). Fifth, we investigate how the algorithm performs when some of the conditions required for the theoretical guarantees are violated. Namely, we investigate the algorithm's performance when it is given (1) a decomposition that is not a OBD, and (2) a data graph that is not an Erdős-Rényi random graph.<fnref n="1" /> Finally, our implementation is publicly available.<fnref n="2" />

<footnote n="1">¹ Graphs where each edge is included, independent of all other edges, with probability p.</footnote>

<footnote n="2">² For the code see: https://dtai.cs.kuleuven.be/software/gs-srl</footnote>
## 2 Related work

Data in domains, such as the Semantic Web, social networks, citation networks, biology, relational learning and geography among others, are often naturally represented as a graph. Therefore, analyzing large graphs or networks is a very active research field. Furthermore, problems related to the analysis (searching, counting, etc.) of subgraphs in large network graphs attract significant attention. Next, we present some of the more prominent lines of research concerned with subgraph analysis and their relation to this paper.

Many of the following approaches employ subgraph identification, either by *subgraph isomorphism* or *subgraph homomorphism* as a subtask. Typically, they employ existing state-of-the-art algorithms for this task (e.g., Cordella et al (2004); Ullmann (1976)).

One line of work looks at defining features of graphs. Graph particles are a promising approach for *graph characterization and comparison* (Pržulj 2007; Shervashidze et al 2009; Bordino et al 2008). Graph particles are usually small subgraphs that contain up to five nodes, and are sometimes called graphlets. The distribution of graph particles can be used as a “fingerprint”, a characteristic invariant, of a graph. Then, the similarity between pairs of graphs can be measured by comparing their distributions of graph particles. *Graph querying* (Di Natale et al 2010; Giugno and Shasha 2002) focuses on searching for a graph in a database of graphs. Here, most research focuses on identifying graph features (e.g., paths, walks, subtrees, small subgraphs, etc.) for filtering and indexing in order to enable faster searching and querying of graphs.

*Frequent pattern mining* in graph databases is another active topic of research (Yan and Han 2002; Inokuchi et al 2003). This usually entails finding frequent subgraphs, which are also called *motifs* in biological networks (Kashtan et al 2004; Wernicke 2005). Motifs are subgraphs that appear more often than expected according to a null model. Although subgraph isomorphism is an integral subtask in this problem, the main research focus is determining the significance of subgraphs, which is measured by their support in the data. Usually, these approaches only consider small subgraphs since the goal is to find all subgraphs that meet the support threshold. Some approaches (e.g., Baskerville et al (2007)) use approximate methods to scale to larger subgraphs.

Finally, general graph sampling approaches explore how to obtain a representative sample of a graph for a specific purpose (Leskovec and Faloutsos 2006). Some approaches, like the work of Zou and Holder (2010), rely on sampling for performing subgraph analysis. Zou and Holder (2010) first sample from a single large graph and then mine the sample for frequent subgraphs. They empirically assess the viability of various graph sampling approaches for frequent subgraph mining.

The two key differentiating factors of our work are that we focus on approximating the frequency of a subgraph in larger graph (i.e., the inverse task of Zou and Holder (2010)), and we consider relatively large subgraphs (up to 15 nodes). Perhaps the most closely related work is that of Bordino et al (2008) who approximately count the subgraphs for graph characterization. But again, they only consider small subgraphs.

Within SRL, there has been work on using approximations to scale up both parameter estimation and structure learning. The early work of Kok and Domingos
(2005) recognized the potential of sampling to improve the scalability of structure learning for Markov logic networks. Their approach offered the ability to both subsample the data, as well as the number of true groundings of a clause. More recently, Venugopal et al (2015) proposed an approach that improves the scalability of parameter learning in Markov logic networks by approximately counting the number of satisfied groundings of a clause. This problem is equivalent to counting paths in a graph. Venugopal et al (2015)'s approach is specific to Markov logic and is focused on the parameter learning task. The fast approximate counting (FACT) approach proposed by Das et al (2016) uses a graph database to scale learning and lifted inference. In contrast to this paper's sampling-based scheme, their approach is based on message passing and looks at the in and out degrees of a node relative to the maximum possible degree. Das et al (2016)'s approach does not offer any theoretical guarantees whereas this paper inherits the original guarantees from the work of Fürer and Kasiviswanathan (2014) (under certain conditions). Furthermore, Das et al (2016) focus solely on estimating a pattern's count, whereas we discuss computing additional statistics (see Section 8).

# 3 Background

Throughout this paper, we use the term *domain graph* to refer to the graph that represents the full dataset. We use the term *pattern graph* to refer to a small subgraph that we are evaluating (e.g., to estimate its frequency in the domain graph). Next, we introduce the necessary background information and terminology used throughout this paper.

## 3.1 Graph representation

A labeled graph is a tuple $G = (V, E, \Sigma, \lambda)$, where $V$ is a set of vertices (also called nodes), $E \subseteq \{\{u, v\} | u, v \in V\}$ is a set of edges, $\Sigma$ is a set of labels and $\lambda: (V \cup E) \rightarrow \Sigma$ is a function that assigns labels to nodes and edges. For a graph $G$, we refer to its set of vertices with $V(G)$, to its set of edges with $E(G)$, to its alphabet with $\Sigma_G$ and to its labeling function with $\lambda_G$. For a graph $G$ and a vertex $v \in V(G)$ the set of *neighbors* of $v$ in $G$ is $N_G(v) = \{u \in V(G) | \{u, v\} \in E(G)\}$. A graph $H$ is a *subgraph* of a graph $G$ if $V(H) \subseteq V(G)$, $E(H) \subseteq E(G)$ and $\forall x \in V(H) \cup E(H) : \lambda_H(x) = \lambda_G(x)$. The label set $\Sigma$ can contain attribute-value pairs as well. More precisely, for a node $v$ with label *attribute = value* we denote the *value* of $v$ with *value($v$)*. We also allow for wildcards `#` to be used as *value*. This means that the node can take on values, but we do not care in this particular case which value it is. We will often omit `#` in *attribute = #*, and denote the label of the node with *attribute*.

Consider the labeled graph in Figure 1. The nodes in the graph denote people (men or women) together with their friendship or marriage relations to other people, and their satisfaction with their salary. Note that labels for people, salary and satisfaction are attribute-value pairs (e.g., *satisfaction = high* for *mike*).
<figure><img src="image_2.png" /><figcaption>Fig. 1: An example of a labeled domain graph depicting a number of *people* (*men* or *women*) together with their *friendship* or *marriage* relations to other people and their *satisfaction* with their *salary*.</figcaption></figure>

## 3.2 Pattern matching

There are multiple ways of performing the matching between two graphs. A homomorphism from a graph $H$ to a graph $D$ is a mapping $\varphi : V(H) \rightarrow V(D)$ such that (i) $\forall v \in V(H), \lambda_H(v) = \lambda_D(\varphi(v))$, (ii) $\forall \{u, v\} \in E(H), \{\varphi(u), \varphi(v)\} \in E(D) \wedge \lambda_H(\{u, v\}) = \lambda_D(\{\varphi(u), \varphi(v)\})$. With $emb_{hom}(H, D)$ we denote the set of all homomorphisms between $H$ and $D$. A subgraph isomorphism from $H$ to $D$ is a homomorphism $\varphi$ from $H$ to $D$ such that $\forall u, v \in V(H), \varphi(u) \neq \varphi(v)$. We denote with $emb_{iso}(H, D)$ the set of all subgraph isomorphisms between $H$ and $D$. Whether homomorphism or subgraph isomorphism is more appropriate to use as a matching operator depends on the application. Our results are generic and apply to both cases. We use $emb(H, D)$ to refer to $emb_x(H, D)$ for any $x \in \{\text{hom, iso}\}$ and refer to its elements as embeddings.

A pattern over $\Sigma$ is a graph $P$ with $\Sigma_P = 2^\Sigma$ (i.e., the labels on vertices and edges of $P$ are subsets of $\Sigma$). We use the term values to refer to the labels of the images of pattern vertices. In particular, for an embedding $\varphi$ of $P$ in $D$ and a set $S \subseteq V(P)$, we define $val_{S:D}(\varphi) = \lambda_D \circ \varphi|_S = \{(u, \lambda_D(v)) \mid u \in S \wedge (u, v) \in \varphi\}$. We define the range of values for $S$ of $P$ in $D$ by $Range(S: P, D) = \{val_{S:D}(\varphi) \mid \varphi \in \text{emb}(P, D)\}$. We also define the multiset of values for $S$ of $P$ in $D$ by $Val(S: P, D) = \{val_{S:D}(\varphi)\}_{\varphi \in \text{emb}(P, D)}$. So the essential difference between $Range(S: P, D)$ and $Val(S: P, D)$ is that the first is a normal set (containing each value once) while the second is a multiset possibly containing the same value several times. For the special case $S = V(P)$, we define $val_D(\varphi) = val_{V(P):D}(\varphi)$, $Range(P, D) = Range(V(P): P, D)$, and $Val(P, D) = Val(V(P): P, D)$.
# 4 Decomposition of graphs

Decomposing graphs permits applying a divide-and-conquer strategy to complex problems and it represents a preprocessing step for many graph theory techniques. In this paper we use this strategy for obtaining the count of a pattern's embeddings in a graph. The general idea is to decompose the vertices involved in a complex pattern into multiple partitions. The matching task then becomes simpler as only the vertices in a single partition need to be matched simultaneously, with the current partition's matching being conditioned on the assignment given to the vertices in previous partitions. We describe two decompositions: an ordered bipartite decomposition (Fürer and Kasiviswanathan 2014) and an arbitrary decomposition.

## 4.1 Ordered bipartite decompositions

Fürer and Kasiviswanathan (2014) introduced the concept of an ordered bipartite decomposition (OBD), which is a crucial component of our approach. We now define the relevant terminology and describe this decomposition.

An *independent set* is a set of vertices in a graph, no two of which are adjacent. A *graph partition* is a division of a graph's vertices and edges into smaller components with specific properties. A *bipartite graph* $G = (U, V, E, \Sigma, \lambda)$ is a graph whose vertices can be divided into two disjoint independent sets $U$ and $V$ such that every edge in $e \in E$ connects a vertex in $U$ to one in $V$. A bipartite graph can be seen as coloring the nodes in $U$ one color, and the nodes in $V$ another color, where each edge connects two vertices with different colors. For example, representing the relationship of football players to their clubs can be modeled as a bipartite graph: at a specific point in a time, each player is associated with one club.

Informally, an OBD decomposes a graph pattern into an implicit sequence of bipartite graphs by numerically labeling all the vertices in a graph pattern. The labeling must ensure that (1) each edge connects two vertices with different labels, and (2) all of a vertex's neighbors with a higher label have the same label. Now we formally define an OBD.

**Definition 1 (Ordered bipartite decomposition (Fürer and Kasiviswanathan 2014))** An ordered bipartite decomposition of a graph $P = (V_P, E_P)$ is a sequence $V_1, ..., V_l$ of subsets of $V_P$ such that:

1. $V_1, ..., V_l$ form a partition of $V_P$

2. Each of the $V_i$ (for $i \in [l] = 1, ..., l$) is an independent set in $P$

3. $\forall v \exists j$ such that $v \in V_i$ implies $N_P(v) \subseteq (\bigcup_{k<i} V_k) \cup V_j$

Fürer and Kasiviswanathan (2014) defined an OBD's *size* as the number of partition classes in it and its *width* as the number of vertices in the largest partition class.

Fürer and Kasiviswanathan (2014) proved that a number of graph classes (e.g., trees, grids, etc.) will always have an OBD. However, not all graphs will have an OBD. Consider the two graphs in Figure 2. Graph (a) is a triangle graph. There exists a trivial decomposition $\{v1\}, \{v2\}, \{v3\}$ satisfying Properties 1 and 2 of Definition 1. However, it is not possible to satisfy Property 3 for this graph,
hence, there is no OBD for it. The OBD for graph (b) is $\{v1\}, \{v2\}, \{v3, v4\}$ where each partition is labeled with, for example, I, II, and III. The reader can check that all three properties hold. For example, for the vertex $v2$ labeled as II, all the neighbors ($v3$ and $v4$) are labeled identically with a higher label III.

<figure><img src="image_3.png" /><figcaption>Fig. 2: An example of a) a non-bipartite triangle graph for which there does not exist an OBD, b) a graph that has the following OBD: $\{v1\}\{v2\}\{v3, v4\}$ (where $\{v1\}$ is the first partition, $\{v2\}$ is second partition, and $\{v3, v4\}$ is the third partition).</figcaption></figure>

## 4.2 Arbitrary decompositions

As many potentially interesting patterns may not have an OBD, we also want to investigate how the Fürer and Kasiviswanathan (2014) approach works when the algorithm is provided with a decomposition that is not an OBD. Here, we will still need a decomposition of the graph pattern and in this case, we will work with an arbitrary decomposition.

**Definition 2 (Arbitrary decomposition)** An arbitrary decomposition of a graph $P = (V_P, E_P)$ is a sequence $V_1, ..., V_l$ of subsets of $V_P$ such that:

1. $V_1, ..., V_l$ form a partition of $V_P$

2. Each of the $V_i$ (for $i \in [l] = 1, ..., l$) is an independent set in $P$

# 5 Counting the number of embeddings in graphs

The problem of collecting statistics about the occurrences of a graph pattern in a network can be formalized as:

- **Given:** a domain graph $D$, a pattern $P$, and a statistic of interest $f$

- **Find:** $f(P, D)$

The basic procedure that is used in all the approaches in this paper is to find all extensions to a partial embedding $\varphi$. For that purpose we use the procedure in Algorithm 1. The function receives a domain graph $D$, a pattern graph $P$, a partial embedding $\varphi$ and a set of vertices $S \subseteq P$ for which we want to find the images. The function loops through all the vertices $v \in S$ and tries to find images of $v$ which are then stored in the $C_v$ variable. It then constructs the set
of all possible embeddings such that each embedding is the concatenation of the
partial embedding $\varphi$ and one of the combination of images found for vertices in
$S$. The test $\varphi' \in \text{emb}(P[\varphi'^{-1}(D)], D)$ checks if the edges and labels are preserved
by $\varphi'$. If subgraph isomorphism is used as the matching operator (rather than
homomorphism), one must additionally check that no two vertices receive the
same image.

<table><caption>**Algorithm 1** Extend partial embeddings.</caption>
  <tr>
    <td>1:</td>
    <td>function <b>PARTIALEXT</b>(<i>D</i>, <i>P</i>, <i>&phi;</i>, <i>S</i>)</td>
    <td rowspan="7">▷ Loop through vertices in the partition</td>
  </tr>
  <tr>
    <td>2:</td>
    <td>for all <i>v</i> &isin; <i>S</i> do</td>
  </tr>
  <tr>
    <td>3:</td>
    <td>if &exist;(<i>x</i>, <i>y</i>) &isin; <i>&phi;</i>: <i>x</i> &isin; <i>N</i><sub><i>P</i></sub>(<i>v</i>) then</td>
  </tr>
  <tr>
    <td>4:</td>
    <td><i>C</i><sub><i>v</i></sub> &leftarrow; {<i>u</i> &isin; <i>N</i><sub><i>D</i></sub>(<i>y</i>) | &lambda;(<i>u</i>) &isin; &lambda;(<i>v</i>)}</td>
  </tr>
  <tr>
    <td>5:</td>
    <td>else</td>
  </tr>
  <tr>
    <td>6:</td>
    <td><i>C</i><sub><i>v</i></sub> &leftarrow; {<i>u</i> &isin; <i>V</i>(<i>D</i>) | &lambda;(<i>u</i>) &isin; &lambda;(<i>v</i>)}</td>
  </tr>
  <tr>
    <td>7:</td>
    <td>return {<i>&phi;</i>' | <i>&phi;</i>' = <i>&phi;</i> &cup; {(<i>v</i>, <i>u</i><sub><i>v</i></sub>)} <sub><i>v</i> &isin; <i>S</i> &and; &forall;<i>v</i> &isin; <i>S</i>: <i>u</i><sub><i>v</i></sub> &isin; <i>C</i><sub><i>v</i></sub> &and; <i>&phi;</i>' &isin; <i>emb</i>(<i>P</i>[<i>&phi;</i>'<sup>-1</sup>(<i>D</i>)], <i>D</i>)}</td>
  </tr>
</table>

Next, we introduce an exhaustive approach and a sampling approach for finding
the number of embeddings of a pattern in a domain graph.

5.1 Exhaustive approach

Algorithm 2 shows an exhaustive algorithm for computing all embeddings of a
pattern. The algorithm receives a domain graph *D*, a pattern graph *P*, and a fixed
ordering *O* on the pattern vertices as input. We use *O*<sub>*i*</sub> to denote the *i*-th vertex
in the ordering *O*, that is, the *i*-th vertex that is assigned an image.

Starting from an empty embedding, the algorithm works as follows. The PARTIALEXT function finds all extensions to a partial embedding $\varphi$ by assigning an image to the $O_{|\varphi|+1}$th vertex. Then, in line 7 these images are appended to the embedding. The test $\varphi' \in \text{emb}(P[\varphi'^{-1}(D)], D)$ checks if the edges and labels are preserved by $\varphi'$. If subgraph isomorphism is used as the matching operator (rather than homomorphism), one must additionally check that no two vertices receive the same image.

When choosing an ordering $\mathcal{O}$, a typical heuristic is to first select images for
vertices that have a small number of possible images. We require that every vertex,
except the first one, is adjacent to at least one other vertex which precedes it in
the ordering. A large literature exists on heuristics and various forms of lookahead
(Ullmann 1976), which is outside the scope of our paper.

<table><caption>**Algorithm 2** Exhaustive approach.</caption>
  <tbody>
    <tr>
      <td>1:</td>
      <td>function ALLEMB(graph D, pattern P, ordering O)</td>
      <td rowspan="2">▷ Return a set of embeddings</td>
    </tr>
    <tr>
      <td>2:</td>
      <td>return ALLEMB(D, P, O, {})</td>
    </tr>
    <tr>
      <td>3:</td>
      <td>function ALLEMB(D, P, O, &phi;)</td>
      <td rowspan="5">▷ If the pattern is fully embedded</td>
    </tr>
    <tr>
      <td>4:</td>
      <td>if |&phi;| = |V(P)| then</td>
    </tr>
    <tr>
      <td>5:</td>
      <td>return {&phi;}</td>
    </tr>
    <tr>
      <td>6:</td>
      <td>else</td>
    </tr>
    <tr>
      <td>7:</td>
      <td>return &cup;<sub>&phi;' &isin; PARTIALEXT(D, P, &phi;, {O<sub>|&phi;|+1</sub>})</sub> ALLEMB(D, P, O, &phi;')</td>
    </tr>
  </tbody>
</table>
## 5.2 Random vertex sampling

Algorithm 3 presents an approximate algorithm for computing all embeddings of a pattern based on sampling. It receives the same input as the exhaustive approach (i.e., function ALLEMB in Algorithm 2). At a high level, this algorithm samples the first part of the embedding randomly and then exhaustively computes all of its completions. This algorithm differs from the exhaustive algorithm in that it first finds an image for the first vertex $O_1$ by randomly sampling, with replacement, from the set of vertices with a matching label instead of exhaustively iterating over all of them.

<caption>**Algorithm 3 Random vertex sampling approach.**</caption>

1: **function** RNDV1EMB(graph *D*, pattern *P*, ordering *O*)
2: &nbsp;&nbsp;&nbsp;&nbsp;*M* ← {}
3: &nbsp;&nbsp;&nbsp;&nbsp;*R* ← {*u* ∈ *V*(*D*) | λ(*u*) = λ(*O*₁)} ▷ Find images for the first vertex of *P*
4: &nbsp;&nbsp;&nbsp;&nbsp;**while** not timeout **do**
5: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Select a random *u* ∈ *R*
6: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;*M* ← RNDV1EMB(*D*, *P*, *O*, {*(*O*₁, *u*)})
7: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;*M* ← *M* ∪ {*M*}
8: &nbsp;&nbsp;&nbsp;&nbsp;**return** *M*
9: **function** RNDV1EMB(*D*, *P*, *O*, *φ*) ▷ Find extensions for the partial embedding
10: &nbsp;&nbsp;&nbsp;&nbsp;**if** timeout **then**
11: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**return** {}
12: &nbsp;&nbsp;&nbsp;&nbsp;**else if** |*φ*| = |*V*(*P*)| **then**
13: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**return** {*φ*}
14: &nbsp;&nbsp;&nbsp;&nbsp;**else**
15: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;*M* ← {}
16: &nbsp;&nbsp;&nbsp;&nbsp;**for all** *φ'* ∈ PARTIALEXT(*D*, *P*, *φ*, {*|*φ*|+1}) **do**
17: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;*M'* ← RNDV1EMB(*D*, *P*, *O*, *φ'*)
18: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;*M* = *M* ∪ *M'*
19: &nbsp;&nbsp;&nbsp;&nbsp;**return** *M*

One can observe that every embedding has the same probability of being en-
countered by this algorithm. However, several embeddings may be included in the
sample at the same iteration, and hence these may not be independent. Therefore,
in line 7 we form a sample by a set of sets of embeddings (indicating these depen-
dencies explicitly) rather than just merging all sets *M* into one large unstructured
sample *M*. This algorithm is simple but has the potential disadvantage that, for
large patterns, it may spend all its time finding embeddings for the first sampled
vertex.

Similar strategies have previously been proposed, such as in the work on tri-
angle counting (c.g., Jowhari and Ghodsi (2005)).

# 6 The Fürer-Kasiviswanathan algorithm

The algorithms we introduced in the previous section represent the baseline algo-
rithms for finding the pattern embeddings. However, they can be very costly for
larger patterns and domain graphs. Next, we propose an extension of a theoreti-
cal sampling idea suggested by Fürer and Kasiviswanathan (2014) for unlabeled
undirected graphs.

Figure 3 depicts one iteration of the Fürer and Kasiviswanathan algorithm. Each iteration of the algorithm returns a single embedding and an estimate for the pattern's total number of embeddings. The key intuition behind the algorithm is that overall problem of estimating the number of embeddings of a pattern can be broken down into the problem of estimating the number of embeddings for each partition $V_i$ of the OBD. A single iteration sequentially goes through all of the OBD's partitions and builds up both an embedding of the pattern and an estimate of the pattern's total number of embeddings. For each partition, it finds all possible extensions of the current partial embedding by matching the vertices in the partition. It then randomly samples one of these extensions and multiplies its estimate of the pattern's number of embeddings by the number of possible extensions found for the current partition.

When searching for an OBD, the crucial step is to ensure that each partition is small. Since all vertices in a partition must be simultaneously matched, it is easier to find a legal extension when fewer vertices must be matched. If the decomposition provided as input is not an ordered bipartite decomposition, the sampling algorithm will still converge to the correct value, but the convergence speed guaranteed by Fürer and Kasiviswanathan (2014) no longer holds.

<figure><img src="image_4.png" /><figcaption>Fig. 3: An illustration of one iteration of the Fürer and Kasiviswanathan approach. The panel labeled I shows how one node is selected from $N_1$ total possibilities to match the vertex in the first partition. The panel labeled II shows how one node is selected from $N_2$ total possibilities to match the vertex in the second partition, and the new count is $N_1 \times N_2$. Finally, the panel labeled III shows how two nodes are selected from $N_3$ total possibilities to match the two vertices in the third partition, and the final count is $N_1 \times N_2 \times N_3$.</figcaption></figure>

Algorithm 4 shows our extension of the Fürer and Kasiviswanathan approach to labeled graphs. The algorithm consists of two parts: the control loop (lines 2-7), and the single attempt procedure (lines 9-18). It receives the same parameters as Algorithm 3, except that it receives a decomposition $V$ instead of an ordering $O$. While the time limit is not exceeded, the control loop randomly samples, with replacement, a path in the search space that possibly leads to an embedding. The
<caption>**Algorithm 4** Fürer and Kasiviswanathan approach.</caption>

1: **function** FK-V(domain graph D, pattern subgraph P, decomposition V)
2:     **M** ← {}
3:     **while** not timeout **do**
4:         **M** ← FK-V(1, D, P, V, 1, {}
5:         **if** **M** ≠ ∅ **then**
6:             **M** ← **M** ∪ {M}
7:             **return** **M**
8:         
9:         **function** FK-V(W, D, P, V, i, φ)
10:         **if** i > |V| **then**
11:             **return** ({(φ, W)})
12:         **else if** timeout **then**
13:             **return** {}
14:         **C** ← PARTIALEXT(D, P, φ, Vᵢ)
15:         **if** **C** = {} **then**
16:             **return** {}
17:         **Select** φ' **randomly from** **C**
18:         **return** FK-V(W × |C|, D, P, V, i + 1, φ')
    ▷ If all partitions embedded
    ▷ Find all extensions for the partition
    ▷ Randomly select one extension

function FK-V returns a set of pairs (φ, W) where φ is the embedding found and W is a positive number equal to the inverse probability that an attempt will find that particular embedding. The expected value of this W (assuming implicitly 0 for the unsuccessful attempts) is the total number of embeddings.

A call to FK-V identifies C, the set of all possible extensions to the current partial embedding. It randomly samples a set of images for a partition Vᵢ and extends the partial embedding accordingly. While it only samples one element from C, it uses the size of C as an estimate for the number of subpattern represented by Vᵢ. It then recurses, with both the extended partial embedding φ' and an updated count estimate of W × |C| (which becomes 0 when the attempt fails).

For example, if $V_i = \{w1, w2, w3\}$, each vertex $w \in V_i$ has a set of candidate matches $C_w$ in D. Suppose these are $C_{w1} = \{a1, a2\}$, $C_{w2} = \{a3\}$ and $C_{w3} = \{a1, a4, a5\}$. Then each element of the Cartesian product $C_{w1} \times C_{w2} \times C_{w3}$ is a valid match for the vertices in $V_i$ if the matching operator is homomorphism. In this case, $|C|$ would equal 6. In the case of subgraph isomorphism, embeddings mapping several vertices to the same image (here $\{(w1, a1), (w2, a3), (w3, a1)\}$) must be eliminated, so $|C|$ would equal 5.

## 7 Computing ordered bipartite decompositions

When estimating the number of embeddings for pattern P, Algorithm 4 requires an ordered bipartite decomposition (OBD) for P. For efficiency, the size of each partition class in an OBD is important. Having more vertices in a single partition requires matching multiple vertices simultaneously, which is more computationally expensive. In particular, the size of set C in line 14 of Algorithm 4 may be exponential in $V_i$. This implies that we want to minimize the number of vertices in each partition class. This also suggests that we prefer having OBDs that are as large as possible (i.e., have as many partition classes as possible). Fürer and Kasiviswanathan (2014) assume that a pattern's OBD is given. Thus, we propose two different search strategies for automatically finding an OBD.
## 7.1 Exact search

Algorithm 5 presents a level-wise search for finding the largest possible OBD (i.e., the one with the most partition classes). First, it checks if *P* contains a triangle, because then no OBD exists. The maximal OBD size for pattern *P* is |*V*(*P*)|, which is the number of vertices in *P*. Starting with the largest possible size *s* = |*V*(*P*)|, the algorithm performs a search through the space of possible OBDs. In each iteration, it checks whether an OBD of size *s* exists by generating all possible candidate OBDs of size *s*. This entails enumerating all ways of dividing *V*(*P*) into *s* partition classes, and then building one candidate for each possible ordering of the partition classes. It then applies the VALID_OBD function to check if each candidate is a valid OBD (i.e., it satisfies the requirements of Definition 1). As soon as it finds a valid OBD, it is returned and the search terminates. If no OBD of size *s* is found, then *s* is decremented, and the search proceeds. If an OBD is found, its size is maximized because all larger decompositions were evaluated and none were valid OBDs.

<caption>**Algorithm 5** Level-wise exact OBD search</caption>

1: **function** FINDLARGESTOBD(graph *P*)
   **return** Largest OBD for *P* or ∅ if no OBD exists
2:   **if** CONTAINS_TRIANGLE(*P*) **then**
3:     **return** ∅
4:   **s** ← |*V*(*P*)|  
5:   **while** *s* > 0 **do**
6:     **for all** ordered partitions *V* of *V*(*P*) of size *s* **do**
7:       **if** VALID_OBD(*V*) **then**
8:         **return** *V*
9:       **s** ← *s* − 1
10:      **return** ∅

## 7.2 Greedy search

Algorithm 6 presents a depth-first search for finding an ODB for a pattern *P*. It greedily searches the most promising path and returns the first OBD found. Because it performs backtracking, it is guaranteed to find an OBD if one exists. However, its greedy nature means that it is not guaranteed to find the largest valid OBD.

The main algorithm loops through each node *n* in *P*, and builds a partial OBD (*n*) consisting of a single partition class containing *n*. (If there exists an OBD, then there exists one with a singleton as first partition class.) Next, it calls GREEDYOBD which recursively extends the partial OBD by adding one partition class at a time.

If GREEDYOBD receives a valid and complete OBD (i.e., it contains all the nodes in *P*), it terminates and returns the current OBD. Otherwise, it generates a set of candidate extensions for *V*. To do so, it first builds a graph *G'* = (*R*, *L*), which has one node for each vertex in *P* that has not been assigned a partition class in the partial OBD *V*. It has one edge {*x*, *y*} ∈ *L* for each *x*, *y* ∈ *R* such that
$x$ and $y$ share a neighbor in $P$ that already appears in $\mathcal{V}$. This implies that $x$ and $y$ must be in the same partition class. Because a node $x \in R$ may neighbor multiple nodes in $\mathcal{V}$, the set of candidate extensions $\mathcal{C}$ consists of each connected component in $G'$. Next, line 14 checks if any candidate partition requires assigning the same label to neighboring vertices in $P$, which would violate the definition of an OBD, meaning no complete OBD can be constructed by extending this partial OBD. In this case, the algorithm backtracks. Otherwise, it creates a set of candidate extensions by concatenating each element (i.e., set of vertices) in $\mathcal{C}$ as a partition class to the end of $\mathcal{V}$. It sorts the set of candidates from smallest to largest on the basis of the size of the newest partition class. This heuristic is an attempt to favor OBDs with small partition classes. It then recursively calls GREEDYOBD on each candidate in this order.

<caption>**Algorithm 6** Depth-first greedy OBD search.</caption>

1: **function** GREEDYFINDOBD(graph $P$)
   **return** OBD for $P$ or $\emptyset$ if no OBD exists
2:   **for** all vertices $v \in V(P)$ **do**
3:     $result$ = GREEDYOBD({$v$}, $P$)
4:     **if** $result \neq \emptyset$ **then**
5:       **return** $result$
6:     **return** $\emptyset$
7: **function** GREEDYOBD(partial OBD $\mathcal{V}$, graph $P$)
8:   **if** COMPLETE($\mathcal{V}$) && VALID_OBD($\mathcal{V}$) **then**
9:     **return** $\mathcal{V}$
10:    **else**
11:     ▷ Vertices not yet in partition class
12:     $L$ = {$x, y \in R | \exists V \in \mathcal{V} : N_P(x) \cap V \neq \emptyset \land N_P(y) \cap V \neq \emptyset$}
13:     $C$ = ConnectedComponents($R, L$)
14:     **if** $\exists C \in \mathcal{C}, \exists x, y \in V(C) : \{x, y\} \in E(P)$ **then**
15:      ▷ $x$ and $y$ must be apart ⇒ no OBD possible
16:      **return** $\emptyset$
17:    **for** all $C \in \mathcal{P}$ **do**
18:      $result$ = GREEDYOBD($\mathcal{V} \cdot V(C), P$)
19:      **if** $result \neq \emptyset$ **then**
20:        **return** $result$
21:    **return** $\emptyset$

## 7.3 Empirical evaluation of OBD search

We empirically compare the run time performance of the two OBD search algorithms. We generate all simple connected graphs with between five and eight vertices and use both Algorithm 5 and Algorithm 6 to find an OBD for each generated pattern. Table 1 reports the mean, median and maximum run time in seconds for both algorithms on each pattern size. While the depth-first approach does not guarantee finding an optimal solution, on average it is several orders of magnitude faster than the level-wise approach. As it scales better than the exact approach,
our empirical evaluation in Section 9 employs the depth-first search approach to find an OBD for each evaluated pattern.

<table><caption>Table 1: This table reports the mean, median and maximum run times in seconds for the exact and depth-first search algorithms for finding an OBD on all simple connected graphs with between five and eight vertices.</caption><thead><tr><th rowspan="2">|V|</th><th colspan="3">Exact</th><th colspan="3">Greedy</th></tr><tr><th>Mean</th><th>Median</th><th>Max</th><th>Mean</th><th>Median</th><th>Max</th></tr></thead><tbody><tr><td>5</td><td>0.0014</td><td>&lt;0.0001</td><td>0.0096</td><td>0.0002</td><td>0.0003</td><td>0.0004</td></tr><tr><td>6</td><td>0.0106</td><td>&lt;0.0001</td><td>0.1611</td><td>0.0004</td><td>0.0004</td><td>0.0010</td></tr><tr><td>7</td><td>0.0618</td><td>&lt;0.0001</td><td>2.8392</td><td>0.0005</td><td>0.0005</td><td>0.0025</td></tr><tr><td>8</td><td>0.3894</td><td>&lt;0.0001</td><td>54.9054</td><td>0.0008</td><td>0.0007</td><td>0.0084</td></tr></tbody></table>

# 8 Computing statistics

The algorithms described in the previous sections all generate a sample of embeddings. In this section we will discuss how we can compute a number of statistics from such samples, and we will study their accuracy.

## 8.1 Common statistics

Let $P$ be a pattern and $D$ a domain graph. We first define the statistics we will consider.

**Count.** The count statistic just returns $|Val(P, D)|$, the number of embeddings of $P$ in $D$.

**Mean and variance** Let $t: Range(P, D) \to \mathbb{R}$ be a function mapping values of embeddings on real numbers. Then, $\mu_t(X) = \sum_{x \in X} t(x)/|X|$ and $var_t(X) = \sum_{x \in X} (t(x) - \mu_t(X))^2/|X|$.

**Confusion matrix.** Some vertices of the pattern could be labeled by finite sets. In that case, it may be interesting to compare the relative frequencies with which these vertices map on vertices with certain labels.

More formally, let $S \subseteq V(P)$ be a set of vertices of $P$. For $x \in Range(S: P, D)$ a possible value of the vertices in $S$ and $y \in Val(P, D)$ the value of a particular embedding, we define $I_{S \to x}(y) = 1$ if $y|_S = x$ and $I_{S \to x}(y) = 0$ if $y|_S \neq x$. Next, for a set of values $X$, $freq_{S \to x}(X) = \mu_{I_{S \to x}}$ and hence $freq_{S \to x}(P, D)$ is the fraction of embeddings of $P$ in $D$ for which the vertices in $S$ get the value $x$. We define $freq_S(P, D) = \{(x, freq_{S \to x}(P, D) | x \in Range(S: P, D))\}$. $freq_S(P, D)$ is a confusion matrix (or tensor), as for each possible combination of values for the vertices in $S$ it gives the observed (relative) frequency.

## 8.2 Estimators

We now discuss how to compute statistics from samples. In the above, we defined the elements of the confusion matrix statistic using a mean statistic, so we can
limit our discussion to *count*, *μ* and *var* here. In the empirical evaluation section, we will focus more on the confusion matrix statistic.

Let us first consider the random sampling. It returns a set $\mathcal{M} = \{M_i\}_{i=1}^n$ of groups of embeddings with $n$ the number of main loop iterations of the algorithm. Each group $M_i = \{\varphi_{i,j}\}_{j=1}^{n_i}$ of embeddings is generated together during one iteration of random and shares the same image for the first pattern vertex. We can estimate count($P, D$) by count<sup>rnd</sup>($P, D$) = $|V(P)| \sum_{i=1}^n n_i/n$. We can estimate $\mu_t(P, D)$ by $\hat{\mu}_t^{rnd} = \sum_{i=1}^n \sum_{j=1}^{n_i} t(val_D(\varphi_{i,j}))/\sum_{i=1}^n n_i$. This is a consistent estimator (it converges to $\mu_t(P, D)$ when $n$ gets large) but even though count<sup>rnd</sup> and $\sum_{i=1}^n \sum_{j=1}^{n_i} t(val_D(\varphi_{i,j}))$ are unbiased estimators of the count and the sum, $\hat{\mu}_t^{rnd}$ is not an unbiased estimator itself. For instance, suppose that the first vertex has two possible images $u_1$ and $u_2$, $u_1$ is in three embeddings with $t$-values 1, 2 and 3 and $u_2$ is in one embedding with $t$-value 6. Then, a single iteration of the sampling algorithm will either produce values 1, 2 and 3 leading to an estimation $(1+2+3)/3 = 2$ or a single value 6, leading to an estimate of 6. On average, we expect an estimation of 4, while the average of all values is 3.

Next, we consider the FK-V algorithm. It returns a set $M = \{(\varphi_i, w_i)\}_{i=1}^n$ of independently sampled weighted embeddings. Fürer and Kasiviswanathan (2014) showed that count<sup>fk</sup>($P, D$) = $\sum_{i=1}^n w_i/n$ is an unbiased estimator. For $\mu_t$ we can derive an unbiased estimator too (but omit details due to lack of space). Let $W = \sum_{i=1}^n w_i$ and $W_2 = \sum_{i=1}^n w_i^2$. Then, $\hat{\mu}_t^{fk}(P, D) = \sum_{i=1}^n w_i t(val_D(\varphi_i))/W$ is an unbiased estimator of $\mu_t(P, D)$. Moreover, we can estimate the square error from the sample with $E[\hat{\mu}_t^{fk}(P, D) - \mu_t(P, D)] = (W_2/W.(W^2 - W_2)) \sum_{i=1}^n w_i(t(val_D(\varphi_i)) - \hat{\mu}_t^{fk})^2$. It is straightforward to show that these estimators remain unbiased even if Fürer-Kasiviswanathan is run with a decomposition which is not an OBD. The main difference is that in that case the error converges more slowly to zero as the sample size increases.

# 9 Experiments and results

The goal of the empirical evaluation is to address the following questions:

**Q1:** How well do the sampling algorithms approximate the true embedding statistics obtained using the exhaustive approach?

**Q2:** How does having an OBD versus an arbitrary decomposition (AD) affect the performance of the Fürer-Kasiviswanathan algorithm?

**Q3:** How does the method of determining an arbitrary decomposition affect the performance of the Fürer-Kasiviswanathan algorithm?

**Q4:** How does the Fürer-Kasiviswanathan algorithm perform on random and non-random graphs?

**Q5:** How does using the Fürer-Kasiviswanathan algorithm to compute approximate counts affect weight learning and inference in statistical relational learning?

**Q6:** How well does the Fürer-Kasiviswanathan algorithm perform on longer patterns?

The first two questions attempt to assess the real-world viability of the algorithm described in this paper, while the last three questions attempt to provide more insight into the Fürer and Kasiviswanathan approach.
The majority of the empirical evaluation focuses on comparing the performance of the following five approaches:

**Exhaustive:** This is the approach of Algorithm 2.

**Random:** This is the approach of Algorithm 3 that randomly samples an initial vertex in the domain graph from which to start the search for the embeddings of a pattern.

**FACT:** This is the approach proposed by Das et al (2016). This method approximately counts the satisfied instances by using the summary statistics of in- and out-degrees of the graph representing the database.

**FK-OBD:** This is the extended Fürer and Kasiviswanathan approach outlined in Algorithm 4 that uses depth-first search to find an OBD for a pattern. It only applies to generated patterns that have an OBD.

**FK-AD:** This is the extended Fürer and Kasiviswanathan approach outlined in Algorithm 4 that works with an arbitrary decomposition of a pattern. To construct an arbitrary decomposition, we simply flatten the OBD, which means that each partition contains exactly one node.

Next, we describe our datasets, the experimental methodology, and then we present and discuss the results.

## 9.1 Datasets

We evaluate the algorithms using 12 synthetic datasets and seven real-world datasets. Our real-world datasets are initially stored in a relational format. We use a standard transformation (e.g., Richards and Mooney (1992)) to convert these datasets into a graph form. The general idea is that graphs consist of constants connected through their relations.

### 9.1.1 Synthetic datasets

We created a series of synthetic datasets that contain the following node types: *person*, *salary*, *satisfaction*, *man*, *woman*, *married* and *friends*. We create two types of graphs by varying the way friend nodes are connected within the graph. The first denoted with *S* adds the nodes according to the power law distribution of the Barabási-Albert algorithm (Barabasi and Albert 1999). The second, denoted with *R*, employs the Erdős-Rény random graph model, and adds them randomly. Finally, each person receives a salary and a satisfaction, 40% of the people are male, the rest are female, and 70% of the people are married. Table 2 shows the synthetic graphs and their characteristics.

### 9.1.2 Real-world datasets

We consider a standard bibliography dataset, three well-known relational learning datasets, and three datasets from the Stanford Network Analysis Platform (SNAP) dataset collection.<fnref n="3" /> For the real-world datasets, we transformed directed or labeled edges into 2-paths with an intermediate labeled vertex such that our graphs only

<footnote n="3">³ https://snap.stanford.edu/data/</footnote>
<table><caption>Table 2: Characteristics of the synthetic graphs, which are sorted in descending order by density. The maximum degree (MaxDeg) is the maximum number of outgoing edges of the nodes, and the density is calculated as $Density = \frac{2m}{n(n-1)}$, where *m* is the number of edges and *n* is the number of nodes in the graph.</caption><thead><tr><th rowspan="2">#Nodes</th><th rowspan="2">#Edges</th><th rowspan="2">#Friends</th><th colspan="2">Power Law</th><th colspan="2">Random</th><th rowspan="2">Density</th></tr><tr><th>Name</th><th>MaxDeg</th><th>Name</th><th>MaxDeg</th></tr></thead><tbody><tr><td>143,445</td><td>236,890</td><td>99,945</td><td>S1</td><td>551</td><td>R1</td><td>32</td><td>1.56 × 10<sup>-5</sup></td></tr><tr><td>24,140</td><td>43,280</td><td>19,790</td><td>S2</td><td>383</td><td>R2</td><td>29</td><td>3.13 × 10<sup>-5</sup></td></tr><tr><td>14,295</td><td>23,590</td><td>9,945</td><td>S3</td><td>261</td><td>R3</td><td>26</td><td>6.26 × 10<sup>-5</sup></td></tr><tr><td>93,485</td><td>136,970</td><td>49,985</td><td>S4</td><td>158</td><td>R4</td><td>26</td><td>1.61 × 10<sup>-4</sup></td></tr><tr><td>243,290</td><td>436,580</td><td>199,790</td><td>S5</td><td>174</td><td>R5</td><td>36</td><td>2.31 × 10<sup>-4</sup></td></tr><tr><td>9,335</td><td>13,670</td><td>4,985</td><td>S6</td><td>131</td><td>R6</td><td>24</td><td>3.12 × 10<sup>-4</sup></td></tr></tbody></table>

have vertex labels. Next, we describe the datasets in more detail. Table 3 lists the properties of the real-world datasets we used.

**DBLP dataset**<fnref n="4" /> is a bibliography database from the ArnetMiner repository (Tang et al 2008). We only consider data from 2001 and 2002, to ensure that the Exhaustive approach is computationally feasible. The DBLP domain graph consists of *paper*, *coauthored*, *reference* and *citations* nodes.

<table><caption>Table 3: Characteristics of the real-world graphs. The average degree represents the average number of outgoing edges of the nodes, the maximum degree is the maximum number of outgoing edges of the nodes, and the density is calculated as $Density = \frac{2m}{n(n-1)}$, where *m* is the number of edges and *n* is the number of nodes in the graph. Note that the number of nodes in our SNAP datasets is larger than specified on their website because in our representation both objects (i.e., users, e-mails) and their attributes appear as nodes.</caption><thead><tr><th>Dataset</th><th>#Nodes</th><th>#Edges</th><th>Average Degree</th><th>MaxDeg</th><th>Density</th></tr></thead><tbody><tr><td>DBLP</td><td>393,230</td><td>447,650</td><td>2.28</td><td>1,036</td><td>5.7e-06</td></tr><tr><td>Yeast</td><td>16,233</td><td>18,355</td><td>2.26</td><td>124</td><td>0.00014</td></tr><tr><td>IMDB</td><td>548</td><td>782</td><td>2.85</td><td>29</td><td>0.005</td></tr><tr><td>WEBKB</td><td>1,477</td><td>61,766</td><td>83.63</td><td>286</td><td>0.057</td></tr><tr><td>Facebook (SNAP)</td><td>28,057</td><td>112,252</td><td>8.00</td><td>1,051</td><td>0.0003</td></tr><tr><td>Enron e-mails (SNAP)</td><td>36,692</td><td>183,831</td><td>10.02</td><td>1,383</td><td>0.0003</td></tr><tr><td>Amazon (SNAP)</td><td>334,863</td><td>925,872</td><td>5.53</td><td>549</td><td>1.65e-05</td></tr></tbody></table>

The following three datasets are commonly used in relational learning

**Yeast** is about protein-protein interactions in yeast (Mewes et al 2000; Davis et al 2005). This paper uses a modified version (Van Haaren et al 2015) from the Alchemy repository, such that all self-referential links are removed.<fnref n="5" /> The graph has the following node types: *protein*, *location*, *function*, *phenotype*, *class*, *enzymes*, and *interaction*.

<footnote n="4">⁴ http://www.informatik.uni-trier.de/ley/db/index.html</footnote>

<footnote n="5">⁵ http://alchemycs.washington.edu</footnote>
<figure><img src="image_5.png" /><figcaption>Fig. 4: Pattern graph $P_{sat-sal}$ for our example. Shaded nodes are the target nodes of the pattern.</figcaption></figure>

**IMDB** consists of *movies* which are described by a *genre*, the movie's *directors*, and the *actors* who appeared in the movie. Actors and directors have a specific *gender*.

**WEBKB** consists of labeled *pages* and their hyperlinks from computer science departments of four universities.

The following three datasets come from SNAP:

**Facebook** consists of *users* from Facebook connected through *friendship* relations and where each user has a number of anonymized features such as *birthday*, *gender*, *education type*, *education degree*, *hometown*, and *language*.

**Enron** is a communication network consisting of email addresses such that if an address *i* sent at least one email to address *j*, the graph contains an undirected edge from *i* to *j*.

**Amazon** is a network consisting of products from Amazon such that if a product *i* is frequently co-purchased with product *j*, the graph contains an undirected edge from *i* to *j*.

## 9.2 Experimental methodology

Next, we introduce our experimental methodology which involves generating and evaluating patterns.

### 9.2.1 Generating pattern graphs

**Synthetic pattern graph** For the synthetic dataset we only consider one pattern graph, which is depicted in Figure 4. This pattern is a slight modification of a phenomena described by Ariely (2008).

**Generating pattern graphs for real-world data** We want to evaluate our approach on a variety of patterns. However, even with a restriction on the pattern size, it is computationally too demanding to evaluate all possible candidate patterns. Hence, we generate patterns of increasing size in a level-wise manner, and evaluate a subset of patterns of each size.

Each pattern is an undirected graph. We vary the pattern size from 4 to 10 nodes. Given a set $C_L$ of patterns of size $L$, we generate a candidate set $C_{L+1}$ of patterns of size $L+1$ by enumerating all valid extensions of each pattern in $C_L$. Each pattern is extended by first adding a new node to it. We consider adding all node types. Then, new candidates are created for all ways of adding up to
$m = n \cdot (n - 1)/2$ edges between the new node and $n$ existing nodes in the pattern. We remove duplicates from the candidate set by checking for isomorphism between pairs of candidate patterns.

We focus our evaluation on a sample of patterns where each selected pattern is neither too frequent nor too infrequent as these tend to be more interesting patterns for analysis. Thus, we randomly sample patterns from $C_{L+1}$ and run FK-OBD or FK-AD for one hour on the larger datasets, and for five minutes for smaller datasets. FK-OBD is used if a pattern has an OBD, and FK-AD if it does not. We continue sampling until we find 100 patterns per each pattern size that meet the following criterion:

$$ \sqrt{M} - 3 \cdot N_{\text{std}} \le N_{\text{avg}} \le M + 3 \cdot N_{\text{std}}, $$

where $M$ denotes the number of nodes in the data graph, $N_{\text{avg}}$ denotes the average and $N_{\text{std}}$ the standard deviation of the estimated number of embeddings after the one-hour run.

### 9.2.2 Pattern evaluation

We evaluate all algorithms on the selected patterns for each size. We run each algorithm for a maximum of $T$ time units and record intermediate results after each $K$ time units. For DBLP, Enron, and Amazon, we set $T$ to ten hours and $K$ to five minutes. For Yeast, IMDB, WEBKB, Facebook, and the synthetic datasets, $T$ is ten minutes and $K$ is five seconds. This means that we record 120 intermediate results for each pattern, regardless of the dataset or method.

### 9.2.3 Evaluation tasks and metrics

We consider two different evaluation tasks.

**Frequency estimation** This corresponds to the standard pattern mining task, where the support of a pattern in the data is calculated.

**Parameter estimation** Parameter learning for statistical relational (SRL) (Getoor and Taskar 2007) formalisms such as Markov logic networks (MLNs) (Richardson and Domingos 2006), Bayesian logic programs (BLPs) (Kersting et al 2000), and logical Bayesian networks (LBNs) (Fierens et al 2004), is a widely studied problem (e.g., Venugopal et al (2015); Das et al (2016); Ravkic et al (2015)). Addressing this task requires estimating the frequency of (many) patterns in the data. In this paper, we explore the suitability of using the proposed approaches for parameter estimation in LBNs.

For frequency estimation of a pattern, we report the relative error of a pattern's estimated frequency (using a sampling algorithm) up until time $t$, denoted as $\#\hat{emb}_i(t)$, versus its true frequency as computed by the exact approach. The relative error is given by the following formula:

$$ RelErr_i(t) = \frac{|\#emb_i - \#\hat{emb}_i(t)|}{\#emb_i} \qquad (1) $$

Parameter estimation in LBNs requires estimating many conditional probability distributions (CPD) of the form $P(Y | X \leftarrow Context)$. In these CPDs, $Y$ is
a single random variable that can take on multiple values, X is a set of random variables and Context defines a set of conditions that specifies when the CPD is valid (e.g., there only exists a random variable for a student's grade in a course C if she has taken course C.).<fnref n="6" /> The CPD captures the conditional probability of Y for each joint assignment of values to X, when Context holds in the data. Thus, in a constructed pattern we must specify which nodes represent each of X, Y, and Context. For example, in Figure 4 the shaded nodes represent X, Y (target) nodes, while other nodes represent Context. For example, when creating the patterns for the DBLP dataset, X, Y are nodes labeled with citations. Nodes labeled as references, dir, coauthored, citations = low, citations = med, citations = high represent Context. For the patterns in the Yeast dataset, X, Y can be nodes labeled with function, location, class, enzyme or phenotype. Each of these has a range of values. Nodes labeled with interaction and predicate = value where predicate ∈ {function, location, class, enzyme, phenotype} and value ∈ range(predicate) are Context. This procedure is similar for the other datasets.

To make the results less susceptible to fluctuations arising from infrequently occurring combinations, we employ the standard Laplace smoothing of our estimates. Furthermore, we aggregate each value combination that has a relative frequency of < 1% in the Exhaustive approach into a single “default” row. In probabilistic models, this is known as using a default table representation for the local probability distributions (Friedman and Goldzsmidt 1996).

In order to compare the quality of the estimated and the true distribution for a CPD we calculate the Kullback-Leibler divergence (KLD).

$$KLD(P||Q) = \sum_{i} p(i) \times \log_{2} \left( \frac{p(i)}{q(i)} \right) \quad (2)$$

The estimated distribution (P) is obtained by one of the sampling algorithms, and the true distribution (Q) is obtained by the Exhaustive approach.

## 9.3 Experiments and results for Q1 and Q2

Ultimately, one of the primary goals of sampling is to achieve similar performance at a fraction of the run time required by an Exhaustive approach. We consider two possible ways to vary the run time of the sampling approaches: by comparing them to the Exhaustive approach's and FACT's run times. In both scenarios, we focus on estimating the counts. Hence, we only consider patterns where we know the true count (i.e., those patterns where the Exhaustive approach finished within the time limit). Table 4 provides for each pattern size the percentage of patterns where the Exhaustive approach does not finish. This analysis does not consider the Amazon and Enron datasets as the Exhaustive approach never completes within the time limit. Table 5 shows the average run times for the Exhaustive approach and FACT for each pattern size.

First, we allow each sampling approach to run for a percentage of the Exhaustive approach's total run time. For example, if the Exhaustive approach finished in

<footnote n="6">⁶ Note that this is a slight simplification as LBN uses first-order logic to perform parameter tying across multiple random variables.</footnote>
<table><caption>Table 4: For each pattern size, the percentage of patterns where the Exhaustive approach failed to finish within the time limit on every dataset.</caption><thead><tr><th rowspan="2">Dataset</th><th colspan="10">Pattern Size</th></tr><tr><th>4</th><th>5</th><th>6</th><th>7</th><th>8</th><th>9</th><th>10</th><th>Average</th></tr></thead><tbody><tr><td>Yeast (10 minutes)</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr><tr><td>WEBKB</td><td>0.0</td><td>9.1</td><td>30.0</td><td>36.7</td><td>58.1</td><td>77.0</td><td>80.0</td><td>41.5</td></tr><tr><td>IMDB</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr><tr><td>Facebook</td><td>10.7</td><td>47.0</td><td>69.0</td><td>94.0</td><td>96.0</td><td>96.0</td><td>98.1</td><td>73.1</td></tr><tr><td>DBLP (10 hours)</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td><td>3.0</td><td>1.0</td><td>10.0</td><td>2.0</td></tr><tr><td>Amazon</td><td>100.0</td><td>100.0</td><td>100.0</td><td>100.0</td><td>100.0</td><td>100.0</td><td>100.0</td><td>100.0</td></tr><tr><td>Enron</td><td>100.0</td><td>100.0</td><td>100.0</td><td>100.0</td><td>100.0</td><td>99.0</td><td>100.0</td><td>99.9</td></tr></tbody></table>

5 hours for a pattern, then 10% represents estimated results for FK-OBD and FK-AD after 30 minutes. Moreover, we focus on more challenging patterns by omitting any pattern where FK-AD had less than 10% relative error on the estimated number of embeddings in the first interval (that is, the first 5 seconds for YEAST, IMDB, FACEBOOK and WEBKB, and the first 5 minutes for DBLP). Figure 5 shows how the average relative error (averaged over all pattern sizes) varies as a function of run time for the FK-OBD, FK-AD and Random approaches. The plots show that the FK-OBD and FK-AD perform similarly and have a better convergence than Random for all the datasets. Figure 6 shows similar plots but for three specific pattern sizes: 4, 7 and 10.<fnref n="7" /> Regardless of the pattern size, both FK-OBD and FK-AD converge very quickly, typically in around 10% of the run time needed by the Exhaustive approach, to a stable estimate. There is a slight decrease in accuracy as the pattern size increases. Finally, regardless of the pattern size, both FK approaches are better than Random, and the gap between them and Random grows as the pattern size increases.

Second, we measure FACT's run time for each pattern and then allow each sampling algorithm (FK-OBD, FK-AD and Random) to run for the same amount of time. Table 6 shows the average relative errors of the sampling approaches and FACT on four real-world datasets. FACT exceeded the amount of available memory<fnref n="8" /> on DBLP, Amazon, and Enron. These results also show that FK-AD and FK-OBD have similar performance. On three datasets, they perform better than Random, and on one dataset their performance is equivalent. FK-OBD, FK-AD, and Random are consistently much more accurate than FACT. Table 7 shows the average relative error for each of these four approaches on a per pattern size basis on each dataset for this setting. On the Yeast, WEBKB, and IMDB datasets, FK-OBD and FK-AD result in substantially better performance than FACT, with reductions in the average relative error between 57.5% to 98% on pattern sizes up to length seven. On larger pattern sizes, there are reductions of between 25% and 94%. On Facebook, particularly for pattern sizes nine and ten, the results are closer, but FK-OBD and FK-AD are always better on average. Hence, both FK approaches seem to offer superior performance than FACT for estimating counts.

<footnote n="7">⁷ We omit FACT on these plots to declutter them and because its results would simply be a point.</footnote>

<footnote n="8">⁸ The experiments were run on a machine with 10Gb of RAM.</footnote>
<table><caption>Table 5: The average run time in seconds for each pattern size for FACT and the Exhaustive approach.</caption><thead><tr><th rowspan="2">Dataset</th><th rowspan="2">Method</th><th colspan="10">Pattern Size</th></tr><tr><th>4</th><th>5</th><th>6</th><th>7</th><th>8</th><th>9</th><th>10</th><th>Average</th></tr></thead><tbody><tr><td rowspan="2">Yeast</td><td>Exhaustive</td><td>3.50</td><td>5.24</td><td>10.10</td><td>5.12</td><td>10.27</td><td>21.54</td><td>24.95</td><td>11.53</td></tr><tr><td>FACT</td><td>16.50</td><td>17.77</td><td>14.69</td><td>14.14</td><td>17.45</td><td>14.31</td><td>17.78</td><td>16.09</td></tr><tr><td rowspan="2">WEBKB</td><td>Exhaustive</td><td>1.08</td><td>96.32</td><td>321.44</td><td>338.03</td><td>505.39</td><td>561.20</td><td>565.10</td><td>341.22</td></tr><tr><td>FACT</td><td>7.00</td><td>7.23</td><td>7.22</td><td>7.37</td><td>7.48</td><td>7.42</td><td>7.83</td><td>7.36</td></tr><tr><td rowspan="2">IMDB</td><td>Exhaustive</td><td>0.00</td><td>0.03</td><td>1.59</td><td>15.18</td><td>25.60</td><td>25.37</td><td>16.86</td><td>12.09</td></tr><tr><td>FACT</td><td>4.34</td><td>4.01</td><td>3.99</td><td>3.99</td><td>4.10</td><td>4.16</td><td>4.01</td><td>4.09</td></tr><tr><td rowspan="2">Facebook</td><td>Exhaustive</td><td>69.90</td><td>447.36</td><td>657.18</td><td>614.36</td><td>625.59</td><td>609.16</td><td>608.44</td><td>518.86</td></tr><tr><td>FACT</td><td>34.00</td><td>34.02</td><td>35.88</td><td>31.30</td><td>41.69</td><td>27.39</td><td>24.91</td><td>32.74</td></tr></tbody></table>

<figure><img src="image_6.png" /><figcaption>Fig. 5: The average relative error as a function of the percentage of the Exhaustive approach's run time for the FK-OBD, FK-AD and Random approaches on the Yeast, IMDB, Facebook, WEBKB and DBLP datasets. The results are averaged over all pattern sizes.</figcaption></figure>

Next, we focus on the parameter estimation task. We only consider the Yeast, IMDB, DBLP, Facebook, and WEBKB datasets because Amazon and Enron only contain links (no attributes), so the conditional probability tables would be less interesting. Again, we omit those patterns where the Exhaustive approach failed to finish within the time limit as we lack the ground truth answers for them.

Table 8 presents the KLD averaged over all pattern sizes and all patterns. Here, we see that FK-OBD and FK-AD obtain highly similar results. On the Yeast, Facebook, WEBKB, and DBLP datasets, FK-OBD and FK-AD do substantially
<figure><img src="image_7.png" /><figcaption>Fig. 6: The average relative error for a specific pattern size in relation to the percentage of the Exhaustive approach's run time for FK-OBD, FK-AD and Random on the Yeast, IMDB, Facebook, WEBKB and DBLP datasets. The plots show results for pattern sizes 4 (left column), 7 (middle column), and 10 (right column).</figcaption></figure>

<table><caption>Table 6: The average relative error and its standard deviation for the FK-OBD, FK-AD, FACT and Random approaches. The average is over all pattern sizes on the Yeast, WEBKB, IMDB, and Facebook datasets. The best results are shown in bold.</caption><thead><tr><th>Dataset</th><th>FK-OBD</th><th>FK-AD</th><th>FACT</th><th>Random</th></tr></thead><tbody><tr><td>Yeast</td><td><strong>0.048±0.061</strong></td><td>0.054±0.072</td><td>0.917±0.191</td><td>0.184±0.190</td></tr><tr><td>WEBKB</td><td>0.456±0.376</td><td><strong>0.432±0.382</strong></td><td>0.872±0.245</td><td>0.592±0.332</td></tr><tr><td>IMDB</td><td>0.124±0.011</td><td><strong>0.077±0.017</strong></td><td>0.861±0.268</td><td>0.298±0.061</td></tr><tr><td>Facebook</td><td><strong>0.605±0.375</strong></td><td>0.607±0.376</td><td>0.999±0.002</td><td>0.606±0.452</td></tr></tbody></table>

better (i.e., at least a 45% reduction in the KLD) than the Random approach. On the IMDB dataset Random does the best, but all three approaches have similar results. Table 9 shows the average KLD for each pattern size. The KLD is relatively stable for all approaches with increases becoming apparent around pattern size seven or eight. At this size, there is spike that is most pronounced for Yeast, WEBKB, and Facebook. For the latter two datasets, this may be due to the fact that we exclude a large number of longer patterns as we lack ground truth for
<table><caption>Table 7: The average relative error for each pattern size for FK-OBD, FK-AD, FACT and Random on the Yeast, WEBKB, IMDB, and Facebook datasets. FK-OBD, FK-AD, and Random's were permitted the same amount of run time as FACT took for each pattern.</caption><thead><tr><th rowspan="2">Dataset</th><th rowspan="2">Method</th><th colspan="10">Pattern Size</th><th rowspan="2">Average</th></tr><tr><th>4</th><th>5</th><th>6</th><th>7</th><th>8</th><th>9</th><th>10</th></tr></thead><tbody><tr><td rowspan="4">Yeast</td><td>FK-OBD</td><td>0.023</td><td>0.029</td><td>0.046</td><td>0.029</td><td>0.058</td><td>0.067</td><td>0.083</td><td>0.048</td></tr><tr><td>FK-AD</td><td>0.026</td><td>0.032</td><td>0.047</td><td>0.037</td><td>0.072</td><td>0.081</td><td>0.082</td><td>0.054</td></tr><tr><td>FACT</td><td>0.887</td><td>0.942</td><td>0.949</td><td>0.956</td><td>0.956</td><td>0.858</td><td>0.873</td><td>0.917</td></tr><tr><td>Random</td><td>0.069</td><td>0.153</td><td>0.232</td><td>0.171</td><td>0.194</td><td>0.210</td><td>0.263</td><td>0.184</td></tr><tr><td rowspan="4">WEBKB</td><td>FK-OBD</td><td>0.012</td><td>0.128</td><td>0.257</td><td>0.378</td><td>0.512</td><td>0.621</td><td>0.637</td><td>0.456</td></tr><tr><td>FK-AD</td><td>0.030</td><td>0.145</td><td>0.290</td><td>0.376</td><td>0.430</td><td>0.569</td><td>0.637</td><td>0.432</td></tr><tr><td>FACT</td><td>0.953</td><td>0.849</td><td>0.824</td><td>0.890</td><td>0.873</td><td>0.873</td><td>0.873</td><td>0.872</td></tr><tr><td>Random</td><td>0.036</td><td>0.281</td><td>0.503</td><td>0.576</td><td>0.658</td><td>0.704</td><td>0.630</td><td>0.592</td></tr><tr><td rowspan="4">IMDB</td><td>FK-OBD</td><td>0.006</td><td>0.011</td><td>0.034</td><td>0.089</td><td>0.218</td><td>0.206</td><td>0.188</td><td>0.124</td></tr><tr><td>FK-AD</td><td>0.007</td><td>0.014</td><td>0.023</td><td>0.063</td><td>0.098</td><td>0.128</td><td>0.135</td><td>0.077</td></tr><tr><td>FACT</td><td>0.825</td><td>0.814</td><td>0.840</td><td>0.831</td><td>0.927</td><td>0.913</td><td>0.838</td><td>0.861</td></tr><tr><td>Random</td><td>0.018</td><td>0.056</td><td>0.118</td><td>0.275</td><td>0.414</td><td>0.486</td><td>0.458</td><td>0.298</td></tr><tr><td rowspan="4">Facebook</td><td>FK-OBD</td><td>0.007</td><td>0.209</td><td>0.405</td><td>0.658</td><td>0.792</td><td>0.887</td><td>0.933</td><td>0.606</td></tr><tr><td>FK-AD</td><td>0.008</td><td>0.208</td><td>0.406</td><td>0.662</td><td>0.794</td><td>0.890</td><td>0.933</td><td>0.607</td></tr><tr><td>FACT</td><td>0.996</td><td>0.999</td><td>0.999</td><td>0.999</td><td>0.999</td><td>0.999</td><td>0.999</td><td>0.999</td></tr><tr><td>Random</td><td>0.014</td><td>0.247</td><td>0.438</td><td>0.673</td><td>0.789</td><td>0.848</td><td>0.901</td><td>0.605</td></tr></tbody></table>

them (i.e., because the Exhaustive approach did not finish in its allotted time). Based on this results, it seems that FK-OBD and FK-AD can be used to learn accurate parameters for a LBN.

<table><caption>Table 8: The average KLD over all pattern sizes and its standard deviation for the FK-OBD, FK-AD and Random approaches on the Yeast, Facebook, WEBKB, IMDB and DBLP dataset. The best results are shown in bold.</caption><thead><tr><th>Dataset</th><th>FK-OBD</th><th>FK-AD</th><th>Random</th></tr></thead><tbody><tr><td>Yeast</td><td>0.053 ± 0.170</td><td><strong>0.049 ± 0.165</strong></td><td>0.283 ± 0.277</td></tr><tr><td>Facebook</td><td><strong>0.237 ± 0.256</strong></td><td>0.243 ± 0.252</td><td>0.331 ± 0.322</td></tr><tr><td>WEBKB</td><td><strong>0.088 ± 0.176</strong></td><td>0.090 ± 0.189</td><td>0.163 ± 0.241</td></tr><tr><td>IMDB</td><td>0.006 ± 0.027</td><td>0.004 ± 0.012</td><td><strong>0.002 ± 0.006</strong></td></tr><tr><td>DBLP</td><td>0.003 ± 0.023</td><td><strong>0.001 ± 0.011</strong></td><td>0.006 ± 0.031</td></tr></tbody></table>

Finally, the Exhaustive approach exceeded the given time limit for all patterns on Amazon and Enron, for a significant number of patterns on WEBKB and Facebook, and some patterns on DBLP. As we lack ground truth for these datasets, we explore the convergence behavior of each sampling approach on these datasets in the following way. We treat each approach's final estimate as the ground truth. Then we show how the approach converges to this estimate as a function of run time. Here, we focus on the average relative error and consider only the patterns were the Exhaustive approach failed to finish in the time limit. Figure 7 and Figure 8 show the results for this experiment. For Amazon, Enron and Facebook,
<table><caption>Table 9: The average KLD for each pattern size for the FK-ODB, FK-AD, and Random approaches on the Yeast, WEBKB, IMDB, Facebook and DBLP datasets. The entries with all zeros represent numbers that are smaller than $10^{-4}$.</caption><thead><tr><th rowspan="2">Dataset</th><th rowspan="2">Method</th><th colspan="10">Pattern Size</th><th rowspan="2">Average</th></tr><tr><th>4</th><th>5</th><th>6</th><th>7</th><th>8</th><th>9</th><th>10</th></tr></thead><tbody><tr><td rowspan="3">Yeast</td><td>FK-OBD</td><td>0.011</td><td>0.006</td><td>0.003</td><td>0.004</td><td>0.026</td><td>0.176</td><td>0.181</td><td>0.053</td></tr><tr><td>FK-AD</td><td>0.011</td><td>0.006</td><td>0.003</td><td>0.004</td><td>0.019</td><td>0.169</td><td>0.176</td><td>0.049</td></tr><tr><td>Random</td><td>0.033</td><td>0.182</td><td>0.155</td><td>0.346</td><td>0.485</td><td>0.466</td><td>0.559</td><td>0.283</td></tr><tr><td rowspan="3">WEBKB</td><td>FK-OBD</td><td>0.001</td><td>0.028</td><td>0.042</td><td>0.139</td><td>0.105</td><td>0.115</td><td>0.159</td><td>0.088</td></tr><tr><td>FK-AD</td><td>0.001</td><td>0.027</td><td>0.044</td><td>0.070</td><td>0.093</td><td>0.123</td><td>0.167</td><td>0.090</td></tr><tr><td>Random</td><td>0.001</td><td>0.024</td><td>0.088</td><td>0.158</td><td>0.233</td><td>0.306</td><td>0.311</td><td>0.163</td></tr><tr><td rowspan="3">IMDB</td><td>FK-OBD</td><td>0.001</td><td>0.001</td><td>0.001</td><td>0.004</td><td>0.009</td><td>0.013</td><td>0.006</td><td>0.006</td></tr><tr><td>FK-AD</td><td>0.001</td><td>0.001</td><td>0.001</td><td>0.004</td><td>0.009</td><td>0.006</td><td>0.005</td><td>0.004</td></tr><tr><td>Random</td><td>0.000</td><td>0.000</td><td>0.000</td><td>0.002</td><td>0.004</td><td>0.004</td><td>0.003</td><td>0.002</td></tr><tr><td rowspan="3">Facebook</td><td>FK-OBD</td><td>0.003</td><td>0.015</td><td>0.106</td><td>0.235</td><td>0.343</td><td>0.508</td><td>0.563</td><td>0.237</td></tr><tr><td>FK-AD</td><td>0.003</td><td>0.015</td><td>0.118</td><td>0.241</td><td>0.366</td><td>0.448</td><td>0.553</td><td>0.243</td></tr><tr><td>Random</td><td>0.005</td><td>0.087</td><td>0.312</td><td>0.549</td><td>0.622</td><td>0.656</td><td>0.780</td><td>0.331</td></tr><tr><td rowspan="3">DBLP</td><td>FK-OBD</td><td>0.000</td><td>0.000</td><td>0.000</td><td>0.000</td><td>0.001</td><td>0.003</td><td>0.008</td><td>0.003</td></tr><tr><td>FK-AD</td><td>0.000</td><td>0.000</td><td>0.000</td><td>0.000</td><td>0.000</td><td>0.002</td><td>0.004</td><td>0.001</td></tr><tr><td>Random</td><td>0.000</td><td>0.000</td><td>0.001</td><td>0.000</td><td>0.001</td><td>0.005</td><td>0.022</td><td>0.006</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>Fig. 7: The average relative error as a function of the run time (in minutes) for the FK-OBD, FK-AD, and Random approach on the Amazon, Enron, and DBLP datasets.</figcaption></figure>

FK-ODB and FK-AD display nearly identical convergence behavior. On the DBLP and WEBKB datasets, FK-AD converges faster than FK-ODB. They typically converge relatively quickly. Random also converges quickly, but typically with worse estimates than FK-AD and FK-OBD.

## 9.4 Experiments and Results for Q3

The performance of FK-ODB and FK-AD are very similar. One explanation is that FK-AD uses the OBD of the pattern to construct the AD. To evaluate if this is the case, we construct random arbitrary decompositions in the following way:
<figure><img src="image_9.png" /><figcaption>Fig. 8: The average relative error as a function of the run time (in seconds) for the FK-OBD, FK-AD, and Random approach on the Facebook and WEBKB datasets.</figcaption></figure>

<table><caption>Table 10: The average relative error and its standard deviation for the FK-AD and FK-AD-10 on the Yeast, Facebook, WEBKB, and IMDB datasets. The average is over all pattern sizes. The best results are shown in bold.</caption><thead><tr><th>Dataset</th><th>FK-AD</th><th>FK-AD-10</th></tr></thead><tbody><tr><td>Yeast</td><td>0.016 ± 0.027</td><td><strong>0.007 ± 0.009</strong></td></tr><tr><td>Facebook</td><td>0.048 ± 0.444</td><td><strong>0.046 ± 0.445</strong></td></tr><tr><td>WEBKB</td><td>0.032 ± 0.077</td><td><strong>0.007 ± 0.012</strong></td></tr><tr><td>IMDB</td><td>0.027 ± 0.074</td><td><strong>0.022 ± 0.072</strong></td></tr></tbody></table>

1. Randomly select a first node in the ordering.

2. Randomly select a node such that it is adjacent to at least one previous node in the existing partial order and add it to the order.

3. Repeat step 2 until all nodes are in the order.

For each pattern, we generate ten such decompositions and compute the counts using each one and report the average over these ten runs. We denote the results for this approach as FK-AD-10. Due to the computational demands of this experiment, we only consider the following four datasets: Yeast, Facebook, WEBKB and IMDB.

Table 10 shows the average relative error for FK-AD and FK-AD-10 averaged over all pattern sizes. The results for both approaches are quite similar, with FK-AD-10 performing slightly better. Figure 9 shows the results for each pattern size. These results show that, regardless of the pattern size, FK-OBD, FK-AD and FK-AD-10 tend to perform similarly. This gives some evidence that the Fürer-Kasiviswanathan approach is robust when a non-OBD decomposition is used to estimate the counts.

Figures 10, 11, 12, and 13 show the average relative error as a function of the percentage of the Exhaustive approach's run time for various pattern sizes on the Facebook, WEBKB, Yeast, and IMDB datasets. Again, these results only consider the patterns where the Exhaustive approach finished within the allotted time so
<figure><img src="image_10.png" /><figcaption>Fig. 9: The average relative error for the FK-OBD, FK-AD, FK-AD-10 and Random approaches on each pattern size for the Yeast, Facebook, WEBKB and IMDB datasets.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 10: The average relative error for a specific pattern size as a function of the percentage of the Exhaustive approach's run time for the FK-OBD, FK-AD, FK-AD-10 and Random approaches on the Facebook dataset.</figcaption></figure>

we know the true count. On all datasets and pattern sizes FK-AD-10 displays a similar convergence behavior to FK-OBD and FK-AD. Like in Table 10, we can see that FK-AD-10 yields similar or slightly better results than FK-OBD and FK-AD.
<figure><img src="image_12.png" /><figcaption>Fig. 11: The average relative error for a specific pattern size as a function of the percentage of the Exhaustive approach's run time for the FK-OBD, FK-AD, FK-AD-10 and Random approaches on the WEBKB dataset.</figcaption></figure>

## 9.5 Experiments and Results for Q4

The theoretical analysis of Fürer and Kasiviswanathan only applies to random graphs. To explore in a controlled setting how the algorithm performs when this condition does and does not hold, we use the synthetic datasets consisting of random and power law graphs. Figure 14 and Figure 15 show the performance of each algorithm for the KLD and average relative error, respectively. On these datasets, the Exhaustive approach always finished within the time limit. On average, FK-OBD has a better KLD and smaller relative errors than both FK-AD and Random. Figure 16 compares the average relative error for FK-OBD and FACT. FACT performs much worse than all the other methods on these datasets. These results provide some evidence that FK-OBD and FK-AD will still perform well even on non-random graphs where the theoretical analysis does not hold.

## 9.6 Experiments and Results for Q5

Next, we evaluate our approach in the context of weight learning and inference. We consider hand-crafted models for the WEBKB and Yeast datasets. In the WEBKB dataset, the goal is to predict label of a web page, which is a standard task. The model consists of four hand-crafted patterns. In the Yeast dataset, we predict the location of a protein. For this problem, the model consists of five hand-crafted patterns. In both models, the patterns were derived from first-order logical formulas where the target predicate (PageClass and Location) appeared exactly once. We
<figure><img src="image_13.png" /><figcaption>Fig. 12: The average relative error for a specific pattern size as a function of the percentage of the Exhaustive approach's run time for the FK-OBD, FK-AD, FK-AD-10 and Random approaches on the Yeast dataset.</figcaption></figure>

train a logistic regression model whose features are the counts associated with each pattern to predict the truth value of each grounding of the target predicate. Our models can be viewed as being analogous to an MLN because a discriminatively trained MLN model that only contains non-recursive clauses (which corresponds to our setting) is equivalent to a logistic regression model (Huynh and Mooney 2008).

In the WEBKB dataset, we restrict ourselves to considering 50 words (Davis and Domingos 2009). Both the WEBKB and Yeast datasets are divided into four folds, and hence we perform four-fold cross validation. Both weight learning (i.e., training the model) and inference (i.e., making the predictions on the test set) require computing counts. We compare the predictive and run time performance of using the Exhaustive approach to obtain the exact counts for both learning and inference versus using FK-OBD to obtain the approximate counts for both learning and inference. For FK-OBD, we consider the count obtained for a pattern using four different time limits: 1 second (FK-OBD-1s), 10 seconds (FK-ODB-10s), 1 minute (FK-OBD-1m), and 5 minutes (FK-OBD-5m).
<figure><img src="image_14.png" /><figcaption>Fig. 13: The average relative error for a specific pattern size as a function of the percentage of the Exhaustive approach's run time for the FK-OBD, FK-AD, FK-AD-10 and Random approaches on the IMDB dataset.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>Fig. 14: The KLD averaged over all sampling iterations on the synthetic datasets. The datasets are sorted by increasing graph density. The error bars represent the standard deviation.</figcaption></figure>

<caption>Fig. 15: The relative error averaged over all iterations of the sampling algorithms on the synthetic datasets. The datasets are sorted by increasing graph density. The error bars represent the standard deviation.</caption>
Table 11 shows the area under the ROC curve (AUC-ROC) for all approaches together with the speedup factor of the FK-ODB approaches over the Exhaustive approach for the counting time needed for learning and inference. On both datasets, the sweet spot seems to be running FK-ODB for 1 minute. There is essentially no loss of accuracy, and it achieves a 100 time speedup on WEBKB and a 10 time speedup on Yeast. Running FK-OBD for a longer time period only marginal improves predictive performance. Similarly, running FK-OBD for 10 second results in a small decline in predictive performance on WEBKB and a larger one on Yeast, albeit with a large saving in run time cost. Unsurprisingly, on these datasets, running FK-ODB for 1 second degrades the predictive performance substantially.

<table><thead><tr><th>Dataset</th><th>Metrics</th><th>Exhaustive</th><th>FK-ODB-1s</th><th>FK-ODB-10s</th><th>FK-ODB-1m</th><th>FK-ODB-5m</th></tr></thead><tbody><tr><td rowspan="3">WEBKB</td><td>Speedup Factor: Learning</td><td>N/A</td><td>12916.7</td><td>1291.7</td><td>215.3</td><td>43.1</td></tr><tr><td>Speedup Factor: Inference</td><td>N/A</td><td>9642.6</td><td>964.3</td><td>160.7</td><td>32.1</td></tr><tr><td>AUC-ROC</td><td>0.99± 0.014</td><td>0.74±0.033</td><td>0.93±0.028</td><td>0.97±0.015</td><td>0.98±0.013</td></tr><tr><td rowspan="3">Yeast</td><td>Speedup Factor: Learning</td><td>N/A</td><td>585.0</td><td>58.5</td><td>9.8</td><td>2.2</td></tr><tr><td>Speedup Factor: Inference</td><td>N/A</td><td>585.7</td><td>58.6</td><td>9.8</td><td>2.1</td></tr><tr><td>AUC-ROC</td><td>1.00±0.000</td><td>0.64±0.005</td><td>0.83±0.010</td><td>0.97±0.006</td><td>0.99±0.005</td></tr></tbody></table>

Table 11: The average AUC-ROC and its standard deviation when the counts required for learning and inference were obtained by the Exhaustive approach and FK-ODB with a time limit of 1 second (FK-OBD-1s), 10 seconds (FK-ODB-10s), 1 minute (FK-OBD-1m), and 5 minutes (FK-OBD-5m). Additionally, the speedup factor of the FK-ODB approaches over the Exhaustive approach for the counting time needed during both learning and inference is shown. The speedup factor for the Exhaustive approach over itself is not applicable (N/A).

## 9.7 Experiments and Results for Q6

Finally, we explore how the approach will handle longer patterns. Employing the same pattern generation and selection approach outlined in Subsection 9.2.1, we generate patterns containing from 11 up to 15 nodes for the DBLP and Yeast datasets. Here, we only focus on the FK-ODB and FK-AD approaches.

Figures 17 and 18 present results for this setting. Each plot shows the average relative error of FK-ODB and FK-AD when those approaches are given a fixed percentage of the run time required for the Exhaustive approach to complete. On DBLP, regardless of the pattern size, both approaches arrive at a relatively stable estimate using only a small fraction of the time needed by the Exhaustive approach. Given more time, the estimates improve slightly. On Yeast, longer pattern sizes require slightly more time to converge. On this dataset, there is not much improvement when running beyond 20 to 40 percent of the Exhaustive approach's run time. On both datasets, FK-OBD and FK-AD perform similarly.
<figure><img src="image_16.png" /><figcaption>Fig. 16: The average relative error for FK-OBD and FACT on the synthetic datasets. The datasets are sorted by increasing graph density.</figcaption></figure>

## 9.8 Discussion

Now we will revisit the experimental questions posed at the beginning of this section in light of the presented results. In terms of **Q1**, across the various datasets and pattern sizes considered, we see that the FK-OBD and FK-AD approaches result in much more accurate counts than Random and FACT when each algorithm has the same run time. Furthermore, both FK-ODB and FK-AD converge rather quickly to counts that are close to the true count found by Exhaustive search. Finally, we see that FK-OBD and FK-AD are able to learn highly accurate CPDs for LBNs on almost all datasets and pattern sizes.

<figure><img src="image_17.png" /><figcaption>Fig. 17: The average relative error for a specific pattern size as a function of the percentage of the Exhaustive approach's run time for the FK-OBD, FK-AD, and Random approaches on the DBLP dataset. The focus is longer pattern sizes.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>Fig. 18: The average relative error for a specific pattern size as a function of the percentage of the Exhaustive approach's run time for the FK-OBD, FK-AD, and Random approaches on the YEAST dataset. The focus is longer pattern sizes.</figcaption></figure>

## 9.9 Discussion

Finally, we revisit the experimental questions posed at the beginning of this section in light of the presented results. In terms of **Q1**, across the various datasets and pattern sizes considered, we see that the FK-OBD and FK-AD approaches result in much more accurate counts than Random and FACT when each algorithm has the same run time. Furthermore, both FK-ODB and FK-AD converge rather quickly to counts that are close to the true count found by Exhaustive search.

For **Q2**, we can see that FK-OBD and FK-AD tend to perform similarly. Thus even when the decomposition is not an OBD the approach seems to converge in practice. One possible reason is that because the size of the AD is larger, it is faster to match. Hence, it makes more observations than FK-OBD which needs to calculate all possible combinations of values for nodes in each partition of an OBD. The larger the partition, the longer it takes to finish one iteration of sampling. This result holds for both approaches we used to construct an arbitrary decomposition (**Q3**).

The synthetic data allow us to address **Q4**. Here, we find that FK-OBD and FK-AD perform well on non-random graphs where the theoretical guarantees of the original Fürer-Kasiviswanathan algorithm are not valid. These results, in conjunction with our findings on the real-world data, indicate that FK is a viable strategy for non-random graphs.

To address **Q5**, we considered using FK-OBD to obtain the counts needed for learning and inference. We see that FK-ODB achieves significant gains in run time performance for both tasks. These gains, which are between one and two orders of magnitude, come only with a very small drop in predictive performance.
By looking at results on patterns up to length 15, we see that the behavior of FK-ODB and FK-AD are consistent with how they perform on shorter patterns. These results thus provide some evidence for a positive answer to **Q6**.

Finally, we hope that this work can benefit and be exploited by the broader community. Therefore, we have released our code and the patterns used.<fnref n="9" />

# 10 Conclusions

This paper explored the problem of estimating the number of embeddings of a pattern in a graph, which is an important task in many areas like (approximate) pattern mining and parameter estimation for statistical relational models. Specifically, the paper presented several different sampling algorithms, studied how to obtain statistics from the embeddings they generated, and performed an extensive empirical comparison.

Empirically, we observed that a strategy based on the theoretical algorithm by Fürer and Kasiviswanathan (2014) performed well, even in conditions where the requirements for its theoretical guarantees were violated. Namely, we found that the algorithm still exhibited good performance when given (1) a decomposition that is not an OBD, and (2) a data graph that is not an Erdős-Rényi random graph. We found that the FK-based approaches offered several benefits. One, they outperformed several competitors. Two, they converged to a good estimate in the early stages of sampling. Three, they were able to produce accurate estimates for longer patterns.

<footnote n="9">⁹ https://dtai.cs.kuleuven.be/software/gs-srl</footnote>
Acknowledgments

IR was supported by the KU Leuven Research Fund (OT/11/051). MZ was partially supported by the KU Leuven Research Fund (OT/11/051) and the Slovenian Research Agency (P2-0103). JD is partially supported by the KU Leuven Research Fund (OT/11/051, C14/17/070, C22/15/015, C32/17/036) and FWO-Vlaanderen (G.0356.12, SBO-150033).

References

Ariely D (2008) Predictably irrational: The hidden forces that shape our decisions. Harper Collins

Barabasi AL, Albert R (1999) Emergence of scaling in random networks. Science 286(5439):509-512

Baskerville K, Grassberger P, Paczuski M (2007) Graph animals, subgraph sampling, and motif search in large networks. Physical Review E 76(3):036,107

Bordino I, Donato D, Gionis A, Leonardi S (2008) Mining large networks with subgraph counting. In: Proceedings of the 2008 IEEE International Conference on Data Mining (ICDM), pp 737-742

Cordella L, Foggia P, Sansone C, Vento M (2004) A (sub)graph isomorphism algorithm for matching large graphs. IEEE Transactions on Pattern Analysis and Machine Intelligence 26(10):1367-1372

Das M, Wu Y, Khot T, Kersting K, Natarajan S (2016) Scaling lifted probabilistic inference and learning via graph databases. In: Proceedings of the 2016 SIAM International Conference on Data Mining (SDM), pp 738-746

Davis J, Domingos P (2009) Deep transfer via second-order Markov logic. In: Proceedings of the Twenty-Sixth International Conference on Machine Learning (ICML), pp 217-224

Davis J, Burnside E, Dutra IC, Page D, Costa VS (2005) An integrated approach to learning Bayesian networks of rules. In: Proceedings of the Sixteenth European Conference on Machine Learning (ECML), pp 84-95

Di Natale R, Ferro A, Giugno R, Mongiovi M, Pulvirenti A, Shasha D (2010) SING: Subgraph search in non-homogeneous graphs. BMC Bioinformatics 11(1):96

Fierens D, Blockeel H, Ramon J, Bruynooghe M (2004) Logical Bayesian networks. In: Proceedings of the Third International Workshop on Multi-Relational Data Mining (MRDM), pp 19-30

Friedman N, Goldzsmidt M (1996) Learning Bayesian networks with local structure. In: Proceedings of the Twelfth Annual Conference on Uncertainty in Artificial Intelligence (UAI), pp 252-262

Fürer M, Kasiviswanathan SP (2014) Approximately counting embeddings into random graphs. Combinatorics, Probability & Computing 23(6):1028-1056

Getoor L, Taskar B (2007) Introduction to statistical relational learning. MIT press

Giugno R, Shasha D (2002) GraphGrep: A fast and universal method for querying graphs. In: Proceedings of the Sixteenth International Conference on Pattern Recognition (ICPR), pp 112-115

Huynh T, Mooney R (2008) Discriminative structure and parameter learning for Markov logic networks. In: Proceedings of the Twenty-Fifth International Conference on Machine Learning, pp 416-423

Inokuchi A, Washio T, Motoda H (2003) Complete mining of frequent patterns from graphs: Mining graph data. Machine Learning 50(3):321-354

Jowhari H, Ghodsi M (2005) New streaming algorithms for counting triangles in graphs. In: Proceedings of the Eleventh International Conference on Computing and Combinatorics (COCOON), pp 710-716

Kashtan N, Itzkovitz S, Milo R, Alon U (2004) Efficient sampling algorithm for estimating subgraph concentrations and detecting network motifs. Bioinformatics 20(11):1746-1758

Kersting K, De Raedt L, Kramer S (2000) Interpreting Bayesian logic programs. In: Proceedings of the AAAI-2000 Workshop on Learning Statistical Models from Relational Data, pp 29-35
Kok S, Domingos P (2005) Learning the structure of Markov logic networks. In: Proceedings of the Twenty-Second International Conference on Machine Learning (ICML), pp 441–448

Leskovec J, Faloutsos C (2006) Sampling from large graphs. In: Proceedings of the Twelfth ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (KDD), pp 631–636

Mewes HW, Frishman D, Gruber C, Geier B, Haase D, Kaps A, Lemcke K, Mannhaupt G, Pfeiffer F, Schüller C, Stocker S, Weil B (2000) MIPS: A database for genomes and protein sequences. Nucleic Acids Research 28(1):37–40

Pržulj N (2007) Biological network comparison using graphlet degree distribution. Bioinformatics 23(2):177–183

Ravkic I, Ramon J, Davis J (2015) Learning relational dependency networks in hybrid domains. Machine Learning 100(2-3):217–254

Richards BL, Mooney RJ (1992) Learning relations by pathfinding. In: Proceedings of the Tenth National Conference on Artificial Intelligence (AAAI), pp 50–55

Richardson M, Domingos P (2006) Markov logic networks. Machine learning 62(1-2):107–136

Shervashidze N, Vishwanathan S, Petri T, Mehlhorn K, Borgwardt K (2009) Efficient graphlet kernels for large graph comparison. In: Proceedings of the Twelfth International Conference on Artificial Intelligence and Statistics (AISTATS), pp 488–495

Tang J, Zhang J, Yao L, Li J, Zhang L, Su Z (2008) Arnetminer: Extraction and mining of academic social networks. In: Proceedings of the Fourteenth ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (KDD), pp 990–998

Ullmann JR (1976) An algorithm for subgraph isomorphism. Journal of the ACM 23(1):31–42

Van Haaren J, Kolobov A, Davis J (2015) TODTLER: Two-order-deep transfer learning. In: Proceedings of the Twenty-Ninth AAAI Conference on Artificial Intelligence, pp 3007–3015

Venugopal D, Sarkhel S, Gogate V (2015) Just count the satisfied groundings: Scalable local-search and sampling based inference in MLNs. In: Proceedings of the Twenty-Ninth AAAI Conference on Artificial Intelligence, pp 3606–3612

Wernicke S (2005) A faster algorithm for detecting network motifs. In: Proceedings of the Fifth International Workshop on Algorithms in Bioinformatics (WABI), pp 165–177

Yan X, Han J (2002) gSpan: Graph-based substructure pattern mining. In: Proceedings of the 2002 IEEE International Conference on Data Mining (ICDM), pp 721–724

Zou R, Holder LB (2010) Frequent subgraph mining on a single large graph using sampling techniques. In: Proceedings of the Eighth Workshop on Mining and Learning with Graphs (MLG), pp 171–178