Add FixedDistribution helper

[ricardoV94]

Some users on discourse have come up with models that rely on fixed distributions (i.e., RVs that should not be treated as standard variables that you want to infer):

https://discourse.pymc.io/t/custom-sampling-step/7243
https://discourse.pymc.io/t/two-stage-bayesian-regression-enforcing-a-fixed-distribution-not-just-hierarchical-regression/10036
https://discourse.pymc.io/t/using-posterior-as-likelihood/10742 (this one is not as clear but in DM with the user it became clear she wanted to do something similar)

In theory one could simply use .dist:

import pymc as pm

with pm.Model() m:
    fixed_normal_dist = pm.Deterministic("fixed_normal_dist", pm.Normal.dist())
    noise = pm.Gamma("noise", alpha=2, beta=1)
    pm.Normal("llike", mu=fixed_normal_dist, sigma=noise, observed=[0, 1, 2])

    idata = pm.sample()

Right now we simply error out if we find RVs in the logp graph (except for SimulatorRVs)

ValueError: Random variables detected in the logp graph: [normal_rv{0, (0, 0), floatX, False}.out].

We could easily remove that check... but I don't think NUTS would appreciate that (@ColCarroll, @junpenglao, @aseyboldt?). Instead, what sounds safer is to use a different, but fixed, draw of fixed_normal_dist in every NUTS step. In that case, what I have been proposing is some hack like this:

import pymc as pm

with pm.Model() m:
    # Use Flat so we don't add any term to the model logp
    fixed_normal_dist = pm.Flat("fixed_normal_dist")
    # In every iteration CustomFixedStep assigns a new Normal draw to `fixed_normal_dist`
    # Care must be taken with seeding in multiprocessing
    step = CustomFixedStep(vars={fixed_normal_dist: Normal.dist()})

    noise = pm.Gamma("noise", alpha=2, beta=1)
    pm.Normal("llike", mu=fixed_normal_dist, sigma=noise, observed=[0, 1, 2])

    idata = pm.sample(step=step)

Does it make sense to add a helper pm.FixedDist(name, rv) that does something similar behind the scenes? It could obviously go into pymc_experimental first, but it's easier to bring up the discussion here.

Count the logp of the fixed dist or not?

Also it's still a bit fuzzy for me if we always/never want to consider the FixedDistribution logp in the rest of model, or if in some cases we do and others we don't... We could of course, give the choice to the user...

If we wanted to count the example above, would look like:

with pm.Model() m:

    fixed_normal_dist = pm.Normal("fixed_normal_dist")
    # In every iteration CustomFixedStep assigns a new Normal draw to `fixed_normal_dist`,
    # without any concern for the model joint_logp.
    step = CustomFixedStep(vars=[fixed_normal_dist])
    noise = pm.Gamma("noise", alpha=2, beta=1)
    pm.Normal("llike", mu=fixed_normal_dist, sigma=noise, observed=[0, 1, 2])

    idata = pm.sample(step=step)

[ricardoV94]

Another use case: https://discourse.pymc.io/t/can-pm-draw-be-used-inside-a-pymc-model-or-will-it-break-links-between-rvs-and-observed-likelihood-variable/12069