tests/testthat_backup/test-exact-Kakadu.R

In certifiedwaif/blma: Bayesian Linear Model Averaging

context("exact Kakadu")

library(parallel)
library(tictoc)

cores <- detectCores()

# modelprior = uniform
test_that("Kakadu produces correct results BIC", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results BIC")
	result <- blma(vy, mX, prior="BIC", modelprior="uniform", cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.11956375866122319,0.43598166701152036,0.03004602756871903,
		0.37141951236536819,0.81868011537400154,0.16833752435869168,
		0.03221165979213939,0.04298271555297178,0.02617847083497719,
		0.52531028106054500,0.92513375271585307,0.99815130365593652,
		0.02454673388614505,0.08100948178569978,0.08169622092074036,
		0.62985862997075337,0.03271752259200703,0.54749994288839865,
		0.99999999999905875,0.26634597085271527,0.99999886673686156,
		0.04953096133476929
	), tolerance = 1e-8)
})

test_that("Kakadu produces correct results ZE", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results ZE")
	result <- blma(vy, mX, prior="ZE", modelprior="uniform", cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.20355684150991538,0.47244398575065949,0.07485696352142855,
		0.42017464791967363,0.86114935243204360,0.26671059509879097,
		0.08889221213490464,0.11092303936894550,0.07192086903966989,
		0.77777906132245778,0.93733237674906777,0.99937746373248582,
		0.06601135712668889,0.19907890141491377,0.18511189624094807,
		0.75297589180094671,0.08527867659504840,0.74927591041652553,
		0.99999999999929812,0.44106509651296133,0.99999873524470562,
		0.13221495773128050
	), tolerance = 1e-8)
})

test_that("Kakadu produces correct results liang_g1", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results liang_g1")
	result <- blma(vy, mX, prior="liang_g1", modelprior="uniform", cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.3469442223401997,0.5033763713275147,0.1710391998737254,0.4686716902393477,
		0.9041282693203453,0.4182787535192226,0.2152694787243199,0.2365794398160191,
		0.1709195101828257,0.9081285313534126,0.9457607536048467,0.9996987392264256,
		0.1584196650685204,0.3866126303418119,0.3524291472817372,0.8329753936885851,
		0.1953702526810372,0.8655055789962616,0.9999999999992009,0.6255640528176040,
		0.9999982103316476,0.2911572010612115
	), tolerance = 1e-8)
})


test_that("Kakadu produces correct results liang_g2", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results liang_g2")
	result <- blma(vy, mX, prior="liang_g2", modelprior="uniform", cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.3469442223402000,0.5033763713275141,0.1710391998737261,0.4686716902393479,
		0.9041282693203457,0.4182787535192236,0.2152694787243207,0.2365794398160196,
		0.1709195101828260,0.9081285313534141,0.9457607536048475,0.9996987392264269,
		0.1584196650685208,0.3866126303418126,0.3524291472817367,0.8329753936885853,
		0.1953702526810377,0.8655055789962616,0.9999999999992023,0.6255640528176055,
		0.9999982103316490,0.2911572010612116
	), tolerance = 1e-8)
})

# We get NANs. This simply doesn't work. Move on for now.
test_that("Kakadu produces correct results liang_g_n_appell", {
	skip("We were unable to get this case to work. It's too numerically difficult.")
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results liang_g_n_appell")
	result <- blma(vy, mX, prior="liang_g_n_appell", modelprior="uniform", cores=1)
	toc()
	expect_equal(result$vinclusion_prob, c(
	), tolerance = 1e-8)
})

test_that("Kakadu produces correct results liang_g_n_approx", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results liang_g_n_approx")
	result <- blma(vy, mX, prior="liang_g_n_approx", modelprior="uniform", cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.3296136691294210,0.4997609497180022,0.1579846499019243,0.4631387408540441,
		0.9006995511063154,0.4005475083511888,0.1982925439416997,0.2209076900143787,
		0.1575545973403185,0.8998253808608440,0.9453194636814199,0.9996900361221982,
		0.1456947072049949,0.3654752994812290,0.3329298710932005,0.8268178993825527,
		0.1806657458627090,0.8573731512895385,0.9999999999992192,0.6083132163536824,
		0.9999983508359184,0.2714475290242928
	), tolerance = 1e-8)
})

