2.1 Answer set programming

In this paper, we consider ASP (Brewka et al. Reference Brewka, Eiter and Trusczynski2011). An answer set program contains disjunctive rules of the form $ h_1 ; \dots ; h_m \,{:\!-}\, \, b_1, \dots, b_n.$ , where the disjunction of atoms $h_i$ s is called head and the conjunction of literals $b_i$ s is called body. Rules with no atoms in the head are called constraints, while rules with no literals in the body and a single atom in the head are called facts. Choice rules are a particular type of rules where a single head is enclosed in curly braces as in $\{a\} \,{:\!-}\, b_1, \dots, b_n$ , whose meaning is that $a$ may or may not hold if the body is true. They are syntactic sugar for $a;not\_a \,{:\!-}\, b_1,\dots, b_n$ where $not\_a$ is an atom for a new predicate not appearing elsewhere in the program. We allow aggregates (Alviano and Faber Reference Alviano and Faber2018) in the body, one of the key features of ASP that allows representing expressive relations between the objects of the domain of interest. An example of a rule $r0$ containing an aggregate is $v(A)\,{:\!-}\, \#count\{X:b(X)\}=A$ . Here, variable $A$ is unified with the integer resulting from the evaluation of the aggregate $\#count\{X:b(X)\}$ , which counts the elements $X$ such that $b(X)$ is true.

The semantics of ASP is based on the concept of a stable model (Gelfond and Lifschitz Reference Gelfond and Lifschitz1988). Every answer set program has zero or more stable models, also called answer sets. An interpretation $I$ for a program $P$ is a subset of its Herbrand base. The reduct of a program $P$ with respect to an interpretation $I$ is the set of rules of the grounding of $P$ with the body true in $I$ . An interpretation $I$ is a stable model or, equivalently, an answer set, of $P$ if it is a minimal model (under set inclusion) of the reduct of $P$ w.r.t. $I$ . We denote with $AS(P)$ the set of answer sets of an answer set program $P$ and with $|AS(P)|$ its cardinality. If $AS(P) = \emptyset$ we say that $P$ is unsatisfiable. If we consider rule $r0$ and two additional rules, $\{b(0)\}$ and $b(1)$ , the resulting program has 2 answer sets: $\{\{b(1), b(0), v(2)\},\{b(1), v(1)\}\}$ . Finally, we will also use the concept of projective solutions of an answer set program $P$ onto a set of atoms $B$ (defined as in (Gebser et al. Reference Gebser, Kaufmann, Schaub, Van Hoeve and Hooker2009)): $AS_{B}(P) = \{A \cap B \mid A \in AS(P)\}$ .

2.2 Probabilistic answer set programming

Uncertainty in Logic Programming can be represented with discrete Boolean probabilistic facts of the form $\Pi ::a$ where $\Pi \in [0,1]$ and $a$ is an atom that does not appear in the head of rules. These are considered independent: this assumption may seem restrictive, but, in practice, the same expressivity of Bayesian networks can be achieved by means of rules and extra atoms (Riguzzi Reference Riguzzi2022). Every probabilistic fact corresponds to a Boolean random variable. With probabilistic facts, a normal logic program becomes a probabilistic logic program. One of the most widely adopted semantics in this context is the Distribution Semantics (DS) (Sato Reference Sato1995). Under the DS, each of the $n$ probabilistic facts can be included or not in a world $w$ , generating $2^n$ worlds. Every world is a normal logic program. The DS requires that every world has a unique model. The probability of a world $w$ is defined as $ P(w) = \prod _{f_i \in w} \Pi _i \cdot \prod _{f_i \not \in w} (1 - \Pi _i).$

If we extend an answer set program with probabilistic facts, we obtain a probabilistic answer set program that we interpret under the CS (Cozman and Mauá Reference Cozman and Maua2016, Reference Cozman and Maua2020). In the following, when we write “probabilistic answer set program,” we assume that the program follows the CS. Similarly to the DS, the CS defines a probability distribution over the worlds. However, every world (which is an answer set program in this case) may have 0 or more stable models. A query $q$ is a conjunction of ground literals. The probability $P(q)$ of a query $q$ lies in the range $[\underline{P}(q), \overline{P}(q)]$ where

(1) \begin{equation} \begin{split} \underline{P}(q) &= \sum _{w_i \, \text{such that} \, \forall m \in AS(w_i), \, m \models q} P(w_i), \\[5pt] \overline{P}(q) &= \sum _{w_i \, \text{such that} \, \exists m \in AS(w_i), \, m \models q} P(w_i). \end{split} \end{equation}

That is, the lower probability is given by the sum of the probabilities of the worlds where the query is true in every answer, while the upper probability is given by the sum of the probabilities of the worlds where the query is true in at least one answer set. Here, as usual, we assume that all the probabilistic facts are independent and that they cannot appear as heads of rules. In other words, the lower (upper) probability is the sum of the probabilities of the worlds from which the query is a cautious (brave) consequence. Conditional inference means computing upper and lower bounds for the probability of a query $q$ given evidence $e$ , which is usually given as a conjunction of ground literals. The formulas are by Cozman and Mauá (Reference Cozman and Maua2017):

