Questions about defining a custom node and rules example.

[Flawless1202]

Hi,
Thanks for the detailed tutorial on defining a custom node and rules. I followed this example, but encountered some confusion regarding the necessity of the rules.

To investigate this, I commented some of the rules and discovered that the following rule alone seems to be necessary:

@rule MyBernoulli(:π, Marginalisation) (q_out :: Any,) = begin
    p = mean(q_out)
    return Beta(one(p) + p, 2one(p) - p)
end

Could you please provide an explanation for this? And how can I determine which rules are essential for the custom node? @ThijsvdLaar @bvdmitri @alstat and others. Thanks a lot!

p.s.
The full code is as follows:

using RxInfer
using Random

# struct for Bernoulli distribution with success probability π
struct MyBernoulli{T <: Real} <: ContinuousUnivariateDistribution
    π :: T
end

# for simplicity, let's also specify the mean of the distribution
Distributions.mean(d::MyBernoulli) = d.π

@node MyBernoulli Stochastic [out, (π, aliases = [p])]

@rule MyBernoulli(:π, Marginalisation) (q_out :: Any,) = begin
    p = mean(q_out)
    return Beta(one(p) + p, 2one(p) - p)
end

rng = MersenneTwister(42)
n = 500
π_real = 0.75
distribution = Bernoulli(π_real)

dataset = float.(rand(rng, distribution, n))

@model function coin_model_mybernoulli(y)
    # We endow θ parameter of our model with some prior
    π ~ Beta(4.0, 8.0)
    # We assume that outcome of each coin flip is governed by the MyBernoulli distribution
    for i in eachindex(y)
        y[i] ~ MyBernoulli(π)
    end
end

result_mybernoulli = infer(
    model = coin_model_mybernoulli(),
    data  = (y = dataset, )
)

@model function coin_model(y)
    p ~ Beta(4.0, 8.0)
    for i in eachindex(y)
        y[i] ~ Bernoulli(p)
    end
end

result_bernoulli = infer(
    model = coin_model(),
    data  = (y = dataset, ),
)

if !(result_bernoulli.posteriors[:p] == result_mybernoulli.posteriors[:π])
    error("Results are not identical")
else
    println("Results are identical 🎉🎉🎉")
end

And the output is:

Results are identical 🎉🎉🎉

[Raphael-Tresor]

Hi,
Thanks for your question Flawless1202.
The example proposes multiple rules to compute the message toward the variable out (the result of the Bernoulli) and multiple rules to compute the message toward the variable \pi (the parameter of the Bernoulli). Only the rules needed for the computation will be called when running inference, and the others will be ignored. In your case, you infer a posterior for the variable \pi, so only a rule toward \pi is required. Indeed, messages pass from the data variable (in your case, y ) to the hidden variable (in your case, \pi).

In general, necessary rules depend on both the algorithm used to infer and the model in which the node is used. When all variables attached to a node have been updated except one, RxInfer will always try to update the last one by calling a rule of this node. Note that feeding a data variable with a new value is an update. This implies that after feeding each data variable a new value, messages will fly through the graph, updating once each variable. This way, you can infer a posterior with the belief propagation algorithm, which computes Baye's posterior.

To determine the required rules, I advise you to draw the graphical representation of your model like in the Getting started documentation page. If you want to run the belief propagation algorithm, you must write a rule for any message that flows out of a node to reach an hidden variable. I wrote an abstract example below and hope that it will help. Data variables are noted $d_i$ and hidden variables $v_i$.

model
$$p( v_1, v_2, v_3, d_1, d_2, d_3,) = \underbrace{p(v_3 \mid v_2, v_1)} \underbrace{p(v_1 \mid d_1, d_2)} \underbrace{p(v_ 2\mid d_3 )}$$
----------------------------node C-------node A-----node B-

graph

For instance, for the node A, you only need
@rule NodeA(:v_1, Marginalisation) (d_1 :: PointMass, d_2 :: PointMass),
and for node C you only need
@rule NodeC(:v_3 Marginalisation) (v_1 :: Any, v_2 :: Any).

For more information, please read the documentation about inference and don't hesitate to come back to us.

[Flawless1202]

Hi @Raphael-Tresor, thanks for your patient reply. I think I have got it. And I have one more question about the constraint specification. Also for this example, I draw the factor graph with 10 observations y as follows:

As the tutorial says:

In the example that we will show later on, we solely use sum-product message passing. Variational message passing requires us to set the local constraints in our model, something which is out of scope of this tutorial.

Why there are no constraints set, but the codes still work well? Could you please give me an explain?

[bvdmitri]

Why there are no constraints set, but the codes still work well? Could you please give me an explain?

We have a section about constraints specification here. Generally, the explanation is simple: sometimes exact inference is computationally feasible, and sometimes it is not (within a reasonable time frame). If exact inference is not possible, we may impose extra constraints on the inference to simplify computation at the cost of accuracy. It depends on the model specification. In this particular case, exact inference is possible with the provided rules, so constraints are not required.