runehaubo / lmerTestR

Repository for the R-package lmerTest
48 stars 9 forks source link

update for lmertest models fails when inside a function #3

Closed singmann closed 6 years ago

singmann commented 6 years ago

Hi Rune,

While preparing the necessary fixes for the new lmerTest version for afex I noticed that one set of tests fails due to an apparent bug in lmertest. update does not work inside a function when passed an additional data argument.

library("lmerTest")

myupdate <- function(m, ...) {
  update(m, ...)
}

fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
update(fm1, data = sleepstudy) # works
myupdate(fm1, data = sleepstudy) # fails

fm2 <- lme4::lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
update(fm2, data = sleepstudy) # works
myupdate(fm2, data = sleepstudy) # works

As long as this is the case, my all_fit function fails and this crashes quite a few tests.

runehaubo commented 6 years ago

Argh, managing how things are evaluated in various environments is not exactly my favourite activity! But thanks for the bug report!

It should be fixed now in the eval-update-in-functions-branch. Could you perhaps check if that solves the problem with afex before I merge into master? I get:

> library(afex, quietly = TRUE)
************
Welcome to afex. For support visit: http://afex.singmann.science/
- Functions for ANOVAs: aov_car(), aov_ez(), and aov_4()
- Methods for calculating p-values with mixed(): 'KR', 'S', 'LRT', and 'PB'
- 'afex_aov' and 'mixed' objects can be passed to emmeans() for follow-up tests
- Get and set global package options with: afex_options()
- Set orthogonal sum-to-zero contrasts globally: set_sum_contrasts()
- For example analyses see: browseVignettes("afex")
************

Attaching package: ‘afex’

The following object is masked from ‘package:lme4’:

    lmer

> library(lmerTest)

Attaching package: ‘lmerTest’

The following object is masked from ‘package:lme4’:

    lmer

The following object is masked from ‘package:stats’:

    step

> require(optimx)
Loading required package: optimx
> data("sleepstudy", package="lme4")
> fm1 <- lmerTest::lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
> res <- all_fit(fm1)
bobyqa. : [OK]
Nelder_Mead. : [OK]
optimx.nlminb : [OK]
optimx.L-BFGS-B : [OK]
nloptwrap.NLOPT_LN_NELDERMEAD : [OK]
nloptwrap.NLOPT_LN_BOBYQA : [OK]
nmkbw. : [OK]

Cheers Rune

singmann commented 6 years ago

Thanks so much, that works for me as well.

And I completely agree. Evaluation in different environments is a never ending source of bugs and problems. I believe I once read a rant by Ben Bolker on R-devel where he suggested to source for all this is the way model.frame is implemented. I just try to avoid this as much as possible. However, I have one pending pull request on lme4 that also appears to suffer from this problem. This is really one of the few things in R that just suck (or I am too stupid).

One more thing. Do you have a concrete plan when to submit to CRAN? I am about to (after running the reverse dependency checks once more), but still get a note with the current CRAN version of lmerTest due to as_lmerModLmerTest not being exported. Maybe we can sync our submission and then just explain that this will solve the problem.

singmann commented 6 years ago

Okay, that was a little premature celebration. Now I get a bug in the lmerTest_anova function you wrote for afex. I am investigating this now to see if I can isolate it and then post again.

singmann commented 6 years ago

Got the bug:

library("lme4")

myfit <- function(formula, data) {
  lme4::lmer(formula = formula, data = data)
}

fm2 <- myfit(Reaction ~ Days + (Days|Subject), sleepstudy)
lmerTest::as_lmerModLmerTest(fm2)
# Error in if (!is.function(devfun) || names(formals(devfun)[1]) != "theta") stop("Unable to extract deviance function from model fit") : 
#   missing value where TRUE/FALSE needed
runehaubo commented 6 years ago

Thanks, I'll have a look.

On 31 March 2018 at 11:28, Henrik Singmann notifications@github.com wrote:

Got the bug:

library("lme4")

myfit <- function(formula, data) { lme4::lmer(formula = formula, data = data) }

fm2 <- myfit(Reaction ~ Days + (Days|Subject), sleepstudy) lmerTest::as_lmerModLmerTest(fm2)

Error in if (!is.function(devfun) || names(formals(devfun)[1]) != "theta") stop("Unable to extract deviance function from model fit") :

missing value where TRUE/FALSE needed

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/runehaubo/lmerTestR/issues/3#issuecomment-377679803, or mute the thread https://github.com/notifications/unsubscribe-auth/AE-crW8ZpTu-AQlSQyteDCXc-AJ3QlKOks5tj0w1gaJpZM4TBphN .

