Nothing
context("data sets [flatlizards]")
tol <- 1e-6
## standard BT model, using the bias-reduced maximum likelihood method:
result <- rep(1, nrow(flatlizards$contests))
BTmodel <- BTm(result, winner, loser, br = TRUE,
data = flatlizards$contests)
## "structured" B-T model: abilities are determined by a linear predictor.
Whiting.model1 <- BTm(result, winner, loser,
~ throat.PC1[..] + throat.PC3[..] +
head.length[..] + SVL[..],
family = binomial, data = flatlizards)
## Equivalently, fit the same model using glmmPQL:
Whiting.model1b <- BTm(result, winner, loser,
~ throat.PC1[..] + throat.PC3[..] +
head.length[..] + SVL[..] + (1|..),
sigma = 0, sigma.fixed = TRUE, data = flatlizards)
## Same predictor but with a normally distributed error
Whiting.model2 <- BTm(result, winner, loser,
~ throat.PC1[..] + throat.PC3[..] +
head.length[..] + SVL[..] + (1|..),
data = flatlizards)
## Now use probit rather than logit as the link function:
Whiting.model3 <- BTm(result, winner, loser,
~ throat.PC1[..] + throat.PC3[..] +
head.length[..] + SVL[..] + (1|..),
family = binomial(link = "probit"), data = flatlizards)
test_that("standard model as expected on flatlizards", {
# check standard model
# ignore family: mode of initialize changes between R versions
res <- summary(BTmodel)
res$family <- NULL
expect_known_value(res,
file = test_path("outputs/flatlizards-BTmodel.rds"),
tol = tol)
# check structured model against Table 1 of Whiting et al. (2006)
# (for coefficients of covariates only, not separate lizard effects)
cf <- coef(summary(Whiting.model1))[-(1:2),]
expect_equal(unname(round(cf[, "Estimate"], 2)),
c(-0.09, 0.34, -1.13, 0.19))
expect_equal(unname(round(cf[, "Std. Error"], 2)),
c(0.03, 0.11, 0.49, 0.1))
# reported Z stat appear to be Chi-squared stat
expect_equal(unname(round(cf[, "z value"]^2, 1)),
c(10.3, 9.5, 5.2, 3.6), tol = 1e-1)
expect_equal(unname(signif(cf[, "Pr(>|z|)"], 1)),
c(0.001, 0.002, 0.02, 0.06))
# check equiv glmmPQL against Table 1 of Whiting et al. (2006)
# (for coefficients of covariates only, not separate lizard effects)
cf <- coef(summary(Whiting.model1b))[-(1:2),]
expect_equal(unname(round(cf[, "Estimate"], 2)),
c(-0.09, 0.34, -1.13, 0.19))
expect_equal(unname(round(cf[, "Std. Error"], 2)),
c(0.03, 0.11, 0.49, 0.1))
# reported Z stat appear to be Chi-squared stat
expect_equal(unname(round(cf[, "z value"]^2, 1)),
c(10.3, 9.5, 5.2, 3.6), tol = 1e-1)
expect_equal(unname(signif(cf[, "Pr(>|z|)"], 1)),
c(0.001, 0.002, 0.02, 0.06))
})
test_that("GLMM models as expected on flatlizards", {
## The estimated coefficients (of throat.PC1, throat.PC3,
## head.length and SVL are not changed substantially by
## the recognition of an error term in the model
cf <- coef(summary(Whiting.model1b))[-(1:2),]
cf2 <- summary(Whiting.model2)$fixef[-(1:2),]
expect_equal(cf[, "Estimate"],
cf2[, "Estimate"], tol = 0.5)
## but the estimated
## standard errors are larger, as expected. The main conclusions from
## Whiting et al. (2006) are unaffected.
expect_true(all(cf2[, "Std. Error"] > cf[, "Std. Error"]))
## Modulo the usual scale change between logit and probit, the results
## are (as expected) very similar to Whiting.model2.
cf3 <- summary(Whiting.model3)$fixef[-(1:2),]
expect_equal(unname(cf2[, "Estimate"]/cf3[, "Estimate"]),
rep(1.6, 4), tol = 0.1)
## drop lizard 996as coef not estimable !! should be 96
abilities <- BTabilities(Whiting.model3)[-55,]
expect_known_value(abilities,
file = test_path("outputs/flatlizards-abilities.rds"),
tol = tol)
resids <- residuals(Whiting.model3, "grouped")
expect_known_value(resids,
file = test_path("outputs/flatlizards-residuals.rds"),
tol = tol)
})
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.