get.ppd.2S {BayesPIM}R Documentation

Posterior Predictive Cumulative Incidence Function

Description

Computes the posterior predictive cumulative incidence function (CIF) from a bayes.2S prevalence-incidence mixture model. The function can return quantiles corresponding to user-specified percentiles (i.e., time points at which the cumulative probability reaches certain thresholds), or vice versa (percentiles at user-specified quantiles). Additionally, it allows for marginal or conditional CIFs of either the mixture population (including prevalence as a point-mass at time zero) or the non-prevalent (healthy) subpopulation.

Usage

get.ppd.2S(
  mod,
  fix_Z.X = NULL,
  fix_Z.W = NULL,
  pst.samples = 1000,
  perc = seq(0, 1, 0.01),
  type = "x",
  ppd.type = "percentiles",
  quant = NULL
)

Arguments

mod

A fitted prevalence-incidence mixture model of class bayes.2S.

fix_Z.X

Either NULL for a marginal CIF or a numeric vector of length ncol(Z.X) to request a conditional CIF. Numeric entries fix those covariates at the given value, whereas NA entries are integrated out. See Details.

fix_Z.W

Same as fix_Z.X but for the prevalence model covariates; must be NULL to obtain a marginal CIF.

pst.samples

Integer; number of posterior samples to draw when computing the posterior predictive CIF. Must not exceed the total available posterior samples in mod. Larger values can improve precision but increase computation time.

perc

A numeric vector of cumulative probabilities (i.e., percentiles in (0,1)) for which time points are returned when ppd.type = "quantiles".

type

Character; "xstar" for the mixture CIF (prevalence + incidence), or "x" for the non-prevalent (healthy) population CIF.

ppd.type

Character; "percentiles" to return cumulative probabilities at times quant, or "quantiles" to return time points at cumulative probabilities perc.

quant

A numeric vector of time points for which the function returns cumulative probabilities when ppd.type = "percentiles".

Details

For a prevalence-incidence mixture model, some fraction of the population may already have experienced the event (prevalent cases) at baseline, while the remaining (healthy) fraction has not. This function estimates the CIF in two ways:

You may request a marginal CIF by setting both fix_Z.X = NULL and fix_Z.W = NULL, thus integrating over all covariates. Alternatively, a conditional CIF can be obtained by partially or fully specifying fixed covariate values in fix_Z.X (and optionally fix_Z.W) while integrating out the unspecified covariates (NA entries).

The function operates in two main modes:

Value

A list with some or all of the following elements:

med.cdf
  • If ppd.type = "quantiles", med.cdf contains the median quantiles (time points) across posterior samples for each percentile in perc.

  • If ppd.type = "percentiles", med.cdf contains the median cumulative probabilities across posterior samples for each time point in quant.

med.cdf.ci

A 2-row matrix with the 2.5% and 97.5% posterior quantiles for the estimated med.cdf, reflecting uncertainty.

quant

If ppd.type = "percentiles", this is a copy of the input quant (time points).

perc

If ppd.type = "quantiles", this is a copy of the input perc (cumulative probabilities).

Examples


# Generate data according to the Klausch et al. (2024) PIM
set.seed(2025)
dat <- gen.dat(kappa = 0.7, n = 1e3, theta = 0.2,
               p = 1, p.discrete = 1,
               beta.X = c(0.2, 0.2), beta.W = c(0.2, 0.2),
               v.min = 20, v.max = 30, mean.rc = 80,
               sigma.X = 0.2, mu.X = 5, dist.X = "weibull",
               prob.r  = 1)

# Fit a bayes.2S model (example with moderate ndraws = 2e4)
mod <- bayes.2S(Vobs = dat$Vobs, Z.X = dat$Z, Z.W = dat$Z, r = dat$r,
                kappa = 0.7, update.kappa = FALSE, ndraws = 1e4,
                chains = 2, prop.sd.X = 0.008, parallel = TRUE,
                dist.X = "weibull")

