Nothing
if (.Platform$OS.type != "windows") {
library(lme4)
library(testthat)
source(system.file("testdata/lme-tst-funs.R", package="lme4", mustWork=TRUE))
##-> gSim(), a general simulation function ...
## hand-coded Pearson residuals {for sumFun() }
mypresid <- function(x) {
mu <- fitted(x)
(getME(x,"y") - mu) * sqrt(weights(x)) / sqrt(x@resp$family$variance(mu))
}
## should be equal (up to numerical error) to weights(.,type="working")
workingWeights <- function(mod) mod@resp$weights*(mod@resp$muEta()^2)/mod@resp$variance()
##' Sum of weighted residuals, 4 ways; the last three are identical
sumFun <- function(m) {
wrss1 <- m@devcomp$cmp["wrss"]
wrss2 <- sum(residuals(m,type="pearson")^2)
wrss3 <- sum(m@resp$wtres^2)
## compare to hand-fitted Pearson resids ...
wrss4 <- sum(mypresid(m)^2)
c(wrss1,wrss2,wrss3,wrss4)
}
## The relative "error"/differences of the weights w[] entries
rel.diff <- function(w) abs(1 - w[-1]/w[1])
set.seed(101)
## GAMMA
g0 <- glmer(y~x+(1|block),data=gSim(),family=Gamma)
expect_true(all(rel.diff(sumFun(g0)) < 1e-13))
expect_equal(weights(g0, type = "working"), workingWeights(g0),
tolerance = 1e-4) ## FIXME: why is such a high tolerance required?
## BERNOULLI
g1 <- glmer(y~x+(1|block),data=gSim(family=binomial(),nbinom=1),
family=binomial)
expect_true(all(rel.diff(sumFun(g1)) < 1e-13))
expect_equal(weights(g1, type = "working"), workingWeights(g1),
tolerance = 1e-5) ## FIXME: why is such a high tolerance required?
## POISSON
(n <- nrow(d.P <- gSim(family=poisson())))
g2 <- glmer(y ~ x + (1|block), data = d.P, family=poisson)
g2W <- glmer(y ~ x + (1|block), data = d.P, family=poisson, weights = rep(2,n))
expect_true(all(rel.diff(sumFun(g2 )) < 1e-13))
expect_true(all(rel.diff(sumFun(g2W)) < 1e-13))
## correct
expect_equal(weights(g2, type = "working"), workingWeights(g2),
tolerance = 1e-5) ## FIXME: why is such a high tolerance required?
expect_equal(weights(g2W, type = "working"), workingWeights(g2W),
tolerance = 1e-5) ## FIXME: why is such a high tolerance required?
## non-Bernoulli BINOMIAL
g3 <- glmer(y ~ x + (1|block), data= gSim(family=binomial(), nbinom=10),
family=binomial)
expect_true(all(rel.diff(sumFun(g3)) < 1e-13))
expect_equal(weights(g3, type = "working"), workingWeights(g3),
tolerance = 1e-4) ## FIXME: why is such a high tolerance required?
d.b.2 <- gSim(nperblk = 2, family=binomial())
g.b.2 <- glmer(y ~ x + (1|block), data=d.b.2, family=binomial)
expect_true(all(rel.diff(sumFun(g.b.2 )) < 1e-13))
## Many blocks of only 2 observations each - (but nicely balanced)
## Want this "as" https://github.com/lme4/lme4/issues/47
## (but it "FAILS" survival already):
##
## n2 = n/2 :
n2 <- 2048
if(FALSE)
n2 <- 100 # for building/testing
set.seed(47)
dB2 <- gSim(n2, nperblk = 2, x= rep(0:1, each= n2), family=binomial())
## -- -- --- --------
gB2 <- glmer(y ~ x + (1|block), data=dB2, family=binomial)
expect_true(all(rel.diff(sumFun(gB2)) < 1e-13))
## NB: Finite sample bias of \hat\sigma_1 and \hat\beta_1 ("Intercept")
## tend to zero only slowly for n2 -> Inf, e.g., for
## n2 = 2048, b1 ~= 4.3 (instead of 4); s1 ~= 1.3 (instead of 1)
## FAILS -----
## library(survival)
## (gSurv.B2 <- clogit(y ~ x + strata(block), data=dB2))
## ## --> Error in Surv(rep(1, 200L), y) : Time and status are different lengths
## summary(gSurv.B2)
## (SE.surf <- sqrt(diag(vcov(gSurv.B2))))
g3 <- glmer(y ~ x + (1|block),data=gSim(family=binomial(),nbinom=10),
family=binomial)
expect_equal(var(sumFun(g3)),0)
## check dispersion parameter
## (lowered tolerance to pass checks on my machine -- SCW)
expect_equal(sigma(g0)^2, 0.4888248, tolerance=1e-4)
} ## skip on windows (for speed)
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.