Simulations based noisy likelihood function I have a problem where I have a measured data vector $D$ with Gaussian uncertainties (covariance matrix $\Sigma$). I am now trying to model this data with a generative model with parameters $\phi$. The problem is that my model is simulation based.
For a given $\phi$, and given random seed s, and given number of particles in the simulation n, I can produce the model prediction $M(\phi, s, n)$ and this model prediction can be compared with the data vector D. The issue here is that the model prediction vector will have some uncertainty (due to randomness of the simulation and finite number of particles). We can assume that this uncertainty can be characterized by a Gaussian distribution with known Covariance $S(\phi, s, n)$
The question is what is the best way to write the likelihood and sample the posterior on $\phi$? Can anyone point at some papers that maybe discuss this kind of problem and how to think about it ? ( I am broadly aware of ABC techniques, but it is not obvious how any of those would be an improvement over what we do now)
For the moment the approach we use is basically ensure that n is large enough (to have smaller simulation errors) and then basically we just ignore the randomness in $M(\phi,s,n)$
and our log-likelihood is just $$ log L(\phi) = -\frac{1}{2}(D-M(\phi))^T \Sigma^{-1} (D-M(\phi))$$ which we sample using MCMC. (the resulting likelihood function is somewhat noisy, but the hope is that it all averages out in the resulting posterior on $\phi$)
Are there any alternatives to this ?

If I understand your problem correctly, you have some data $D$ that is assumed to have a Gaussian distribution with an unknown mean $\mu(\phi)$ and known covariance $\Sigma$, namely
$$ D = \mu(\phi) + \mathcal N (0,\Sigma) $$
As well as a simulation output $M(\phi,n)$ which is also assumed to normally distributed with known covariance $S(\phi,n)$
$$ M(\phi,n) = \mu(\phi) + \mathcal N (0, S(\phi,n)) $$
(Where we can think of $\mu(\phi)$ as being equal to $\mu(\phi) = M(\phi,\infty)$ , since $S(\phi,\infty)=0$ )
But notice that this just implies that
$$ D = M(\phi,n) + \mathcal N (0,\Sigma + S(\phi,n)) $$
So your likelihood should be
$$ \log \mathcal L(\phi) = -\frac{1}{2}(D - M(\phi,n))^T(\Sigma + S(\phi,n))^{-1}(D - M(\phi,n) )$$
Now you can sample it using MCMC as you did, or try to optimize it using something like simulated annealing or Bayesian optimization, or any combination of the above (get close to the global minimum and sample around it).