(2) \begin{equation} \begin{split} \underline{P}(q \mid e) &= \frac{\underline{P}(q,e)}{\underline{P}(q,e) + \overline{P}(\lnot q,e)} \\[5pt] \overline{P}(q \mid e) &= \frac{\overline{P}(q,e)}{\overline{P}(q,e) + \underline{P}(\lnot q,e)} \end{split} \end{equation}

Conditional probabilities are undefined if the denominator is 0.

Example 1. The following probabilistic answer set program has two probabilistic facts, $0.3::a$ and $0.4::b$ :

There are $2^2 = 4$ worlds: $w_0 = \{not\, a, not\, b\}$ , $w_1 = \{not\, a, b\}$ , $w_2 = \{a, not\, b\}$ , and $w_3 = \{a,b\}$ (with $not\, f$ we mean that the fact $f$ is absent, not selected), whose probability and answer sets are shown in Table 1. For example, the answer sets for $w_2$ are computed from the following program: $a.\, q0 ; q1 \,{:\!-}\, a.\, q0 \,{:\!-}\, b.$ If we consider the query $q0$ , it is true in all the answer sets of $w_1$ and $w_3$ and in one of those of $w_2$ ; thus, we get $[P(w_1) + P(w_3), P(w_1) + P(w_2) + P(w_3)] = [0.4,0.58]$ as probability bounds. $\square$

Table 1. Worlds, probabilities, and answer sets for Example1

DS and CS can be given an alternative but equivalent definition based on sampling in the following way. We repeatedly sample worlds by sampling every probabilistic fact obtaining a normal logic program or an answer set program. Then, we compute its model(s), and we verify whether the query is true in it (them). The probability of the query under DS is the fraction of worlds in whose model the query is true as the number of sampled worlds tends to infinity. For CS, the lower probability is the fraction of worlds where the query is true in all models, and the upper probability is the fraction of worlds where the query is true in at least one model. We call this the sampling semantics, and it is equivalent to the DS and CS definitions by the law of large numbers.

The complexity of the CS has been thoroughly studied by Mauá and Cozman (Reference Maua and Cozman2020). In particular, they focus on analyzing the cautious reasoning (CR) problem: given a PASP $P$ , a query $q$ , evidence $e$ , and a value $\gamma \in \mathbb{R}$ , the result is positive if $\underline{P}(q \mid e) \gt \gamma$ , negative otherwise (or if $\underline{P}(e) = 0$ ). A summary of their results is shown in Table 2. For instance, in PASP with stratified negation and disjunctions in the head, the complexity of the CR problem is in the class $PP^{\Sigma _2^p}$ , where $PP$ is the class of the problems that can be solved in polynomial time by a probabilistic Turing machine (Gill Reference Gill1977). The complexity is even higher if aggregates are allowed.

Table 2. Tight complexity bounds of the CR problem in PASP from Mauá and Cozman (Reference Maua and Cozman2020). Not stratified denotes programs with stratified negations and $\vee$ disjunction in the head

2.3 ProbLog and hybrid probLog

A ProbLog (De Raedt et al. Reference De RAEDT, Kimmig and Toivonen2007) program is composed by a set of Boolean probabilistic facts as described in Section 2.2 and a set of definite logical rules, and it is interpreted under the DS. The probability of a query $q$ is computed as the sum of the probabilities of the worlds where the query is true: every world has a unique model, so it is a sharp probability value.

Gutmann et al. (Reference Gutmann, Jaeger and De Raedt2011a) proposed Hybrid ProbLog, an extension of ProbLog with continuous facts of the form

\begin{equation*} (X,\phi )::b \end{equation*}

where $b$ is an atom, $X$ is a variable appearing in $b$ , and $\phi$ is a special atom indicating the continuous distribution followed by $X$ . An example of continuous fact is $(X,gaussian(0,1))::a(X)$ , stating that the variable $X$ in $a(X)$ follows a Gaussian distribution with mean 0 and variance 1. A Hybrid ProbLog program $P$ is a pair $(R,T)$ where $R$ is a set of rules and $T = T^c \cup T^d$ is a finite set of continuous ( $T^c$ ) and discrete ( $T^d$ ) probabilistic facts. The value of a continuous random variable $X$ can be compared only with constants through special predicates: $below(X,c_0)$ and $above(X,c_0)$ , with $c_0 \in \mathbb{R}$ , which succeed if the value of $X$ is, respectively, less than or greater than $c_0$ , and $between(X,c_0,c_1)$ , with $c_0, c_1 \in \mathbb{R}$ , $c_0 \lt c_1$ , which succeeds if $X \in [c_0,c_1]$ .

