I’m implementing a brms
custom family and would like to add a log_lik
method as described in the “Post-Processing Custom Family Models” section of the vignette “Define Custom Response Distributions with brms”.
The challenge I am having is that the log likelihood \log p(y | \theta) depends on some “latent variables” (not only the distributional parameters of the model, which can be obtained via brms::get_dpar
). I am unsure as to how I can extract these latent variable parameters from the model in order to compute the log likelihood.
Here is what I have written for the log_lik_latent_lognormal
function so far (where latent_lognormal
is the name of the family):
#' Calculate the pointwise log likelihood
#'
#' See [brms::log_lik()].
#'
#' @param i The index of the observation to calculate the log likelihood of
#' @param prep The result of a call to [`brms::prepare_predictions`]
#' @family postprocess
#' @autoglobal
#' @export
log_lik_latent_lognormal <- function(i, prep) {
mu <- brms::get_dpar(prep, "mu", i = i)
sigma <- brms::get_dpar(prep, "sigma", i = i)
y <- prep$data$Y[i]
obs_t <- prep$data$vreal1[i]
pwindow_width <- prep$data$vreal2[i]
swindow_width <- prep$data$vreal3[i]
# Generate values of the swindow_raw and pwindow_raw
# Really these should be extracted from prep
# How to do this?
swindow_raw <- runif(prep$ndraws)
pwindow_raw <- runif(prep$ndraws)
pwindow <- pwindow_raw * pwindow_width
swindow <- swindow_raw * swindow_width
d <- y - pwindow + swindow
obs_time <- obs_t - pwindow
lpdf <- dlnorm(d, meanlog = mu, sdlog = sigma, log = TRUE)
lcdf <- plnorm(obs_time, meanlog = mu, sdlog = sigma, log.p = TRUE)
return(lpdf - lcdf)
}
The section with the comment is where I’d like to be getting the pwindow_raw
and swindow_raw
draws from the fitted model (via prep
somehow).
I’ve looked a bit a the brms
code (see issue for a log of my thinking [including at one point I think wrongly thinking about integrating these variables out to get a marginal likelihood]). Anyway, in brms:::prepare_predictions.bframel
there are these prepare_predictions_xx
functions for a range of xx
. Are these like special cases where more than the distributional parameters are required? Would my situation be like one of these special cases?
out$fe <- prepare_predictions_fe(x, draws, sdata, ...)
out$sp <- prepare_predictions_sp(x, draws, sdata, ...)
out$cs <- prepare_predictions_cs(x, draws, sdata, ...)
out$sm <- prepare_predictions_sm(x, draws, sdata, ...)
out$gp <- prepare_predictions_gp(x, draws, sdata, ...)
out$re <- prepare_predictions_re(x, sdata, ...)
out$ac <- prepare_predictions_ac(x, draws, sdata, nat_cov = FALSE,
Another idea: is it possible it put within ...
the fit
object, and then rely on the user providing this to any functions requiring the log_lik
?
Thank you! Let me know if I could make any part of this more clear.
Session information
- Operating System: Mac OS
brms
version: 2.21.6