###################
# (1) Provide percentiles, get back quantiles (times)
###################
cif_nonprev <- get.ppd.2S(mod, pst.samples = 1e3, type = "x",
                          ppd.type = "quantiles", perc = seq(0, 1, 0.01))
cif_mix     <- get.ppd.2S(mod, pst.samples = 1e3, type = "xstar",
                          ppd.type = "quantiles", perc = seq(0, 1, 0.01))

# Plot: Non-prevalent stratum CIF vs. mixture CIF (marginal)
oldpar <- par(no.readonly = TRUE)
par(mfrow = c(1,2))

plot(cif_nonprev$med.cdf, cif_nonprev$perc, type = "l", xlim = c(0,300), ylim = c(0,1),
     xlab = "Time", ylab = "Cumulative Incidence")
lines(cif_nonprev$med.cdf.ci[1,], cif_nonprev$perc, lty = 2)
lines(cif_nonprev$med.cdf.ci[2,], cif_nonprev$perc, lty = 2)

plot(cif_mix$med.cdf, cif_mix$perc, type = "l", xlim = c(0,300), ylim = c(0,1),
     xlab = "Time", ylab = "Cumulative Incidence")
lines(cif_mix$med.cdf.ci[1,], cif_mix$perc, lty = 2)
lines(cif_mix$med.cdf.ci[2,], cif_mix$perc, lty = 2)

###################
# (2) Provide quantiles (times), get back percentiles (cumulative probabilities)
###################
cif2_nonprev <- get.ppd.2S(mod, pst.samples = 1e3, type = "x",
                           ppd.type = "percentiles", quant = 1:300)
cif2_mix     <- get.ppd.2S(mod, pst.samples = 1e3, type = "xstar",
                           ppd.type = "percentiles", quant = 1:300)

# Plot: Non-prevalent vs. mixture CIF using times on the x-axis
plot(cif2_nonprev$quant, cif2_nonprev$med.cdf, type = "l", xlim = c(0,300), ylim = c(0,1),
     xlab = "Time", ylab = "Cumulative Incidence")
lines(cif2_nonprev$quant, cif2_nonprev$med.cdf.ci[1,], lty = 2)
lines(cif2_nonprev$quant, cif2_nonprev$med.cdf.ci[2,], lty = 2)

plot(cif2_mix$quant, cif2_mix$med.cdf, type = "l", xlim = c(0,300), ylim = c(0,1),
     xlab = "Time", ylab = "Cumulative Incidence")
lines(cif2_mix$quant, cif2_mix$med.cdf.ci[1,], lty = 2)
lines(cif2_mix$quant, cif2_mix$med.cdf.ci[2,], lty = 2)

###################
# (3) Conditional CIFs by fixing some covariates
###################
cif_mix_m1 <- get.ppd.2S(mod, fix_Z.X = c(-1, NA), pst.samples = 1e3,
                         type = "xstar", ppd.type = "quantiles", perc = seq(0,1,0.01))
cif_mix_0  <- get.ppd.2S(mod, fix_Z.X = c(0, NA),  pst.samples = 1e3,
                         type = "xstar", ppd.type = "quantiles", perc = seq(0,1,0.01))
cif_mix_p1 <- get.ppd.2S(mod, fix_Z.X = c(1, NA),  pst.samples = 1e3,
                         type = "xstar", ppd.type = "quantiles", perc = seq(0,1,0.01))

# Plot: mixture CIF for three different values of the first covariate
par(mfrow = c(1,1))
plot(cif_mix_m1$med.cdf, cif_mix_m1$perc, type = "l", xlim = c(0,300), ylim = c(0,1),
     xlab = "Time", ylab = "Cumulative Incidence", col=1)
lines(cif_mix_0$med.cdf,  cif_mix_m1$perc, col=2)
lines(cif_mix_p1$med.cdf, cif_mix_m1$perc, col=3)

par(oldpar)



[Package BayesPIM version 1.0.0 Index]