R Under development (unstable) (2024-06-30 r86854 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> ## require(agridat)
> ## dat <- gotway.hessianfly
> 
> if (.Platform$OS.type != "windows") {
+ ## don't actually use gotway_hessianfly_fit or gotway_hessianfly_prof,
+ ## so we should be OK even with R< 3.0.1
+ load(system.file("testdata","gotway_hessianfly.rda",package="lme4"))
+ # Block random.  See Glimmix manual, output 1.18.
+ # Note: (Different parameterization)
+ 
+ ## require("lme4.0")
+ ## fit2 <- glmer(cbind(y, n-y) ~ gen + (1|block), data=dat, family=binomial)
+ ## params <- list(fixef=fixef(fit2),theta=getME(fit2,"theta"))
+ ## detach("package:lme4.0")
+ lme4.0fit <- structure(list(fixef = structure(c(1.50345713031203, -0.193853259383803, 
+ -0.540808391060274, -1.43419379979154, -0.203701042949808, -0.978322555343941, 
+ -0.604078624475678, -1.67742449813309, -1.39842466673692, -0.681709344788684, 
+ -1.46295367186169, -1.45908310198959, -3.55285756517073, -2.50731975980307, 
+ -2.08716296677356, -2.96974270029992), .Names = c("(Intercept)", 
+ "genG02", "genG03", "genG04", "genG05", "genG06", "genG07", "genG08", 
+ "genG09", "genG10", "genG11", "genG12", "genG13", "genG14", "genG15", 
+ "genG16")), theta = structure(0.0319087494293615, .Names = "block.(Intercept)")), .Names = c("fixef", 
+ "theta"))
+ 
+ ## start doesn't work because we don't get there
+ library(lme4)
+ m1 <- glmer(cbind(y, n-y) ~ gen + (1|block), data=gotway.hessianfly,
+             family=binomial)
+ m1B <- update(m1,control=glmerControl(optimizer="bobyqa"))
+ max(abs(m1@optinfo$derivs$gradient))  ## 0.0012
+ max(abs(m1B@optinfo$derivs$gradient)) ## 2.03e-5
+ abs(m1@optinfo$derivs$gradient)/abs(m1B@optinfo$derivs$gradient)
+ ## bobyqa gets gradients *at least* 1.64* lower
+ 
+ lme4fit <- list(fixef=fixef(m1),theta=getME(m1,"theta"))
+ 
+ ## hack around slight naming differences
+ lme4fit$theta <- unname(lme4fit$theta)
+ lme4.0fit$theta <- unname(lme4.0fit$theta)
+ ## difference in theta on x86_64-w64-mingw32 (64-bit) with r-devel is  0.000469576
+ stopifnot(all.equal(lme4fit, lme4.0fit, tolerance = 5e-4))
+ 
+ ## Fun stuff: visualize and alternative model
+ 
+ ## library(ggplot2)
+ ## dat$prop <- dat$y/dat$n
+ ## theme_set(theme_bw())
+ ## ggplot(dat,aes(x=gen,y=prop,colour=block))+geom_point(aes(size=n))+
+ ##     geom_line(aes(group=block,colour=block))+
+ ##     geom_smooth(family=binomial,aes(weight=n,colour=block,group=block),method="glm",
+ ##                 alpha=0.1)
+ 
+ ## dat$obs <- factor(seq(nrow(dat)))
+ ## m2 <- glmer(cbind(y, n-y) ~ block+ (1|gen) + (1|obs), data=dat, family=binomial)
+ 
+ }  ## not on windows/CRAN
> 
> proc.time()
   user  system elapsed 
   0.10    0.04    0.14