pwe.logml.post: Log marginal likelihood of a piecewise exponential (PWE)...

In hdbayes: Bayesian Analysis of Generalized Linear Models with Historical Data

Log marginal likelihood of a piecewise exponential (PWE) model under a normal/half-normal prior

Description

Uses bridge sampling to estimate the logarithm of the marginal likelihood of a PWE model under the normal/half-normal prior.

Usage

pwe.logml.post(post.samples, bridge.args = NULL)

Arguments

post.samples	output from pwe.post() giving posterior samples of a PWE model under the normal/half-normal prior, with an attribute called 'data' which includes the list of variables specified in the data block of the Stan program.
bridge.args	a list giving arguments (other than samples, log_posterior, data, lb, and ub) to pass onto bridgesampling::bridge_sampler().

Value

The function returns a list with the following objects

model

"pwe_post"

logml

the estimated logarithm of the marginal likelihood

bs

an object of class bridge or bridge_list containing the output from using bridgesampling::bridge_sampler() to compute the logarithm of the marginal likelihood of the PWE model under the normal/half-normal prior

References

Gronau, Q. F., Singmann, H., and Wagenmakers, E.-J. (2020). bridgesampling: An r package for estimating normalizing constants. Journal of Statistical Software, 92(10).

Examples

if (instantiate::stan_cmdstan_exists()) {
  if(requireNamespace("survival")){
    library(survival)
    data(E1690)
    ## take subset for speed purposes
    E1690 = E1690[1:100, ]
    ## replace 0 failure times with 0.50 days
    E1690$failtime[E1690$failtime == 0] = 0.50/365.25
    E1690$cage = as.numeric(scale(E1690$age))
    data_list = list(currdata = E1690)
    nbreaks = 3
    probs   = 1:nbreaks / nbreaks
    breaks  = as.numeric(
      quantile(E1690[E1690$failcens==1, ]$failtime, probs = probs)
    )
    breaks  = c(0, breaks)
    breaks[length(breaks)] = max(10000, 1000 * breaks[length(breaks)])
    d.post = pwe.post(
      formula = survival::Surv(failtime, failcens) ~ treatment + sex + cage + node_bin,
      data.list = data_list,
      breaks = breaks,
      chains = 1, iter_warmup = 500, iter_sampling = 1000
    )
    pwe.logml.post(
      post.samples = d.post,
      bridge.args = list(silent = TRUE)
    )
  }
}

hdbayes documentation built on Nov. 21, 2025, 1:07 a.m.