runehaubo commented 6 years ago

Hmm. Looking at this, I am beginning to doubt that this is ever going to work... I don't see how any function outside of myfit can know what data is and how to evaluate it:

> data("sleepstudy", package="lme4")
> myfit <- function(formula, data) {
+   fm2 <- lme4::lmer(formula = formula, data = data)
+   fm2
+ }
> fm2 <- myfit(Reaction ~ Days + (Days|Subject), sleepstudy)
> getCall(fm2)
lme4::lmer(formula = Reaction ~ Days + (Days | Subject), data = data)
> environment(myfit)
<environment: R_GlobalEnv>
> environment(formula(fm2))
<environment: R_GlobalEnv>
> environment(fm2)
NULL
> eval(getCall(fm2))
Error in list2env(data) : first argument must be a named list
> fm3 <- lme4::lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
> eval(getCall(fm3))
Linear mixed model fit by REML ['lmerMod']
Formula: Reaction ~ Days + (Days | Subject)
   Data: sleepstudy
REML criterion at convergence: 1743.628
Random effects:
 Groups   Name        Std.Dev. Corr
 Subject  (Intercept) 24.740       
          Days         5.922   0.07
 Residual             25.592       
Number of obs: 180, groups:  Subject, 18
Fixed Effects:
(Intercept)         Days  
     251.41        10.47  

But maybe I've gone blind and lost in the "environmental maze"...

If, instead, we move the call to as_lmerModLmerTest inside the function, it appears to work alright:

> myfun <- function(formula, data) {
+   fm <- lme4::lmer(formula = formula, data = data)
+   fm2 <- lmerTest::as_lmerModLmerTest(fm)
+   anova(fm2)
+ }
> myfun(Reaction ~ Days + (Days|Subject), sleepstudy)
Type III Analysis of Variance Table with Satterthwaite's method
     Sum Sq Mean Sq NumDF DenDF F value    Pr(>F)    
Days  30031   30031     1    17  45.853 3.264e-06 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# Now with lmerTest_anova:
> myfun2 <- function(formula, data) {
+   fm <- lme4::lmer(formula = formula, data = data)
+   lmerTest_anova(fm)
+ }
> myfun2(Reaction ~ Days + (Days|Subject), sleepstudy)
Type III Analysis of Variance Table with Satterthwaite's method
     Sum Sq Mean Sq NumDF DenDF F value    Pr(>F)    
Days  30031   30031     1    17  45.853 3.264e-06 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Anyway, long story short: if you could point me to when and where you see the problem in afex, I'll be happy to take a look at afex+lmerTest (I think at this point I would like to be able to reproduce the original problem myself).

Cheers Rune

runehaubo commented 6 years ago

Actually, the only problem I see when I check afex is

...
1. Error: mixed allows both lme4 and lmerText calls and exports lmerTest::lmer (@test-lmerTest-support.R#26) 
could not find function "skip_if"
1: .handleSimpleError(function (e) 
   {
       e$call <- sys.calls()[(frame + 11):(sys.nframe() - 2)]
       register_expectation(e, frame + 11, sys.nframe() - 2)
       signalCondition(e)
   }, "could not find function \"skip_if\"", quote(skip_if(packageVersion(pkg = "lmerTest") < 
       pkg_version))) at /Users/rhbc/GitHub/afex/package/tests/testthat/test-lmerTest-support.R:26
2: eval(code, test_env)

DONE ===================================================================================

but if you use skip_if_not instead of skip_if afex appears to check out with the new lmerTest :-)

Cheers Rune

singmann commented 6 years ago

Hi Rune,

Sorry for the delay (I was travelling today) and this indeed appears to be my error. When I test it now again, the tests run fine. So the only question open for me is now when to submit or not. Let me know if you have any concrete plans, but I will nevertheless close the issue now.

Cheers, Henrik

runehaubo commented 6 years ago

Hi Henrik,

Thanks for the feedback.

The plan is to submit the new lmerTest "3 weeks from the date of [the] email" which notified maintainers so about March 15th. I see an additional 3 packages no longer failing checks with the newest version of afex, so in a way it would be nice to have the new version of afex on CRAN for the other packages to lean on. If you would rather hold the CRAN submission until ~March 15th (and potentially Depend on lmerTest >= 3.0-0) I completely understand. Let me know what you prefer.

Cheers Rune