Markov Chain Monte Carlo (MCMC) is one of the most used families of algorithms based on sampling. They allow to sample from the posterior distribution when direct sampling from it is infeasible, due to the complexity of the distribution itself. Gibbs sampling is one of these algorithm that has been applied in many situations. In this paper we compare an implementation of Gibbs sampling for Probabilistic Logic Programs on several datasets, in order to better understand its performance. For all the experiments we compute the convergence time, execution time and population standard deviation of the samples.
Probabilistic Logic Programming has been proved e ective in modelling several
real world situations [
Inference in probabilistic logic programs is a computationally hard task. There are two main types of inference. The rst type is exact inference, which is usually based on the representation of the program in a compact form such as decision diagrams, and it is aimed to compute exact answers. However, this type of inference is not always easy to apply due to some limitations, among them, the time required to get an answer. Approximate inference based on sampling may overcome the exact inference limitations but with a cost: the accuracy of the results depends on the number of samples. Moreover, di erent algorithms may have di erent performance on the same dataset.
Among the approximate inference algorithms, Markov Chain Monte Carlo (MCMC) allows to e ciently sample from the posterior distribution. In this ? Copyright c 2020 for this paper by its authors. Use permitted under Creative Commons License Attribution 4.0 International (CC BY 4.0). paper, we focus on one of them, Gibbs sampling, and we test it on several datasets in order to get a quantitative analysis of its performance.
The paper is structured as follows: in Section 2 we introduce the basic concepts of Probabilistic Logic Programming and in Section 3 we analyze approximate inference and Gibbs sampling. Section 4 illustrates the results of the experiments and Section 5 concludes the paper. 2
Here we consider Probabilistic Logic Programs (PLPs) based on the Distribution
Semantics [
An LPAD clause has the form h1 : 1; : : : ; hn : n : b1; : : : ; bm where the head is a set of logical atoms and the body is a set of logical literals. is are real values in the interval [0; 1] that represent the probability that the i-th head is chosen. Moreover, Pi i = 1. If this is not true, i.e. 0 < Pi i < 1 an extra null atom is added to the head with probability 1 Pi i.
A PLP program de nes a probability distribution over worlds. Recall that a substitution is a function mapping variables to terms represented with = fX1=t1; : : : ; Xn=tng, where Xk=tk means that all the occurrences of the variable Xk are replaced by the term tk. A world can be obtained by choosing one atom from the head of a grounding of a LPAD clause. The probability of a query is then computed by considering a joint probability distribution between the query and the worlds and by summing out the worlds.
The previous de nition can only be applied for programs without functions
symbols, since their grounding is nite. In case of function symbols, the de nition
must be extended. See [
Let us introduce an example of LPAD: mistake(X) : 0:1 : good player(X): mistake(X) : 0:05 : f ocused(X): good player(kasparov): f ocused(kasparov): The previous program represents a situation where X makes a mistake with probability 0.1 if X is a good player and nothing happens with probability 0.9. Similarly, if X is f ocused during the match, he will make a mistake with probability 0.05 and nothing happens with probability 0.95. Moreover, we know for sure that kasparov is a good player and he is f ocused during a match. We may be interested, for instance, in the probability that kasparov makes a mistake, denoted with P (mistake(kasparov)). This can be computed as: 0:1 0:05 + 0:1 (1 0:05) + (1 0:1) 0:05 = 0:145.
The probability of the previous example can be computed exactly since it has a nite grounding. Consider now the following example, representing a player that repeatedly toss a coin.
heads(N ) : 0:5; tails(N ) : 0:5: on(0; h) : heads(0): on(0; t) : tails(0): on(s(N ); h) : on(N; t); heads(s(N )): on(s(N ); t) : on(N; t); tails(s(N )): at least once head : on( ; h): If we want to know the probability that the player will get at least once heads (P (at least once head)), exact inference cannot be used since the grounding of the program is in nite. However, using approximate methods based on sampling, we can sample the query a certain number of times and return the number of successes over the number of samples. As the number of samples goes to in nity, the approximated probability reaches value 0.5 (since we need to compute the sum of a geometric series).
In this paper, we utilize the cplint framework [
The computation of the probability of a query, a task called inference, can
be performed both with exact and approximate methods. Exact inference from
discrete programs (programs without continuous random variables) is a
#Pcomplete task [
Approximate inference methods may be implemented as a variation of exact inference algorithms or can be based on sampling. Here we focus on the latter. One of the most used approximate approach is given by Monte Carlo sampling: in a nutshell, the algorithm samples a world by sampling every ground probabilistic fact, checks if the query is true in the world and computes the probability of the query as the fraction of samples where the query is true. This process is repeated for a xed number of steps or until convergence (i.e., the di erence between two consecutive samples is less than a certain threshold).
The cplint framework [
The simplest algorithm used to perform this task is rejection sampling: it queries the evidence in a sample program and, if the query is successful, queries the goal. The probability is then computed and the fraction of successes over the number of samples. Despite its simplicity, rejection sampling may su er from the fact that, if the evidence has low probability, a lot of samples are discarded making the algorithm very slow to converge.
An alternative is represented by the Markov Chain Monte Carlo (MCMC) methods. The main idea behind MCMC is to construct a Markov Chain that has as equilibrium distribution the desired distribution. Then, samples can be drawn from this Markov Chain in order to get the probability of a query. As the number of steps increase, the approximation gets closer to the real distribution. In the limit, MCMC can represent the real posterior distribution. A common operation is to discard the rst samples since they do not represent the real distribution, with an operation called burn-in. There are several algorithms in this family but here we focus on Gibbs sampling.
Gibbs sampling [
To deeply analyze the performance of Gibbs sampling for PLPs, we tested the
algorithm on eight di erent datasets. For each experiment, we plotted three
graphs: one to track the performance in terms of number of samples required to
convergence, one to track the execution time and one to track the population
standard deviation of the samples. For all the experiments the code is available
online at the indicated URLs. All the experiments were conducted on a cluster3
with Intel R Xeon R E5-2630v3 running at 2.40 GHz. The results are averages of
10 runs. Execution times are computed using the SWI-Prolog built-in predicate
statistics/2 with the keyword walltime. The probability of the evidence is
computed using the cplint predicate mc sample/3 with 106 samples. The
predicate samples the query a certain number of times and returns the probability
of success. The number of discarded samples can be set with the option mix/1
while the block size with block/1. Here we set mix to 100 and block variable
between 1 and 5. The experiments are described below:
{ Arithm4: the program represents a random arithmetic function. The goal
is to compute the conditional probability of the value of the function given
that a couple of input output was observed. A characteristic of this program
is that it has an in nite number of explanations. The probability of evidence
is 0.05. Fig. 1 shows that with a relatively small number of samples, the
probability computed with Gibbs sampling with block set to 1 oscillates
between 0.1 and 0.25, as also con rmed by the large variation in the standard
deviation (Fig. 9 Left). For the other values, the oscillation is smaller but
still present.
{ Diabetes5: the program models the probability of insurgence of diabetes
given that some genetic factors are observed. This is an example of
probabilistic constraint logic program [
{ Graph6: in the following experiment we want to test the accuracy of Gibbs
sampling on a Barabasi Albert preferential attachment model. Given a graph
with initially connected nodes, new nodes are added to the graph. The
probability that these new nodes are connected to other nodes is
proportional to the number of edges that the already existent nodes have. We
generated the graph with the python library networkx7 using the function
barabasi albert graph(40,10). For all the edges we set a probability of
0.1. We compute the probability that two nodes are connected given that a
portion of the path has already been observed (probability of the evidence
0.42). The probability of this query can also be computed using exact
inference algorithms. However, when the size of the graph increases, exact
inference may be too expensive in terms of execution time. Fig. 3 shows
that, as for the previous experiments, smaller values of blocks leads to very
di erent probability values. In particular, with block set to 1, there is a gap
of 0.2 between some values. Execution times for all the ve block settings
are equivalent. Standard deviation of the samples (Fig. 10 Left) decreases
as the block number increases. With block set to 3, 4 and 5, the computed
values are almost the same.
{ Hidden Markov Model8 (HMM). The program represents a model of DNA
sequences using an HMM [
known distributions. We want to know the probability that the rst drawn
ball is made of wood given that its color is black. The probability of the
evidence is 0.38. Fig. 6 shows that the probability completely stabilizes after
104 samples only for large values of block. Block set to 1 and 3 seems to
be the values with more instability. However, the standard deviation of the
samples decreases for all the block values, as the number of samples increase
(Fig.11 Right). Execution times are comparable.
{ Pre x parser11. This program models a pre x parser for probabilistic context
free grammars [
0:8 0:6
1 104 0:6
0:8 0:2 0:4 0:8 0:6
1 104 0:6
1 104 Fig. 7. Results for pre x parser experiment. 105 104 6 ) s (m4 e m i T n o i tu 2 c e x E 0 4 ) s ( 3 m e m i T n 2 o it u c e x 1 E 0 2:00 ) s (m1:50 e m i T n 1:00 o it u c e x 0:50 E 0:00 0:2 0:4 0:6
0:8
1 104 0:6 Samples 1 104 0:2 0:4 0:8 0:6
1 104 0:6 Samples 1 104 0:2 0:4 0:6 0:8
1 104 Fig. 10. Standard Deviation for graph (left) and HMM (right) experiments. 106
0:6 Samples
In this paper we tested Gibbs sampling for probabilistic logic programs on eight
di erent datasets. For each dataset we tracked execution time, computed
probability and standard deviation of the samples. We used the cplint framework
and the code for all the experiments can be analyzed trough a web interface
accessible at cplint:eu. Empirical results shows that, when the value of block
increases, the probability and the standard deviation seems to stabilize with a
smaller number of samples. However, the execution time increases as well with
an increment of the block value. As a future work, we plan to apply approximate
reasoning also to abduction [