Markov Chain Monte Carlo with known normalisation

Question

I would like to compute the expectation value $\langle O \rangle = \sum_x P(x) O(x)$ of some random variable over an extremely large sample space that I cannot simply exhaustively go through. Usually I would use Metropolis Monte Carlo for this, but for this particular example I already know the exact normalized probability $P(x)$ with which each sample $x$ occurs.

How can I leverage this? Does this simply mean I do not need to wait for equilibration of my Markov Chain? I feel like knowing the normalization should allow me to use a much more efficient algorithm/method here. Ideally I would like to directly generate samples $x_i$ according to the distribution of $P(x)$ and then compute $1/N \sum_{i=1}^N O(x_i)$ but I do not see an easy way of doing this, apart from MCMC which does not utilise the full information I have available.

Thank you!

Answer 1

Another way of looking at the issue of approximating$$\mathfrak I = \sum_{x\in\mathfrak X} p(x)O(x)$$by stochastic techniques is to aim at adding primarily large values of $p(x)O(x)$. Assuming no information is available about the [location or values of the] largest probabilities over $\mathfrak X$, one could consider a self-avoiding Markov chain by never returning to entries in $\mathfrak X$ already visited and move from $X_t$ to $X_{t+1}$ by choosing among neighbouring entries with probabilities proportional to $p(x)O(x)$ or $\exp\{\alpha p(x)O(x)\}$. With the added perk that since $p(x)O(x)$ is computed for these neighbours, they can all be added to the approximation of $\mathfrak I$ and excluded from future steps. The algorithm could stop when the probability of the visited values is close enough to one.

Here is a toy illustration

where $\mathfrak X$ is an $N\times N$ grid, $O(\cdot)$ is a discretised Normal density (represented by the level set on the above picture) and $\mathfrak X$ is explored by a random exploration until the accumulated probability mass is $0.999$. The white dots in the above picture correspond to the points $x$ that have not been explored.

Here is the attached R code:

#preliminaries
N=5e2;N2=N*N
Y=matrix(0,N,N)#visited points
P=exp(matrix(rnorm(N2),N))
P=P/sum(P)#probability mass function
O=matrix(dnorm((1:N2)%%N,mean=N/2,sd=N/3)+
         dnorm((1:N2)%/%N,mean=N/2,sd=N/3),N)#function O(x)
#stochastic exploration
m=sample(1:N2,1)#current point
t=P[m]*O[m]#targeted expectation
p=P[m]#visited mass
Y[m]=1
while(p<.999){
  b=sample(which(!Y),4)#neighbours
  Y[b]=1
  p=p+sum(P[b])
  t=t+sum(P[b]*O[b])
  m=sample(b,1,prob=P[b]*O[b])}#next value