Changing a parameter that isn't on the chain every sample


Sorry if this is a very easy solution - I searched for a few hours and wasn’t able to find a solution that I understood.

I was wondering if there was a way to introduce uncertainty in a parameter that is being used in the model, without adding it to the chain. For example, suppose I have a linear relationship between two parameters x and y, such that

y = Ax + B, where A = -0.006 +/- 0.001 and B = 1 +/- 0.04. I am sampling x and want to change A and B each chain sample so I can incorporate the uncertainty in my input parameters into the model uncertainty. Is there a way to do this? Thanks.

EDIT: Sorry all, it looks like my phrasing above was confusing. In pseudo-code, here’s what I want:

import numpy as np
import pymc3 as pm

# Dummy data
Y = np.arange(0, 50) + np.random.randn(50)

for sample in int(nburnin + nsamples):
    with pm.Model() as model:
        A = np.random.normal(-0.006, 0.001)
        B = np.random.normal(1.0, 0.04)

        y_true = pm.Data('y', Y)
        x = pm.Flat('x')
        y_syn = A*x + B
       #... augment likelihood to look at difference between y_true and y_syn ...
        trace = pm.sample(...)

Of course the above code wouldn’t work the way it is intended with PyMC3, but this demonstrates my point visually, I think.

I don’t want to use pm.Normal('A', mu=-0.006, sigma=0.001 and pm.Normal('B', mu=1, sigma=0.04) because (my understanding) is that this would sample those parameters to solve them in the MCMC. I just want the values and their uncertainties added into my model solving for x.

it depends on your statistical model. for example, if A is normally distributed, uncertainty is the ‘sigma’ parameter. (Note: I’m not sure if I understand your question properly!)

Not that I know of. But note that 1 chain is just 1 sampling of the model, so you could just perturb A and B between calls to pm.sample something like

for j in range(nchains):
    A += np.random.(...)   # I *think* += operates on the underlying memory rather than creating a new object
    B += np.random.(...)
    x_chains.append(pm.sample(..., chains=1))

1 Like

Can’t you just add a prior to A and B that has the noise distribution you want?

A_noisy = pm.Normal("A_noisy", A, 1)

Where the prior can be any distribution centered around your value A

This has the advantage that your sample can take into consideration gradient information. Otherwise you are in untested territory in PyMC3.

I think this is the closest thing to what I need, thanks. I also had the idea of perturbing the A & B each chain, but I thought I’d ask if there was a way to do it each sample.


Sorry for the late reply. Wouldn’t using A_noisy = pm.Normal("A_noisy", A, 1) allow the sampling to change the distribution of A over the chain? I would want my posterior distribution for A to be the same as its prior.

I am not sure what your goal is. You want a prior that’s fixed but random? The way I would think about it would be to have fixed parameters (if the approximation is not too bad). So you don’t have a fixed A with a special pdf but a distribution with fixed parameters. Say A_mean = 10 and A_std = 1. Then A_noisy = pm.Normal("A_noisy", A_mean, A_std). The latter would change obviously during inference but A itself would be fixed as those are just numbers. Then the prior for A_noisy2 = pm.Normal("A_noisy", A_mean, A_std) in another context is exactly the same as A_noisy1 was. It’s a different way of thinking about it but I feel it achieves the same result, no? If you don’t like to force a Normal prior you can use something else.

Basically, I have two parameters A and B that have known value and uncertainty. I don’t need to sample in the MCMC, since their value is already determined. I would, however, like to incorporate their uncertainties into my model. Ideally, I would accomplish this by changing the A and B parameters each sample, but not updating these parameters each step.

I hope that makes more sense.

I think I understand the idea (although not the rationale). I really don’t know how wise it is to change things that pymc expects to be fixed during sampling.

I still feel that what you want is equivalent in practice to just resetting your prior of A every time you sample, but I might be missing something. How would the generative model for your data look like?

Actually I think you’re correct - resetting the prior each sample would achieve the same thing. I don’t see how to do that using pymc3.sample(), however.