The semantics of Hybrid ProbLog is given in Gutmann et al. Reference Gutmann, Jaeger and De Raedt(2011a) a proof theoretic way. A vector of samples, one for each continuous variable, defines a so-called continuous subprogram, so the joint distribution of the continuous random variables defines a joint distribution (with joint density $f(\mathbf{x})$ ) over continuous subprograms. An interval $I \in \mathbb{R}^n$ , where $n$ is the number of continuous facts, is defined as the cartesian product of an interval for each continuous random variable, and the probability $P(\mathbf{X} \in I)$ can be computed by integrating $f(\mathbf{x})$ over $I$ , that is,

\begin{equation*} P(\mathbf {X} \in I) = \int _{I} f(\mathbf {x}) \, d\mathbf {x}. \end{equation*}

Given a query $q$ , an interval $I$ is called admissible if, for any $\mathbf{x}$ and $\mathbf{y}$ in $I$ and for any truth value of the probabilistic facts, the truth value of the $q$ evaluated in the program obtained by assuming $\mathbf{X}$ takes value $\mathbf{x}$ and $\mathbf{y}$ is the same. In other words, inside an admissible interval, the truth value of a query is not influenced by the values of the continuous random variables; that is, it is always true or false once the value of the discrete facts is fixed. For instance, if we have the continuous fact $(X,gaussian(0,1))::a(X)$ and a rule $q \,{:\!-}\, a(X), between(X,0,2)$ , the interval $[0,3]$ is not admissible for the query $q$ , while $[0.5,1.5]$ is (since for any value of $X \in [0.5,1.5]$ , $between(X,0,2)$ is always true). A partition is a set of intervals, and it is called admissible if every interval is admissible. The probability of a query is defined by extracting its proofs together with the probabilistic facts, continuous facts, and comparison predicates used in the proofs. These proofs are finitely many since the programs do not include function symbols and generalize the infinitely many proofs of a hybrid program (since there can be infinitely many values for continuous variables). Gutmann et al. (Reference Gutmann, Jaeger and De Raedt2011a) proved that an admissible partition exists for each query having a finite number of proofs and the probability of a query does not depend on the admissible partition chosen. They also propose an inference algorithm that first computes all the proofs for a query. Then, according to the continuous facts and comparison predicates involved, it identifies the needed partitions. To avoid counting multiple times the same contribution, the proofs are made disjoint. Lastly, the disjoint proofs are converted into a binary decision diagram (BDD) from which the probability can be computed by traversing it bottom up (i.e., with the algorithm adopted in ProbLog (De Raedt et al. Reference De RAEDT, Kimmig and Toivonen2007)).

Let us clarify it with a simple example.

Example 2. The following Hybrid ProbLog program has one discrete fact and one continuous fact.

Consider the query $q$ . It has two proofs: $a$ and $b(X), above(X,0.3)$ . The probability of $q$ is computed as $P(a) \cdot P(X \gt 0) + P(a) \cdot (1 - P(X \gt 0)) + (1 - P(a)) \cdot P(X \gt 0) = 0.4\cdot 0.5 + 0.4 \cdot 0.5 + 0.6 \cdot 0.5 = 0.7$ . $\square$

The Hybrid ProbLog framework imposes some restrictions on the continuous random variables and how they can be used in programs, namely: (i) comparison predicates must compare variables with numeric constants, (ii) arithmetic expressions involving continuous random variables are not allowed, and (iii) terms inside random variables definitions can be used only inside comparison predicates. These restrictions allow the algorithm for the computation of the probability to be relatively simple.

4.1 Exact inference

4.1.1 Modifying the PASTA solver

We first modified the PASTA solver (Azzolini et al. Reference Azzolini, Bellodi, Riguzzi, Gottlob, Inclezan and Maratea2022),Footnote ² by implementing the conversion of the hybrid program into a regular probabilistic answer set program (Algorithm1). PASTA performs inference on PASP via projected answer set enumeration. In a nutshell, to compute the probability of a query (without loss of generality assuming here it is an atom), the algorithm converts each probabilistic fact into a choice rule. Then, it computes the answer sets projected on the atoms for the probabilistic facts and for the query. In this way, PASTA is able to identify the answer set pertaining to each world. For each world, there can be three possible cases: (i) a single projected answer set with the query in it, denoting that the query is true in every answer set, so this world contributes to both the lower and upper probability; (ii) a single answer set without the query in it: this world does not contribute to any probability; and (iii) two answer sets, one with the query and one without: the world contributes only to the upper probability. There is a fourth implicit and possible case: if a world is not present, this means that the ASP obtained from the PASP by fixing the selected probabilistic facts is unsatisfiable. Thus, in this case, we need to consider the normalization factor of equation (9). PASTA already handles this by keeping track of the sum of the probabilities of the computed worlds. The number of generated answer sets depends on the number of Boolean probabilistic facts and on the number of intervals for the continuous random variables, since every interval requires the introduction of a new Boolean probabilistic fact. Overall, if there are $d$ discrete probabilistic facts and $c$ continuous random variables, the total number of probabilistic facts (after the conversion of the intervals) becomes $T = d + \sum _{i = 1}^{c}(k_i - 1)$ , where $k_i$ is the number of intervals for the $i$ -th continuous fact. Finally, the total number of generated answer sets is bounded above by $2^{T+1}$ , due to the projection on the probabilistic facts. Clearly, generating an exponential number of answer sets is intractable, except for trivial domains.

