Nothing
#' @title Prints a fitted RoBMA object
#'
#' @param x a fitted RoBMA object.
#' @param ... additional arguments.
#'
#'
#' @return \code{print.RoBMA} invisibly returns the print statement.
#'
#' @seealso [RoBMA()]
#' @export
print.RoBMA <- function(x, ...){
cat("Call:\n")
print(x$call)
cat("\nEstimates:\n")
print(stats::coef(x))
}
#' @title Summarize fitted RoBMA object
#'
#' @description \code{summary.RoBMA} creates summary tables for a
#' RoBMA object.
#'
#' @param object a fitted RoBMA object
#' @param type whether to show the overall RoBMA results (\code{"ensemble"}),
#' an overview of the individual models (\code{"models"}), an overview of
#' the individual models MCMC diagnostics (\code{"diagnostics"}), or a detailed summary
#' of the individual models (\code{"individual"}). Can be abbreviated to first letters.
#' @param conditional show the conditional estimates (assuming that the
#' alternative is true). Defaults to \code{FALSE}. Only available for
#' \code{type == "ensemble"}.
#' @param output_scale transform the meta-analytic estimates to a different
#' scale. Defaults to \code{NULL} which returns the same scale as the model was estimated on.
#' @param probs quantiles of the posterior samples to be displayed.
#' Defaults to \code{c(.025, .975)}
#' @param logBF show log of Bayes factors. Defaults to \code{FALSE}.
#' @param BF01 show Bayes factors in support of the null hypotheses. Defaults to
#' \code{FALSE}.
#' @param short_name whether priors names should be shortened to the first
#' (couple) of letters. Defaults to \code{FALSE}.
#' @param remove_spike_0 whether spike prior distributions with location at zero should
#' be omitted from the summary. Defaults to \code{FALSE}.
#' @param ... additional arguments
#'
#' @examples \dontrun{
#' # using the example data from Anderson et al. 2010 and fitting the default model
#' # (note that the model can take a while to fit)
#' fit <- RoBMA(r = Anderson2010$r, n = Anderson2010$n, study_names = Anderson2010$labels)
#'
#' # summary can provide many details about the model
#' summary(fit)
#'
#' # estimates from the conditional models can be obtained with
#' summary(fit, conditional = TRUE)
#'
#' # overview of the models and their prior and posterior probability, marginal likelihood,
#' # and inclusion Bayes factor can be obtained with
#' summary(fit, type = "models")
#'
#' # diagnostics overview, containing the maximum R-hat, minimum ESS, maximum MCMC error, and
#' # maximum MCMC error / sd across parameters for each individual model can be obtained with
#' summary(fit, type = "diagnostics")
#'
#' # summary of individual models and their parameters can be further obtained by
#' summary(fit, type = "individual")
#' }
#'
#' @note See [diagnostics()] for visual convergence checks of the individual models.
#'
#'
#' @return \code{summary.RoBMA} returns a list of tables of class 'BayesTools_table'.
#'
#' @seealso [RoBMA()], [diagnostics()], [check_RoBMA()]
#' @export
summary.RoBMA <- function(object, type = "ensemble", conditional = FALSE,
output_scale = NULL, probs = c(.025, .975), logBF = FALSE, BF01 = FALSE,
short_name = FALSE, remove_spike_0 = FALSE, ...){
BayesTools::check_bool(conditional, "conditional")
BayesTools::check_char(type, "type")
BayesTools::check_char(output_scale, "output_scale", allow_NULL = TRUE)
BayesTools::check_real(probs, "probs", allow_NULL = TRUE, check_length = 0)
BayesTools::check_bool(BF01, "BF01")
BayesTools::check_bool(logBF, "logBF")
BayesTools::check_bool(short_name, "short_name")
BayesTools::check_bool(remove_spike_0, "remove_spike_0")
# apply version changes to RoBMA object
object <- .update_object(object)
# check the scales
if(is.null(output_scale)){
output_scale <- object$add_info[["output_scale"]]
}else if(object$add_info[["output_scale"]] == "y" & .transformation_var(output_scale, estimation = FALSE) != "y"){
stop("Models estimated using the generall effect size scale 'y' / 'none' cannot be transformed to a different effect size scale.")
}else{
output_scale <- .transformation_var(output_scale, estimation = FALSE)
}
# print diagnostics if all models fail to converge
if(!any(.get_model_convergence(object))){
if(substr(type,1,1) != "d")
warning("All models failed to converge. Model diagnostics were printed instead.")
type <- "diagnostics"
}
if(substr(type,1,1) == "e"){
# transform the estimates if needed
if(object$add_info[["output_scale"]] != output_scale){
object <- .transform_posterior(object, object$add_info[["output_scale"]], output_scale)
}
# obtain components overview
components <- BayesTools::ensemble_inference_table(
inference = object$RoBMA[["inference"]],
parameters = names(object$RoBMA[["inference"]])[names(object$RoBMA[["inference"]]) %in% c("Effect", "Heterogeneity", "Bias", "Hierarchical", "Baseline")],
logBF = logBF,
BF01 = BF01,
title = "Components summary:"
)
# obtain estimates tables
estimates <- BayesTools::ensemble_estimates_table(
samples = object$RoBMA[["posteriors"]],
parameters = names(object$RoBMA[["posteriors"]]),
probs = probs,
title = "Model-averaged estimates:",
footnotes = c(.scale_note(object$add_info[["prior_scale"]], output_scale), .note_omega(object)),
warnings = .collect_errors_and_warnings(object)
)
# deal with possibly empty table in case of no alternative models
if(is.null(object$RoBMA[["posteriors_conditional"]])){
estimates_conditional <- data.frame(matrix(nrow = 0, ncol = length(probs) + 2))
colnames(estimates_conditional) <- c("Mean", "Median", probs)
class(estimates_conditional) <- c("BayesTools_table", "BayesTools_ensemble_summary", class(estimates_conditional))
attr(estimates_conditional, "type") <- rep("estimate", ncol(estimates_conditional))
attr(estimates_conditional, "rownames") <- TRUE
attr(estimates_conditional, "title") <- "Conditional estimates:"
attr(estimates_conditional, "footnotes") <- c(.scale_note(object$add_info[["prior_scale"]], output_scale), .note_omega(object))
attr(estimates_conditional, "warnings") <- .collect_errors_and_warnings(object)
}else{
estimates_conditional <- BayesTools::ensemble_estimates_table(
samples = object$RoBMA[["posteriors_conditional"]],
parameters = names(object$RoBMA[["posteriors_conditional"]]),
probs = probs,
title = "Conditional estimates:",
footnotes = c(.scale_note(object$add_info[["prior_scale"]], output_scale), .note_omega(object)),
warnings = .collect_errors_and_warnings(object)
)
}
# create the output object
output <- list(
call = object[["call"]],
title = .object_title(object),
components = components,
estimates = estimates
)
if(conditional){
output$estimates_conditional <- estimates_conditional
}
# add meta-regression summaries
if(is.RoBMA.reg(object)){
# rename the inference components
for(i in seq_along(object$RoBMA[["inference_predictors"]])){
attr(object$RoBMA[["inference_predictors"]][[i]], "parameter_name") <- gsub("(mu) ", "", attr(object$RoBMA[["inference_predictors"]][[i]], "parameter_name"), fixed = TRUE)
}
if(!is.null(object$RoBMA[["inference_predictors"]])){
output$components_predictors <- BayesTools::ensemble_inference_table(
inference = object$RoBMA[["inference_predictors"]],
parameters = names(object$RoBMA[["inference_predictors"]]),
logBF = logBF,
BF01 = BF01,
title = "Meta-regression components summary:"
)
}
if(!is.null(object$RoBMA[["posteriors_predictors"]])){
# obtain estimates tables
output$estimates_predictors <- BayesTools::ensemble_estimates_table(
samples = object$RoBMA[["posteriors_predictors"]],
parameters = names(object$RoBMA[["posteriors_predictors"]]),
probs = probs,
title = "Model-averaged meta-regression estimates:",
formula_prefix = FALSE,
footnotes = .scale_note(object$add_info[["prior_scale"]], output_scale),
warnings = .collect_errors_and_warnings(object)
)
}
# deal with possibly empty table in case of no alternative models
if(conditional){
if(is.null(object$RoBMA[["posteriors_predictors_conditional"]])){
estimates_predictors_conditional <- data.frame(matrix(nrow = 0, ncol = length(probs) + 2))
colnames(estimates_predictors_conditional) <- c("Mean", "Median", probs)
class(estimates_predictors_conditional) <- c("BayesTools_table", "BayesTools_ensemble_summary", class(estimates_predictors_conditional))
attr(estimates_predictors_conditional, "type") <- rep("estimate", ncol(estimates_predictors_conditional))
attr(estimates_predictors_conditional, "rownames") <- TRUE
attr(estimates_predictors_conditional, "title") <- "Conditional meta-regression estimates:"
attr(estimates_predictors_conditional, "footnotes") <- .scale_note(object$add_info[["prior_scale"]], output_scale)
attr(estimates_predictors_conditional, "warnings") <- .collect_errors_and_warnings(object)
}else{
estimates_predictors_conditional <- BayesTools::ensemble_estimates_table(
samples = object$RoBMA[["posteriors_predictors_conditional"]],
parameters = names(object$RoBMA[["posteriors_predictors_conditional"]]),
probs = probs,
title = "Conditional meta-regression estimates:",
formula_prefix = FALSE,
footnotes = .scale_note(object$add_info[["prior_scale"]], output_scale),
warnings = .collect_errors_and_warnings(object)
)
}
output$estimates_predictors_conditional <- estimates_predictors_conditional
}
}
class(output) <- "summary.RoBMA"
attr(output, "type") <- "ensemble"
return(output)
}else if(substr(type,1,1) == "m"){
components <- names(object$RoBMA[["inference"]])[names(object$RoBMA[["inference"]]) %in% c("Effect", "Heterogeneity", "Bias", "Hierarchical", "Baseline")]
parameters <- list()
if(any(components == "Effect")){
if(.is_regression(object)){
parameters[["intercept"]] <- "mu_intercept"
for(i in seq_along(object$add_info[["predictors"]])){
parameters[[object$add_info[["predictors"]][i]]] <- .BayesTools_parameter_name(object$add_info[["predictors"]][i])
}
}else{
parameters[["Effect"]] <- "mu"
}
}
if(any(components == "Heterogeneity")){
parameters[["Heterogeneity"]] <- "tau"
if(.is_multivariate(object)){
parameters[["Hierarchical"]] <- "rho"
}
}
if(any(components == "Bias")){
parameters[["Bias"]] <- c("PET", "PEESE", "omega")
}
if(is.BiBMA(object)){
parameters[["Baseline"]] <- "pi"
}
summary <- BayesTools::ensemble_summary_table(
models = object[["models"]],
parameters = parameters,
title = "Models overview:",
footnotes = NULL,
warnings = .collect_errors_and_warnings(object),
short_name = short_name,
remove_spike_0 = remove_spike_0
)
output <- list(
call = object[["call"]],
title = .object_title(object),
summary = summary
)
class(output) <- "summary.RoBMA"
attr(output, "type") <- "models"
return(output)
}else if(substr(type,1,1) == "d"){
components <- names(object$RoBMA[["inference"]])[names(object$RoBMA[["inference"]]) %in% c("Effect", "Heterogeneity", "Bias", "Hierarchical", "Baseline")]
parameters <- list()
if(any(components == "Effect")){
if(.is_regression(object)){
parameters[["intercept"]] <- "mu_intercept"
for(i in seq_along(object$add_info[["predictors"]])){
parameters[[object$add_info[["predictors"]][i]]] <- .BayesTools_parameter_name(object$add_info[["predictors"]][i])
}
}else{
parameters[["Effect"]] <- "mu"
}
}
if(any(components == "Heterogeneity")){
parameters[["Heterogeneity"]] <- "tau"
if(.is_multivariate(object)){
parameters[["Hierarchical"]] <- "rho"
}
}
if(any(components == "Bias")){
parameters[["Bias"]] <- c("PET", "PEESE", "omega")
}
if(is.BiBMA(object)){
parameters[["Baseline"]] <- "pi"
}
diagnostics <- BayesTools::ensemble_diagnostics_table(
models = object[["models"]],
parameters = parameters,
title = "Diagnostics overview:",
footnotes = NULL,
warnings = .collect_errors_and_warnings(object),
short_name = short_name,
remove_spike_0 = remove_spike_0
)
output <- list(
call = object[["call"]],
title = .object_title(object),
diagnostics = diagnostics
)
class(output) <- "summary.RoBMA"
attr(output, "type") <- "diagnostics"
return(output)
}else if(substr(type, 1, 1) == "i"){
output <- list(
call = object[["call"]],
title = .object_title(object),
models = list()
)
for(i in seq_along(object[["models"]])){
summary <- BayesTools::model_summary_table(
model = object[["models"]][[i]],
remove_parameters = "gamma",
short_name = short_name,
remove_spike_0 = remove_spike_0
)
if(output_scale == "y"){
estimates <- object[["models"]][[i]][["fit_summary"]]
attr(estimates, "warnings") <- object[["models"]][[i]][["warnings"]]
attr(estimates, "title") <- "Parameter estimates:"
}else{
estimates <- object[["models"]][[i]][["fit_summaries"]][[output_scale]]
attr(estimates, "footnotes") <- .scale_note(object[["models"]][[i]][["prior_scale"]], output_scale)
attr(estimates, "warnings") <- object[["models"]][[i]][["warnings"]]
attr(estimates, "title") <- "Parameter estimates:"
}
output[["models"]][[i]] <- list(
summary = summary,
estimates = estimates
)
}
class(output) <- "summary.RoBMA"
attr(output, "type") <- "individual"
return(output)
}else{
stop(paste0("Unknown summary type: '", type, "'."))
}
}
#' @title Prints summary object for RoBMA method
#'
#' @param x a summary of a RoBMA object
#' @param ... additional arguments
#'
#'
#' @return \code{print.summary.RoBMA} invisibly returns the print statement.
#'
#' @seealso [RoBMA()]
#' @export
print.summary.RoBMA <- function(x, ...){
cat("Call:\n")
print(x[["call"]])
cat("\n")
cat(x[["title"]])
if(attr(x, "type") == "ensemble"){
cat("\n")
print(x[["components"]])
if(!is.null(x[["components_predictors"]])){
cat("\n")
print(x[["components_predictors"]])
}
cat("\n")
print(x[["estimates"]])
if(!is.null(x[["estimates_predictors"]])){
cat("\n")
print(x[["estimates_predictors"]])
}
if(!is.null(x[["estimates_conditional"]])){
cat("\n")
print(x[["estimates_conditional"]])
}
if(!is.null(x[["estimates_predictors_conditional"]])){
cat("\n")
print(x[["estimates_predictors_conditional"]])
}
return(invisible())
}else if(attr(x, "type") == "models"){
cat("\n")
print(x[["summary"]])
return(invisible())
}else if(attr(x, "type") == "diagnostics"){
cat("\n")
print(x[["diagnostics"]])
return(invisible())
}else if(attr(x, "type") == "individual"){
for(i in seq_along(x[["models"]])){
if(i > 1){
cat("\n")
}
print(x[["models"]][[i]][["summary"]])
cat("\n")
print(x[["models"]][[i]][["estimates"]])
}
return(invisible())
}
}
#' @title Reports whether x is a RoBMA object
#'
#' @param x an object to test
#'
#'
#' @return returns a boolean.
#'
#' @name is.RoBMA
#' @aliases is.RoBMA.reg is.NoBMA is.NoBMA.reg is.BiBMA
#' @export is.RoBMA
#' @export is.RoBMA.reg
#' @export is.NoBMA
#' @export is.NoBMA.reg
#' @export is.BiBMA
#' @rdname is.RoBMA
is.RoBMA <- function(x){
inherits(x, "RoBMA")
}
#' @rdname is.RoBMA
is.RoBMA.reg <- function(x){
inherits(x, "RoBMA.reg")
}
#' @rdname is.RoBMA
is.NoBMA <- function(x){
inherits(x, "NoBMA")
}
#' @rdname is.RoBMA
is.NoBMA.reg <- function(x){
inherits(x, "NoBMA.reg")
}
#' @rdname is.RoBMA
is.BiBMA <- function(x){
inherits(x, "BiBMA")
}
.object_title <- function(object){
if(is.NoBMA(object)){
return("Bayesian model-averaged meta-analysis (normal-normal model)")
}else if(is.NoBMA.reg(object)){
return("Bayesian model-averaged meta-regression (normal-normal model)")
}else if(is.BiBMA(object)){
return("Bayesian model-averaged meta-analysis (binomial-normal model)")
}else if(is.RoBMA.reg(object)){
return("Robust Bayesian meta-regression")
}else if(is.RoBMA(object)){
return("Robust Bayesian meta-analysis")
}
}
#' @title Interprets results of a RoBMA model.
#'
#' @description \code{interpret} creates a brief textual summary
#' of a fitted RoBMA object.
#'
#' @inheritParams summary.RoBMA
#'
#'
#' @return \code{interpret} returns a character.
#'
#' @export
interpret <- function(object, output_scale = NULL){
if(is.null(output_scale)){
output_scale <- object$add_info[["output_scale"]]
}else if(object$add_info[["output_scale"]] == "y" & .transformation_var(output_scale) != "y"){
stop("Models estimated using the generall effect size scale 'y' / 'none' cannot be transformed to a different effect size scale.")
}else{
output_scale <- .transformation_var(output_scale)
}
if(object$add_info[["output_scale"]] != output_scale){
object <- .transform_posterior(object, object$add_info[["output_scale"]], output_scale)
}
specification <- list()
if(any(names(object$RoBMA[["inference"]]) == "Effect")){
specification[["Effect"]] <- list(
inference = "Effect",
samples = "mu",
inference_name = "effect",
inference_BF_name = "BF_10",
samples_name = .transformation_names(object$add_info[["output_scale"]])
)
}
if(any(names(object$RoBMA[["inference"]]) == "Heterogeneity")){
specification[["Heterogeneity"]] <- list(
inference = "Heterogeneity",
samples = "tau",
inference_name = "heterogeneity",
inference_BF_name = "BF^rf",
samples_name = "tau"
)
}
if(any(names(object$RoBMA[["inference"]]) == "Bias")){
specification[["Bias"]] <- list(
inference = "Bias",
inference_name = "publication bias",
inference_BF_name = "BF_pb"
)
}
text <- BayesTools::interpret(
inference = object$RoBMA[["inference"]],
samples = object$RoBMA[["posteriors"]],
specification = specification,
method = .object_title(object)
)
return(text)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.