checkConsistency error

[Jan-Ole Koslik]

Hi,

running

library(RTMB)
data(InsectSprays)

X <- cbind(1, model.matrix(~ spray - 1, data = InsectSprays))

par <- list(
  beta0 = log(mean(InsectSprays$count)),
  beta = rep(0, length(levels(InsectSprays$spray))),
  log_sigma = log(1)
)
dat <- list(count = InsectSprays$count, X = X)

nll <- function(par) {
  getAll(par, dat)
  count <- OBS(count)
  lambda <- exp(X %*% c(beta0, beta)); ADREPORT(lambda)
  nll <- -sum(dpois(count, lambda, log = TRUE))
  sigma <- exp(log_sigma); ADREPORT(sigma)
  nll - sum(dnorm(beta, 0, sigma, log = TRUE))
}

obj <- MakeADFun(nll, par, random = "beta")
opt <- nlminb(obj$par, obj$fn, obj$gr)

checkConsistency(obj)

gives error 

Error in X %*% c(beta0, beta) :
  requires numeric/complex matrix/vector arguments

Whats going wrong in the automatic simulation here?

Thanks :)

[Kasper Kristensen]

Automatic simulation works under 'some restrictions' specified under ?Simulation (and you're not the first to be confused by this - see https://github.com/kaskr/RTMB/issues/56)

In your case, reordering the likelihood accumulation should fix it:

nll <- function(par) {
  getAll(par, dat)
  count <- OBS(count)

  sigma <- exp(log_sigma); ADREPORT(sigma)
  nll <- - sum(dnorm(beta, 0, sigma, log = TRUE))

  lambda <- exp(X %*% c(beta0, beta)); ADREPORT(lambda)

  nll <- nll - sum(dpois(count, lambda, log = TRUE))
  nll
}

[Jan-Ole Koslik]

Thanks Kasper, does that mean in general it's recommended to first evaluate random effect density and then the observation density?

[Kasper Kristensen]

Yes, density accumulation should be in 'topological order' just like you would order a simulation procedure.