| Title: | Bayesian Model Averaging for Univariate Link Latent Gaussian Models |
|---|---|
| Description: | Bayesian model averaging (BMA) algorithms for univariate link latent Gaussian models (ULLGMs). For detailed information, refer to Steel M.F.J. & Zens G. (2024) "Model Uncertainty in Latent Gaussian Models with Univariate Link Function" <doi:10.48550/arXiv.2406.17318>. The package supports various g-priors and a beta-binomial prior on the model space. It also includes auxiliary functions for visualizing and tabulating BMA results. Currently, it offers an out-of-the-box solution for model averaging of Poisson log-normal (PLN) and binomial logistic-normal (BiL) models. The codebase is designed to be easily extendable to other likelihoods, priors, and link functions. |
| Authors: | Gregor Zens [aut, cre], Mark F.J. Steel [aut] |
| Maintainer: | Gregor Zens <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.1.1 |
| Built: | 2024-11-29 08:44:17 UTC |
| Source: | CRAN |
plotBeta produces a visualization of the estimated posterior means of the coefficients, extracted from ULLGM_BMA results.
plotBeta(x, variable_names = NULL, sort = TRUE)plotBeta(x, variable_names = NULL, sort = TRUE)
x |
The output object of |
variable_names |
A character vector specifying the names of the columns of X. |
sort |
Logical, indicating whether the plot should be sorted by posterior mean. Defaults to |
Returns a 'ggplot2::ggplot' object.
Gregor Zens
# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) plotBeta(results_pln)# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) plotBeta(results_pln)
plotModelSize produces a visualization of the posterior distribution of model size, extracted from ULLGM_BMA results.
plotModelSize(x)plotModelSize(x)
x |
The output object of |
Returns a 'ggplot2::ggplot' object visualizing the posterior distribution of model size.
Gregor Zens
# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) plotModelSize(results_pln)# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) plotModelSize(results_pln)
plotPIP produces a visualization of the posterior inclusion probabilities (PIPs) extracted from ULLGM_BMA results.
plotPIP(x, variable_names = NULL, sort = TRUE)plotPIP(x, variable_names = NULL, sort = TRUE)
x |
The output object of |
variable_names |
A character vector specifying the names of the columns of X. |
sort |
Logical, indicating whether the plot should be sorted by PIP. Defaults to |
Returns a 'ggplot2::ggplot' object.
Gregor Zens
# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) plotPIP(results_pln)# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) plotPIP(results_pln)
ULLGM_BMA Estimation ResultssummaryBMA produces a table with estimated posterior means, standard deviations, and posterior inclusion probabilities (PIPs) for the results of a ULLGM_BMA estimation.
summaryBMA(x, variable_names = NULL, digits = 3, sort = FALSE, type = "pandoc")summaryBMA(x, variable_names = NULL, digits = 3, sort = FALSE, type = "pandoc")
x |
The output object of |
variable_names |
A character vector specifying the names of the columns of X. |
digits |
Number of digits to round the table to. Defaults to 3. |
sort |
Logical, indicating whether the table should be sorted by PIPs. Default is |
type |
A character string indicating the format of the table. Options are |
Returns a 'knitr::kable' object containing the summary table.
Gregor Zens
# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) summaryBMA(results_pln)# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) summaryBMA(results_pln)
ULLGM_BMA Estimation ResultstopModels produces a table of the top n models from a ULLGM_BMA object, sorted by posterior model probabilities.
topModels(x, variable_names = NULL, type = "pandoc", digits = 3, n = 5)topModels(x, variable_names = NULL, type = "pandoc", digits = 3, n = 5)
x |
The output object of |
variable_names |
A character vector specifying the names of the columns of X. |
type |
A character string indicating the format of the table. Options are |
digits |
Number of digits to round the table to. Defaults to 3. |
n |
Number of top models to be returned. Defaults to 5. |
Returns a 'knitr::kable' object containing the table of top models.
Gregor Zens
# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) # Top 5 models topModels(results_pln)# Load package library(LatentBMA) # Example: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) # Top 5 models topModels(results_pln)
tracePlot produces traceplots for selected parameters, extracted from ULLGM_BMA results.
tracePlot(x, parameter = "beta", index = 1)tracePlot(x, parameter = "beta", index = 1)
x |
The output object of |
parameter |
Specifies which parameter should be considered for the traceplot. Options are |
index |
If |
Returns a 'ggplot2::ggplot' object.
Gregor Zens
ULLGM_BMA estimates Bayesian regression models using either a Poisson log-normal (PLN) or binomial logistic-normal (BiL) regression framework. It accounts for model uncertainty via Bayesian model averaging.
ULLGM_BMA(X, y, model = "PLN", gprior = "BRIC", nsave = 10000, nburn = 2000, Ni = NULL, m = NULL, verbose = TRUE)ULLGM_BMA(X, y, model = "PLN", gprior = "BRIC", nsave = 10000, nburn = 2000, Ni = NULL, m = NULL, verbose = TRUE)
X |
A n x p design matrix where n is the number of observations and p is the number of explanatory variables. |
y |
A n x 1 response vector. For PLN and BiL models, this is a count response. |
model |
Indicates the model to be estimated. Options are |
gprior |
Specifies the g-prior to be used. Options under fixed g are |
nsave |
The number of saved posterior samples. Defaults to 10,000. |
nburn |
The number of initial burn-in samples. Defaults to 2,000. |
Ni |
A vector containing the number of trials for each observation when estimating a binomial logistic-normal model. Required if |
m |
The prior expected model size as per the beta-binomial prior of Ley and Steel (2009). Defaults to |
verbose |
Logical indicator of whether progress should be printed during estimation. Default is |
A list containing the inputs and selected posterior simulation outputs, such as posterior chains for the coefficients and inclusion vectors.
All explanatory variables in X are automatically demeaned within the function. All models do automatically include an intercept term.
Liang, F., Paulo, R., Molina, G., Clyde, M. A., & Berger, J. O. (2008). Mixtures of g priors for Bayesian variable selection. Journal of the American Statistical Association, 103(481), 410-423.
Zellner, A., & Siow, A. (1980). Posterior odds ratios for selected regression hypotheses. Trabajos de estadÃstica y de investigación operativa, 31, 585-603.
Ley, E., & Steel, M. F. J. (2009). On the effect of prior assumptions in Bayesian model averaging with applications to growth regression. Journal of Applied Econometrics, 24(4), 651-674.
# Load package library(LatentBMA) # Example 1: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) # Example 2: Estimate a BiL model under a Zellner-Siow prior with m = 5 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) Ni <- rep(50, 100) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rbinom(100, Ni, 1 / (1 + exp(-z))) results_bil <- ULLGM_BMA(X = X, y = y, Ni = Ni, model = "BiL", nsave = 250, nburn = 250, m = 5, gprior = "zellnersiow")# Load package library(LatentBMA) # Example 1: Estimate a PLN model under a BRIC prior with m = p/2 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rpois(100, exp(z)) results_pln <- ULLGM_BMA(X = X, y = y, model = "PLN", nsave = 250, nburn = 250) # Example 2: Estimate a BiL model under a Zellner-Siow prior with m = 5 using simulated data # Note: Use more samples for actual analysis # Note: nsave = 250 and nburn = 250 are for demonstration purposes X <- matrix(rnorm(100*20), 100, 20) Ni <- rep(50, 100) z <- 2 + X %*% c(0.5, -0.5, rep(0, 18)) + rnorm(100, 0, sqrt(0.25)) y <- rbinom(100, Ni, 1 / (1 + exp(-z))) results_bil <- ULLGM_BMA(X = X, y = y, Ni = Ni, model = "BiL", nsave = 250, nburn = 250, m = 5, gprior = "zellnersiow")