Cannot reproduce rstan::optimizing() and rstan::sampling() results with the same seed for my paper

I am writing a paper to accompany an R package that I have written, which calls rstan v2.21.3. I wish to replicate some analysis I did last year and have it as an example in the README file on my Github repo.

So, a few months ago I had used rstan v2.21.3 to draw from the posterior of a complex Bayesian hierarchical model that I had defined in my own .stan program file. I executed rstan::optimizing() and the estimates were then used as initialisation points for the rstan::sampling() routine. I was happy with the convergence diagnostics, posterior outputs, etc.

I have now created an R package that calls rstan v2.21.3 and wish to replicate the analysis.

  • I have ensured that the data, prior specification and seeds in rstan::optimizing() and rstan::sampling() are identical to my MCMC experiment a few months ago.
  • I have also thoroughly checked that the user-defined options are corrently passed to the .stan program.
  • I have ensured that the version of StanHeaders is the same.

Having executed `rstan::optimizing(), this time I receive the error :

Error in chol.default(-H) : 
  the leading minor of order 633 is not positive definite

Additionally, the estimates from rstan::optimizing() called by my package do not agree with the respective estimates.

Further, the computationally intensive rstan::sampling() routine for my example (needs approx 11 hrs to complete) now fails to converge.

Q: is there something happening under the hood with the related packages I am calling like Rcpp, RcppParallel, RcppEigen, etc, that do not allow reproducibility of results? Despite my checks, I really cannot understand how to address this.

Lampros Bouranis