Example 7 (Inference with PASTA.) Consider the discretized program described in Example 5. It is converted into

It has 10 answer sets projected on the atoms $q0/0$ and $a0/0$ , $a1/0$ , and $b/0$ : $AS_{1} = \{\}$ , $AS_{2} = \{a1\}$ , $AS_{3} = \{a2\}$ , $AS_{4} = \{a2,q0\}$ , $AS_{5} = \{b\}$ , $AS_{6} = \{a1,a2\}$ , $AS_{7} = \{a1,b\}$ , $AS_{8} = \{a2,b,q0\}$ , $AS_{9} = \{a1,a2,q0\}$ , and $AS_{10} = \{a1,a2,b,q0\}$ . For example, the world where only $a2$ is true is represented by the answer sets $AS_3$ and $AS_4$ . The query $q0$ is present only in one of the two, so this world contributes only to the upper probability. The world where $a1$ and $a2$ are true and $b$ is false is represented only by the answer set $AS_9$ : the query is true in it so this world contributes to both the lower and upper probability. The world where only $a1$ is true is represented by $AS_2$ , $q0$ is not present in it so it does not contribute to the probability. By applying similar considerations for all the worlds, it is possible to compute the probability of the query $q0$ . $\square$

4.1.2 Modifying the aspcs solver

We also added the conversion from HPASP to PASP to the aspcs solver (Azzolini and Riguzzi Reference Azzolini, Riguzzi, Basili, Lembo, Limongelli and Orlandini2023a), built on top of the aspmc solver (Eiter et al. Reference Eiter, Hecher, Kiesel, Bienvenu, Lakemeyer and Erdem2021; Kiesel et al. Reference Kiesel, Totis and Kimmig2022), which performs inference on Second Level Algebraic Model Counting (2AMC) problems, an extension of AMC (Kimmig et al. Reference Kimmig, Van Den Broeck and De Raedt2017). Some of them are MAP inference (Shterionov et al. Reference Shterionov, Renkens, Vlasselaer, Kimmig, Meert and Janssens2015), decision theory inference in PLP (Van den Broeck et al. Reference Van Den Broeck, Thon, Van Otterlo and De Raedt2010), and probabilistic inference under the smProbLog semantics (Totis et al. Reference Totis, De Raedt and Kimmig2023). More formally, let $X_{in}$ and $X_{out}$ be a partition of the variables in a propositional theory $\Pi$ . Consider two commutative semirings $\mathcal{R}_{in} = (R^{i},\oplus ^i,\otimes ^i,e_{\oplus }^i,e_{\otimes }^i)$ and $\mathcal{R}_{out} = (R^{o},\oplus ^o,\otimes ^o,e_{\oplus }^o,e_{\otimes }^o)$ , connected by a transformation function $f : \mathcal{R}^{i} \to \mathcal{R}^{o}$ and two weight functions, $w_{in}$ and $w_{out}$ , which associate each literal to a weight (i.e., a real number). 2AMC requires computing:

(10) \begin{align} \begin{split} 2AMC(A) =& \bigoplus \nolimits _{I_{out} \in \sigma (X_{out})}^{o} \bigotimes \nolimits ^{o}_{a \in I_{out}} w_{out}(a) \otimes ^{o} \\[5pt] & f( \bigoplus \nolimits _{I_{in} \in \delta (\Pi \mid I_{out})}^{i} \bigotimes \nolimits ^{i}_{b \in I_{in}} w_{in}(b) ) \end{split} \end{align}

