Dear all,

suppose I want to create a distribution for a discrete multivariate random variable that has a random field distribution:

X \sim \frac{1}{Z} \exp\left( \sum_c \phi_c(X_c) \right),

where \phi is a potential function. I believe I cannot really implement this, since I would need to implement the inherited `logp`

function which in this case is not tractable to compute:

```
class RF(Discrete):
def __init__(self, cliques, **kwargs):
# cliques is a list of list of indexes
...
def logp(self, value):
??
```

Thus models like this cannotâ€™ work:

```
with pm.Model() as model:
A = RF('A', cliques)
B = pm.Normal('B', mu=A, sd=10, shape=len(cliques))
```

My question is: is there a way to manually sample from the MRF and then use the samples for the model block below, for instance by creating a custom sampler that inherits from `pm.Metropolis`

and then use this in a compound step? I.e. I maybe sth like this:

```
with model:
step1 = pm.MyCustomMetropolis([A])
step2 = pm.Metropolis([B])
trace = pm.sample(..., step=[step1, step2]...)
```

Thanks and sorry if this sounds confusing,

Simon