Nothing
library(testthat)
library(metaBMA)
test_that("Scheibehenne (2017): model averaging for towels dataset works correctly", {
data(towels, package = "metaBMA")
tmp <- capture_output(eval(towels))
############ TESTING: scale = 1/4
# d.rand ~ dnorm(0, 1/0.3^2)T(0,)
# prec <- 1/(tau*tau) # btw-study precision
# tau ~ dt(0, prior.scaleTau, 1) T(0,)
# sdVec <- c(1/64, 1/32, 1/16, 1/8, 1/4, 1/2, 1, 10, 100) # TEST
# scaleTauVec <- 1/sdVec^2
expect_silent(prior_d <- prior(
family = "norm", param = c(mean = 0, sd = 0.3),
lower = 0, label = "d"
))
expect_silent(prior_d_unrestricted <- prior(family = "norm", param = c(mean = 0, sd = 0.3)))
expect_silent(prior_tau <- prior(
family = "t", param = c(location = 0, scale = 1 / 4, nu = 1),
lower = 0, label = "tau"
))
set.seed(123)
suppressWarnings( # divergent transitions
m_testing <- meta_bma(
y = logOR, SE = SE, data = towels,
d = prior_d, tau = prior_tau,
iter = 1500, logml_iter = 2000, rel.tol = .1, summ = "stan"
)
)
BF_fixed_random <- 2.450923
BFplus0_random <- 4.055704
BFplus0_fixed <- 23.88449
BFincl <- 9.882879
expect_equal(m_testing$BF["fixed_H1", "random_H1"], BF_fixed_random, tolerance = .03)
expect_equal(m_testing$BF["random_H1", "random_H0"], BFplus0_random, tolerance = .01)
expect_equal(m_testing$BF["fixed_H1", "fixed_H0"], BFplus0_fixed, tolerance = .01)
expect_equal(m_testing$inclusion$incl.BF, BFincl, tolerance = .01)
skip_on_cran()
############ ESTIMATION
m.fixed <- 0.2116763
m.rand <- 0.1848269
m.mixed <- 0.2033756
hpd95.fixed <- c(0.063656, 0.3623073)
hpd95.rand <- c(-0.02162486, 0.38411546)
hpd95.mixed <- c(0.03867814, 0.37240522)
suppressWarnings(
meta_est <- meta_bma(
y = logOR, SE = SE, data = towels,
d = prior_d_unrestricted, tau = prior_tau,
iter = 1500, logml_iter = 1000, rel.tol = .1, summ = "stan"
)
)
expect_equal(meta_est$estimates["fixed", "mean"], m.fixed, tolerance = .01)
expect_equal(unname(meta_est$estimates["fixed", c("hpd95_lower", "hpd95_upper")]),
hpd95.fixed,
tolerance = .03
)
expect_equal(meta_est$estimates["random", "mean"], m.rand, tolerance = .01)
expect_equal(unname(meta_est$estimates["random", c("hpd95_lower", "hpd95_upper")]),
hpd95.rand,
tolerance = .03
)
expect_equal(meta_est$estimates["averaged", "mean"], m.mixed, tolerance = .01)
expect_equal(unname(meta_est$estimates["averaged", c("hpd95_lower", "hpd95_upper")]),
hpd95.mixed,
tolerance = .03
)
})
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.