where $\sigma (X_{out})$ are the set of possible assignments to $X_{out}$ and $\delta (\Pi \mid I_{out})$ are the set of possible assignments to the variables of $\Pi$ that satisfy $I_{out}$ . We can cast inference under the CS as a 2AMC problem (Azzolini and Riguzzi Reference Azzolini, Riguzzi, Basili, Lembo, Limongelli and Orlandini2023a) by considering as inner semiring $\mathcal{R}_{in} = (\mathbb{N}^2, +, \cdot, (0,0), (1,1))$ , where $+$ and $\cdot$ are component-wise and $w_{in}$ is a function associating $not\, q$ to $(0,1)$ and all the other literals to $(1,1)$ , where $q$ is the query. The first component $n_u$ of the elements $(n_u,n_l)$ of the semiring counts the models where the query is true, while the second component $n_l$ counts all the models. The transformation function $f(n_u,n_l)$ returns a pair of values $(f_u,f_l)$ such that $f_l = 1$ if $n_l = n_u$ , 0 otherwise, and $f_u = 1$ if $n_u \gt 0$ , 0 otherwise. The outer semiring $\mathcal{R}_{out} = ([0, 1]^2, +, \cdot, (0,0), (1,1))$ is a double probability semiring, where there is a separate probability semiring for each component. $w_{out}$ associates $a$ to $(p,p)$ and $not\ a$ to $(1-p,1-p)$ for every probabilistic fact $p :: a$ , while it associates all the remaining literals to $(1,1)$ . aspmc, and so aspcs, solves 2AMC using knowledge compilation (Darwiche and Marquis Reference Darwiche and Marquis2002) targeting NNF circuits (Darwiche Reference Darwiche2004) where the order of the variables is guided by the treewidth of the program (Eiter et al. Reference Eiter, Hecher, Kiesel, Bienvenu, Lakemeyer and Erdem2021). Differently from PASTA, aspcs does not support aggregates, but we can convert rules and constraints containing them into new rules and constraints without aggregates using standard techniques (Brewka et al. Reference Brewka, Eiter and Trusczynski2011).

4.2 Approximate inference

Approximate inference can be performed by using the definition of the sampling semantics and returning the results after a finite number of samples, similar to what is done for programs under the DS (Kimmig et al. Reference Kimmig, Demoen, De Raedt, Costa and Rocha2011; Riguzzi Reference Riguzzi2013). It can be performed both on the discretized program and directly on the hybrid program.

4.2.3 Approximate algorithm description

Algorithm2 illustrates the pseudocode for the sampling procedure: for the sampling on the discretized program, the algorithm discretizes the program by calling Discretize (Algorithm1). Then, for a number of samples $s$ , it samples a world by including or not every probabilistic fact into the program (function SampleWorld) according to the probability values, computes its answer sets with function ComputeAnswerSets (recall that a world is an answer set program), checks whether the query $q$ is true in every answer set (function QueryInEveryAnswerSet) or in at least one answer set (function QueryInAtLeastOneAnswerSet), and updates the lower and upper probability bounds accordingly. At the end of the $s$ iterations, it returns the ratio between the number of samples contributing to the lower and upper probability and the number of samples taken. The procedure is analogous (but without discretization) in the case of the sampling on the original program. A world is sampled with function SampleVariablesAndTestConstraints: it takes a sample for each continuous random variable, tests that value against each constraint specified in the program, and removes it if the test succeeds; otherwise, it removes the whole rule containing it. The remaining part of the algorithm is the same.

We now present two results regarding the complexity of the sampling algorithm. The first provides a bound on the number of samples needed to obtain an estimate of the upper (or lower) probability of a query within a certain absolute error. The second result provides a bound on the number of samples needed to obtain an estimate within a certain relative error.

Theorem 3 (Absolute Error.) Let $q$ be a query in a hybrid probabilistic answer set program $P$ whose exact lower (upper) probability of success is $p$ . Suppose the sampling algorithm takes $s$ samples, $k$ of which are successful, and returns an estimate $\hat{p} = \frac{k}{s}$ of the lower (upper) probability of success. Let $\epsilon$ and $\delta$ be two numbers in $[0,1]$ . Then, the probability that $\hat{p}$ is within $\epsilon$ of $p$ is at least $1-\delta$ , that is,

\begin{equation*}P(p-\epsilon \leq \hat {p} \leq p+\epsilon ) \geq 1 - \delta \end{equation*}

if

\begin{equation*} s\geq \frac {\epsilon +\frac {1}{2}}{\epsilon ^2 \delta }. \end{equation*}

Algorithm 2. Function Sampling: computation of the probability of a query q by taking s samples in a hybrid probabilistic answer set program P with rules R, continuous probabilistic facts C, and discrete probabilistic factsD. Variable type indicates whether the sampling of the discretized program or the original program is selected.

Proof. We must prove that

(11) \begin{equation} P(p-\epsilon \leq \hat{p}\leq p+\epsilon )\geq 1 - \delta \end{equation}

or, equivalently, that

(12) \begin{equation} P(sp-s\epsilon \leq k\leq sp+s\epsilon ) \geq 1 - \delta . \end{equation}

Since $k$ is a binomially distributed random variable with number of trials $s$ and probability of success $p$ , we have that (Feller Reference Feller1968):

(13) \begin{equation} P(k\geq r_1) \leq \frac{r_1(1-p)}{(r_1-sp)^2} \end{equation}

if $r_1\geq sp$ . Moreover

(14) \begin{equation} P(k\leq r_2) \leq \frac{(s-r_2)p}{(sp-r_2)^2} \end{equation}

if $r_2\leq sp$ . Since $P(k \geq r_1) = 1-P(k \lt r_1)$ , from equation (13), we have

(15) \begin{eqnarray} 1-P(k\lt r_1)&\leq & \frac{r_1(1-p)}{(r_1-sp)^2}\nonumber \\[5pt] P(k\lt r_1)&\geq & 1-\frac{r_1(1-p)}{(r_1-sp)^2} \end{eqnarray}

