IBIS.logreg {UNCOVER}

Logistic regression iterated batch importance sampling

Description

This function uses an Iterated Batch Importance Sampling (IBIS) scheme with batch size one to go from prior to full posterior. We assume a Bayesian logistic regression model.

Usage

IBIS.logreg(
  X,
  y,
  options = IBIS.logreg.opts(),
  prior_mean = rep(0, ncol(X) + 1),
  prior_var = diag(ncol(X) + 1)
)

Arguments

X

Co-variate matrix

y

Binary response vector

options

Additional arguments that can be specified for IBIS.logreg. See IBIS.logreg.opts() for details. Can be ignored.

prior_mean

Mean for the multivariate normal prior used in the SMC sampler. See details. Defaults to the origin.

prior_var

Variance matrix for the multivariate normal prior used in the SMC sampler. See details. Defaults to the identity matrix.

Details

Details of the internal mechanisms of the SMC sampler such as the Metropolis-Hastings MCMC resample move can be found in Emerson and Aslett (2023) and Chopin (2002).

It is never recommended to use anything other than IBIS.logreg.opts to provide the options argument. See examples and IBIS.logreg.opts() for more information.

The prior used for the IBIS procedure will take the form of a multivariate normal, where the parameters can be specified directly by the user. It is however possible to override this default prior distributional form by specifying prior.override=TRUE and providing the relevant prior functions in IBIS.logreg.opts.

Value

An object of class "IBIS", which is a list consisting of:

covariate_matrix

The co-variate matrix provided.

response_vector

The binary response vector provided.

samples

A matrix of samples from the posterior.

log_Bayesian_evidence

An estimate of the log Bayesian evidence (or normalisation constant) of the posterior.

diagnostics

A data frame recording the features of the SMC sampler as the observations were added.

If weighted==TRUE then an additional element of the list (weights) is added detailing the weights of the posterior samples.

References

• Emerson, S.R. and Aslett, L.J.M. (2023). Joint cohort and prediction modelling through graphical structure analysis (to be released)

• Chopin, N. (2002). A sequential particle filter method for static models. Biometrika, 89(3), 539-552, doi:10.1093/biomet/89.3.539

See Also

IBIS.logreg.opts(), print.IBIS(), predict.IBIS(), plot.IBIS()

Examples

require(graphics)
# First we generate a co-variate matrix X and binary response vector y
CM <- matrix(rnorm(200),100,2)
rv <- sample(0:1,100,replace=TRUE)

# Now we can obtain 1000 samples from the posterior from a standard
# multivariate normal prior
out.1 <- IBIS.logreg(X = CM,y = rv)
plot(out.1)
out.1$log_Bayesian_evidence

# We can specify that the samples be weighted
out.1.w <- IBIS.logreg(X = CM,y = rv,
                       options = IBIS.logreg.opts(weighted = TRUE))
out.1.w$weights
plot(out.1.w)

# We can also specify different arguments for a specific prior
out.2 <- IBIS.logreg(X = CM,y = rv,prior_mean = rep(-3,3),
                     prior_var = 0.1*diag(3))
samp.df <- data.frame(rbind(out.1$samples,out.2$samples))
colnames(samp.df) <- paste0("beta[",c(0:2),"]")
GGally::ggpairs(samp.df,
                labeller = "label_parsed",
                ggplot2::aes(color = as.factor(rep(c(1,2),each=1000))),
                upper = list(continuous = GGally::wrap("density")),
                lower = list(continuous = GGally::wrap("points",size=0.5)))
out.2$log_Bayesian_evidence
out.3 <- IBIS.logreg(X = CM,y = rv,prior_mean = rep(3,3),
                     prior_var = 0.1*diag(3))
samp.df <- data.frame(rbind(out.1$samples,out.2$samples,out.3$samples))
colnames(samp.df) <- paste0("beta[",c(0:2),"]")
GGally::ggpairs(samp.df,
                labeller = "label_parsed",
                ggplot2::aes(color = as.factor(rep(c(1,2,3),each=1000))),
                upper = list(continuous = GGally::wrap("density")),
                lower = list(continuous = GGally::wrap("points",size=0.5)))
out.3$log_Bayesian_evidence

# We can also change the prior, for example a multivariate independent
# uniform
rmviu <- function(n,a,b){
return(mapply(FUN = function(min.vec,max.vec,pn){stats::runif(pn,a,b)},
              min.vec=a,max.vec=b,MoreArgs = list(pn = n)))
}
dmviu <- function(x,a,b){
for(ii in 1:ncol(x)){
  x[,ii] <- dunif(x[,ii],a[ii],b[ii])
}
return(apply(x,1,prod))
}

out.4 <- IBIS.logreg(X = CM,y = rv,
                     options = IBIS.logreg.opts(prior.override = TRUE,
                                                rprior = rmviu,
                                                dprior = dmviu,a=rep(0,3),
                                                b=rep(1,3)))
samp.df <- data.frame(rbind(out.1$samples,out.4$samples))
colnames(samp.df) <- paste0("beta[",c(0:2),"]")
GGally::ggpairs(samp.df,
                labeller = "label_parsed",
                ggplot2::aes(color = as.factor(rep(c(1,4),each=1000))),
                upper = list(continuous = GGally::wrap("points",size=0.5)),
                lower = list(continuous = GGally::wrap("points",size=0.5)))
out.4$log_Bayesian_evidence