Following the talk at MCM 2017 about the strange impact of the random generator on the outcome of an MCMC generator, I tried in Montréal airport the following code on the banana target of Haario et al. (1999), copied from Soetaert and Laine and using the MCMC function of the FME package:
library(FME) Banana <- function (x1, x2) { return(x2 - (x1^2+1)) } pmultinorm <- function(vec, mean, Cov) { diff <- vec - mean ex <- -0.5*t(diff) %*% solve(Cov) %*% diff rdet <- sqrt(det(Cov)) power <- -length(diff)*0.5 return((2.*pi)^power / rdet * exp(ex)) } BananaSS <- function (p) { P <- c(p[1], Banana(p[1], p[2])) Cov <- matrix(nr = 2, data = c(1, 0.9, 0.9, 1)) N=1e3 ejd=matrix(0,4,N) RNGkind("Mars") for (t in 1:N){ MCMC <- modMCMC(f = BananaSS, p = c(0, 0.7), jump = diag(nrow = 2, x = 5), niter = 1e3) ejd[1,t]=mean((MCMC$pars[-1,2]-MCMC$pars[1,2])^2)}
since this divergence from the initial condition seemed to reflect the experiment of the speaker at MCM 2017. Unsurprisingly, no difference came from using the different RNGs in R (which may fail to contain those incriminated by the study)…