if $r_1\geq sp$ .

In our case, we have $r_1 = sp+s\epsilon$ (since $r_1 \geq sp$ ) and $r_2 = sp-s\epsilon$ (since $r_2 \leq sp$ ). So

\begin{eqnarray*} && P(p-\epsilon \leq \hat{p}\leq p+\epsilon ) =\\[5pt] & = & P(r_2 \leq k \leq r_1) \\[5pt] & = &P(k\leq r_1)-P(k\lt r_2) \\[5pt] &\geq &1-\frac{r_1(1-p)}{(r_1-sp)^2}-\frac{(s-r_2)p}{(sp-r_2)^2} \qquad \text{(Eq.\,15 and\,14 and since $P(x \leq v) \geq P(x \lt v)$)} \\[5pt] &=&1-\frac{(sp+s\epsilon )(1-p)}{(s\epsilon )^2}-\frac{(s-sp+s\epsilon )p}{(s\epsilon )^2} \quad \text{(by replacing the values of $r_1$ and $r_2$)} \\[5pt] &=&\frac{s^2\epsilon ^2-sp-s\epsilon +sp^2+sp\epsilon -sp+sp^2-sp\epsilon }{s^2\epsilon ^2} \quad \text{(by expanding)} \\[5pt] &=&\frac{s\epsilon ^2-2p-\epsilon +2p^2}{s\epsilon ^2} \quad \text{(by collecting $s$ and simplifying)}\\[5pt] \end{eqnarray*}

and

(16) \begin{eqnarray} \frac{s\epsilon ^2-2p-\epsilon +2p^2}{s\epsilon ^2} & \geq & 1-\delta \nonumber \\[5pt] s\epsilon ^2-2p-\epsilon +2p^2 & \geq & s\epsilon ^2- s\epsilon ^2 \delta \nonumber \\[5pt] s\epsilon ^2\delta & \geq & 2p + \epsilon - 2p^2\nonumber \\[5pt] s & \geq & \frac{2p + \epsilon - 2p^2}{\epsilon ^2\delta } = \frac{2p(1-p) + \epsilon }{\epsilon ^2\delta }. \end{eqnarray}

Since $0\leq 2p (1-p) \leq \frac{1}{2}$ with $p \in [0,1]$ , then equation (16) is implied by

\begin{eqnarray*} s & \geq & \frac{\epsilon +\frac{1}{2}}{\epsilon ^2 \delta } \end{eqnarray*}

Theorem 4 (Relative Error.) Let $q$ be a query in a hybrid probabilistic answer set program $P$ whose exact lower (upper) probability of success is $p$ . Suppose the sampling algorithm takes $s$ samples, $k$ of which are successful, and returns an estimate $\hat{p} = \frac{k}{s}$ of the lower (upper) probability of success. Let $\epsilon$ and $\delta$ be two numbers in $[0,1]$ . Then, the probability that the error between $\hat{p}$ and $p$ is smaller than $p\epsilon$ is at least $1-\delta ^p$ , that is,

\begin{equation*}P(|p-\hat {p}| \leq p\epsilon ) \geq 1 - \delta ^p\end{equation*}

if

\begin{equation*} s\geq \frac {3}{\epsilon ^2}\ln (\frac {1}{\delta }). \end{equation*}

Proof. According to Chernoff’s bound (Mitzenmacher and Upfal Reference Mitzenmacher and Upfal2017, Corollary 4.6), we have that

\begin{equation*} P( |ps-k| \ge ps\epsilon ) \leq 2e^{-\frac {\epsilon ^2ps}{3}}. \end{equation*}

So

\begin{equation*} P( |p - \hat {p}| \geq p\epsilon ) \le 2e^{-\frac {\epsilon ^2ps}{3}} \end{equation*}

and

\begin{equation*} P( |p - \hat {p}| \le p\epsilon ) \geq 1-2e^{-\frac {\epsilon ^2ps}{3}}. \end{equation*}

Then,

\begin{eqnarray*} 1-2e^{-\frac{\epsilon ^2ps}{3}} & \geq & 1-\delta ^p\\[2pt] 2e^{-\frac{\epsilon ^2ps}{3}} & \leq & \delta ^p\\[2pt] \ln (2) -\frac{\epsilon ^2ps}{3} & \leq & p \ln (\delta ) \\[2pt] -\frac{\epsilon ^2ps}{3} & \leq & p \ln (\delta ) - \ln (2) \leq p \ln (\delta ) \\[2pt] s & \geq & \frac{3}{\epsilon ^2}\ln (\frac{1}{\delta }). \end{eqnarray*}

Please note that in Theorem4, we exponentiate $\delta$ to $p$ . This is needed to avoid the appearance of $p$ in the bound since $p$ is unknown. When $p$ is 0, we have no guarantees on the error. When $p$ is 1, the confidence is $1-\delta$ . For $p$ growing from 0–1, the bound provides increasing confidence. In fact, approximating small probabilities is difficult due to the low success probability of the sample.