test_that("Kakadu produces correct results liang_g_n_quad", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results liang_g_n_quad")
	result <- blma(vy, mX, prior="liang_g_n_quad", modelprior="uniform", cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.3203591831238068,0.4977774306550474,0.1512859110442455,0.4601315123732897,
		0.8984723576023463,0.3910244742565688,0.1894701351919556,0.2126145072480971,
		0.1506655077406368,0.8944304522063428,0.9448935925720320,0.9996801335893382,
		0.1391831130834718,0.3538609740786210,0.3223705659773634,0.8229190449064608,
		0.1730543020195790,0.8521810234828430,0.9999999999993161,0.5982779020650832,
		0.9999983946351332,0.2609239705895227
	), tolerance = 1e-8)
})

test_that("Kakadu produces correct results robust_bayarri1", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results robust_bayarri1")
	result <- blma(vy, mX, prior="robust_bayarri1", modelprior="uniform", cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.2645609801935807,0.4872708351002377,0.1119549057643435,0.4428310083494666,
		0.8859403686981394,0.3327165328332443,0.1381842637777420,0.1633673877360185,
		0.1104117927906037,0.8591578715890693,0.9435211879410100,0.9996313610321645,
		0.1012638529130835,0.2837032904837962,0.2586857176889374,0.7998311650018346,
		0.1284988708797292,0.8194918851464008,1.0000000000003251,0.5357989689853954,
		0.9999987783206028,0.1987407247561971
	), tolerance = 1e-8)
})

test_that("kakadu produces correct results robust_bayarri2", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("kakadu produces correct results robust_bayarri2")
	toc()
	result <- blma(vy, mX, prior="robust_bayarri2", modelprior="uniform", cores=cores)
	expect_equal(result$vinclusion_prob, c(
		0.2645609801933496,0.4872708350999790,0.1119549057642774,0.4428310083492054,
		0.8859403686976074,0.3327165328329659,0.1381842637776008,0.1633673877358905,
		0.1104117927905364,0.8591578715885575,0.9435211879405069,0.9996313610316295,
		0.1012638529130122,0.2837032904835066,0.2586857176888843,0.7998311650014662,
		0.1284988708795221,0.8194918851458997,0.9999999999992317,0.5357989689844272,
		0.9999987783204403,0.1987407247563582
	), tolerance = 1e-5)
})

test_that("Kakadu produces correct results zellner_siow_gauss_laguerre", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	tic("Kakadu produces correct results zellner_siow_gauss_laguerre")
	result <- blma(vy, mX, prior="zellner_siow_gauss_laguerre", modelprior="uniform", cores=1)
	toc()
	expect_equal(result$vinclusion_prob, 
c(
0.20578507196599863,0.47307551042974938,0.07612923670551421,
0.42112901060440383,0.86215067202350637,0.26921889275215943,
0.09053504955485733,0.11279139046662487,0.07323093390440316,
0.78214745397319363,0.93756983435089636,0.99939292585211370,
0.06720265691059692,0.20221362639019935,0.18783534784339340,
0.75520354356412289,0.08676136168658173,0.75285716863824470,
0.99999999999929501,0.44505775064088876,0.99999873218671564,
0.13452936396921064
)
	, tolerance = 1e-8)
})

# prior = BIC, modelprior = beta-binomial
test_that("Kakadu produces correct results modelprior beta-binomial", {
	Kakadu <- get_Kakadu()
	vy <- Kakadu$vy
	mX <- Kakadu$mX
	p <- ncol(mX)
	tic("Kakadu produces correct results modelprior beta-binomial")
	result <- blma(vy, mX, prior="BIC", modelprior="beta-binomial", modelpriorvec = c(1, p), cores=cores)
	toc()
	expect_equal(result$vinclusion_prob, c(
		0.056281843799996970,0.306314199561247891,0.007552489928519051,
		0.247616597500152053,0.742709650039889202,0.095631352679238363,
		0.006351015631741534,0.009101941725939111,0.005241750496263586,
		0.153817084200635196,0.907193762314498886,0.984666525605746967,
		0.005113387083452965,0.017754736585571605,0.021245104548298487,
		0.371729547784183711,0.006733387430678084,0.216129983894335492,
		0.999999999998894218,0.110633351931119187,0.999999030703775826,
		0.009702560359548688
	), tolerance = 1e-8)
})