Apologies, I feel that this is a fairly basic question, however, after searching high and low I can find no answer. There is a similar post here, however, it doesn’t seem to answer my question.
This stems from trying to find bayesian point estimates from a poisson regression model where some of the data can be quite sparse. Note that I need to both a point estimate and a credible interval and the relationship between these is not as expected.
I have isolated my difficulty in understanding what is going on to the following;
Given a simple poisson model with improper prior in stan as follows:
poisson_model <- "data{
int<lower=1> n;
int<lower=0> x[n];
}
parameters{
real<lower=0> lambda;
}
model{
x ~ poisson(lambda);
}"
library(rstan)
set.seed(1)
n <- 2
lambda <- 10
x <- rpois(n , lambda)
mod <- stan_model(model_code = poisson_model)
map <- optimizing(mod, data = list(n = n, x = x))
fit_mod <- sampling(mod, data = list(n = n, x = x), iter = 12000, warmup = 2000, chains = 1)
# using just 1 chain as don't have any issues with convergence in such a simple model
map$par
summary(fit_mod)
mcmc <- extract(fit_mod)
mcmc$lambda[which.max(mcmc$lp__)]
I can use this to find the MLE (or MAP if I include a prior) estimator using the optimizing function in rstan. This comes out as expected and is equal to sum(x)/n or 9 in the example above.
However, sampling from the model using the sampling function and looking at the value for lambda that maximizes lp__ this turns out to be equal to (sum(x) + 1)/n or 9.5 in the example above.
My expectation was that the maximum value of lp__ would be at the mode of the posterior of lambda and therefore the MLE of lambda. It turns out that lp__ is maximum at the mean of the posterior of lambda. I am clearly missing something basic.
A further observation is that the posterior distribution of lambda is a Gamma(sum(x) + 1, n) distribution which is equivalent to the prior being a Gamma(alpha, beta) distribution with alpha = 1 and beta approaching 0.
I am happy to do the hard reading if there is a suitable reference or link to a response.