5.1 Dataset t1

In $t1$ we consider instances of increasing size of the program shown in Example4. Every instance of size $n$ has $n/2$ discrete probabilistic facts $d_i$ , $n/2$ continuous random variables $c_i$ with a Gaussian distribution with mean 0 and variance 1, $n/2$ pair of rules $q0 \,{:\!-}\, below(c_i, 0.5), not \, q1$ and $q1 \,{:\!-}\, below(c_i, 0.5), not \, q0$ , $i = \{1,\dots, n/2\}$ , and $n/2$ rules $q0 \,{:\!-}\, below(c_i, 0.7), d_i$ , one for each $i$ for $i=\{1,\dots, n/2\}$ . The query is $q0$ . The goal of this experiment is to analyze the behaviour of the algorithm by increasing the number of discrete and continuous random variables. The instance of size 2 is:

5.2 Dataset t2

In $t2$ we consider a variation of $t1$ , where we fix the number $k$ of discrete probabilistic facts to 2, 5, 8, and 10, and increase the number of continuous random variables starting from 1. Every instance of size $n$ contains $k$ discrete probabilistic facts $d_i$ , $i \in \{0,\dots, k-1\}$ , $n$ continuous random variables $c_i$ , $i \in \{1,\dots, n\}$ with a Gaussian distribution with mean 0 and variance 1, $n$ pair of rules $q0 \,{:\!-}\, below(c_i, 0.5), not \, q1$ and $q1 \,{:\!-}\, below(c_i, 0.5), not \, q0$ , $i = \{1,\dots, n\}$ , and $n$ rules $q0 \,{:\!-}\, below(c_i, 0.7), d_{(i-1) \% k}$ , one for each $i$ for $i=\{1,\dots, n\}$ . The query is $q0$ . Here, when the instance size $n$ is less than the number of discrete probabilistic facts $k$ , $k - n$ probabilistic facts are not relevant for the computation of the probability of the query. The instance of size 2 with $k = 2$ is:

5.3 Dataset t3

In $t3$ , we consider again a variation of $t1$ where we fix the number $k$ of continuous random variables to 2, 5, 8, and 10, and increase the number of discrete probabilistic facts starting from 1. Every instance of size $n$ contains $k$ continuous probabilistic facts $c_i$ , $i \in \{0,\dots, k-1\}$ , with a Gaussian distribution with mean 0 and variance 1, $n$ discrete probabilistic facts $d_i$ , $i \in \{1,\dots, n\}$ , $k$ pair of rules $q0 \,{:\!-}\, below(c_i, 0.5), not \, q1$ and $q1 \,{:\!-}\, below(c_i, 0.5), not \, q0$ , $i = \{1,\dots, n\}$ , and $n$ rules $q0 \,{:\!-}\, below(c_{(i-1) \% k}, 0.7), d_i$ , one for each $i$ for $i=\{1,\dots, n\}$ . The query is $q0$ . The instance of size 3 with $k = 2$ is:

5.4 Dataset t4

In $t4$ , we consider programs with one discrete probabilistic fact $d$ and one continuous random variable $a$ that follows a Gaussian distribution with mean 0 and standard deviation 10. An instance of size $n$ has $n$ pairs of rules $q0 \,{:\!-}\, between(a,lb_i,ub_i), not \, q1$ and $q1 \,{:\!-}\, between(a,lb_i,ub_i), not \, q0$ where $lb_i$ and $ub_i$ are randomly generated (with, $lb_i \lt ub_i$ ), for $j = 1,\dots, n$ , and $n$ rules of the form $q0 \,{:\!-}\, d, between(a,LB_j,UB_j)$ , where the generation of the $LB_j$ and $UB_j$ follows the same process of the previous rule. We set the minimum value of $lb_i$ and $LB_j$ to -30 and, for both rules, the lower and upper bounds for the $\mathit{between}/2$ comparison predicate are uniformly sampled in the range $[u_{i-1}, u_{i-1} + 60/n]$ , where $u_{i-1}$ is the upper bound of the previous rule. The query is $q0$ . Here, the goal is to test the algorithm with an increasing number of intervals to consider. An example of an instance of size 2 is:

5.5 Dataset t5

Lastly, in $t5$ we consider programs of the form shown in Example3: the instance of index $n$ has $n$ people involved, $n$ discrete probabilistic facts, and $n$ $d/1$ and $s/1$ continuous random variables. The remaining part of the program is the same. Example3 shows the instance of index 4 (four people involved). The query is $high\_number\_strokes$ .

Fig 1. Results for $t1$ , $t4$ , and $t5$ (left) and results for aspcs applied to $t3$ (right) with two and five continuous facts.

Fig 2. Results for the experiment $t2$ with a fixed number of discrete facts and an increasing number of continuous variables.

Fig 3. Results for the experiment $t3$ with a fixed number of continuous random variables and an increasing number of discrete facts. The x axis of the left plot is cut at size 20, to keep the results readable.

Table 4. Results for the approximate algorithms based on sampling. The columns contain, from the left, the name of the dataset, the instance index, and the average time required to take $10^2$ , $10^3$ , $10^4$ , $10^5$ , and $10^6$ samples on the original and converted program. O.O.M. and T.O. stand, respectively, for out of memory and timeout

Table 5. Standard deviations for the results listed in Table 4. A dash denotes that there are no results for that particular instance due to either a memory error or a time limit

Table 6. Largest solvable instances of $t4$ by sampling the converted program when reducing the available memory. The columns contain, from the left, the maximum amount of memory, the largest solvable instance together with the number of probabilistic facts (# p.f.) and rules (# rules) obtained via the conversion, and the maximum number of samples that can be taken (max. # samples). Note that we increase the number of samples by starting from $10^2$ and iteratively multiplying the number by 10, up to $10^6$ , so in the last column, we report a range: this means that we get a memory error with the upper bound while we can take the number of samples in the lower bound. Thus, the maximum values of samples lie in the specified range

5.6 Exact inference results

The goal of benchmarking the exact algorithms is threefold: (i) identifying how the number of continuous and discrete probabilistic facts influences the execution time, (ii) assessing the impact on the execution time of an increasing number of intervals, and (iii) comparing knowledge compilation with projected answer set enumeration. Figure 1(a) shows the results of exact inference for the algorithms backed by PASTA (dashed lines) and aspcs (straight lines) on $t1$ , $t4$ , and $t5$ . For $t5$ PASTA is able to solve up to instance 4, while aspcs can solve up to instance 9 (in a few seconds). Similarly for $t1$ , where aspcs doubles the maximum sizes solvable by PASTA. The difference is even more marked for $t4$ : here, PASTA solves up to instance 11 while aspcs up to instance 35. Figures 2 and 3 show the results for $t2$ and $t3$ . Here, PASTA cannot solve instances with around more than 20 probabilistic facts and continuous random variables combined. The performance of aspcs are clearly superior, in particular for $t3$ with 2 and 5 continuous random variables (Figure 3(a)). Note that this is the only plot where we cut the x axis at instance 20, to keep the results readable. To better investigate the behavior of aspcs in this case, we keep increasing the number of discrete facts until we get a memory error, while the number of continuous variables is fixed to 2 and 5. Figure 1(b) shows the results: with 2 continuous variables, we can solve up to instance 110, while with 5 continuous variables, we can only solve up to 75. In general, the results of exact inference are in accordance with the theoretical complexity results. However, as also empirically shown by Azzolini and Riguzzi (Reference Azzolini, Riguzzi, Basili, Lembo, Limongelli and Orlandini2023a), knowledge compilation has a huge impact on the execution times. Overall, as expected, PASTA is slower in all the tests, being based on (projected) answer set enumeration. For all, both PASTA (exact algorithm) and aspcs stop for lack of memory.

5.7 Approximate inference results

The goal of the experiments run with approximate algorithms is: (i) analyzing the impact on the execution time of the number of samples taken, (ii) comparing the approach based on sampling the original program against the one based on sampling the converted program in terms of execution times, and (iii) assessing the memory requirements. Table 4 shows the averages on five runs of the execution time of the approximate algorithm applied to both the discretized and original program with $10^2$ , $10^3$ , $10^4$ , $10^5$ , and $10^6$ samples, for each of the five datasets. Standard deviations are listed in Table 5. For $t2$ and $t3$ , we consider programs with the same number of continuous variables and discrete probabilistic facts. In four of the five tests (all except $t4$ ), sampling the original program is slower than sampling the converted program, sometimes by a significant amount. This is probably due to the fact that sampling a continuous distribution is slower than sampling a Boolean random variable. For example, in instance 100 of $t3$ , sampling the original program is over 12 times slower than sampling the converted program. However, the discretization process increases the number of probabilistic facts and so the required memory: for $t4$ , from the instance 70, taking even 100 samples requires more than 8 GB of RAM. To better assess the memory requirements, we repeat the experiments with the approximate algorithms with 6, 4, 2, and 1 GB of maximum available memory. For all, taking up to $10^5$ samples in the original program is feasible also with only 1 GB of memory. The same considerations hold for sampling the converted program, except for $t4$ . Table 6 shows the largest solvable instances together with the number of probabilistic facts, rules (for the converted program), and samples for each instance: for example, with 1 GB of memory, it is possible to take only up to $10^4$ samples in the instance of size 40. For this dataset, the increasing number of $between/3$ predicates requires an increasing number of rules and probabilistic facts to be included in the program during the conversion, to properly handle all the intervals: the instance of size 70 has 142 probabilistic facts and more than 30,000 rules, which make the generation of the answer sets very expensive and sampling it with only 8 GB of memory is not feasible.