Closed JCruk closed 4 years ago
That looks like maybe simr's heuristics are choosing the wrong test.
Does it work if you manually specify the test? i.e. fixed("A", "???")
.
Is A
numeric or categorical here?
A is categorical here. I tried powerSim(Mod, test = fixed("A", "z"), nsim = 25)
with the same no applicable method... result
How about if you use simr::fixed
in place of fixed
?
That's what I was thinking too - and it was a function masking issue. I usually load simr last, but didn't here. Sorry for the bother!
I'm now getting "Test did not return a p-value" with all interaction terms. I suspect that's due to model specification rather than anything in simr. Is that a correct assumption?
That's a catch-all in the test wrapper, so could be a number of things, but most likely an NA
p-value.
You should be able to bypass the wrapper with simr::lrtest(Mod, "A")
to see what's happening.
Thanks!
I get an error trying to bypass the wrapper:
Error: 'lrtest' is not an exported object from 'namespace:simr
Sorry, triple-double dots: simr:::lrtest(Mod, "A")
.
So there is an NA
simr:::lrtest(Mod "A")
[1] NA
Is there anything I can do about that?
Try drop1(Mod, ~A, test="Chisq")
, which is the test simr's calling from lme4. If that doesn't work, it's probably a problem with your model, but you might get a more helpful error messsge.
That does work. Is there a way to combine the single term deletion with simr so I can plot out a powerCurve?
fixed("A", "lr")
or even just fixed("A")
should be using that drop1
result. If drop1
is working though, maybe the wrapper's having trouble extracting the p-value.
What does the drop1
output look like?
In the model with my actual data (where I get the no p-value errors):
powerSim(Valence_glmer, test = fixed("Category", "lr"), nsim = 10)
,
drop1(Valence_glmer, ~Category, test="Chisq")
returns:
Model:
valence ~ condition + filter_type + Category + condition:filter_type +
condition:Category + Category:filter_type + (1 | initials)
Df AIC LRT Pr(Chi)
<none> 27427
Category 0 27427 1.4767e-08
Aha, drop1
isn't returning a p-value. I suspect that the problem is with the zero degrees of freedom, although I have no idea what could cause that.
It's likely due to something in specification. When I run powerSim
on a "typical" glmer specification I get cannot find valid starting values: please specify some
for every simulation - even though the model runs just fine. So I added in mustart = pmax(data$y)
to the model call, but then I get the no p-value errors.
Have you encountered the cannot find valid starting values: please specify some
error previously?
I've tried a few more things and while my model runs fine in glmer with the specified values for mustart
, it seems something is happening when I pass the model over to powerCurve
.
When I try plot(powerCurve(Mod, test = fixed("A), along = "G", nsim = 10))
,
I get a flat plot at zero and lastResult()$errors shows:
variable lengths differ (found for '(mustart)')
for each simulation
In case it's helpful to be able to reproduce my issue:
## List of required packages
Pkgs <- c("tidyverse","magrittr", "lme4", "lmerTest", "MCMCglmm", "simr")
# Load packages
lapply(Pkgs, require, c = T)
## Build a simulated data set
Sub_Group <- data.frame(Subject = LETTERS[1:20],
Group = c(rep("Y", 10), rep("Z", 10)))
Cat_Item <- data.frame(expand.grid(Category = c(rep("U", 25), rep("V", 25)),
Item = factor(rep(c(1:25),2))))
Cond <- data.frame(Condition = letters[1:5])
Cond_Cat_Item <- merge(Cond, Cat_Item, by = NULL)
Data <- merge(Sub_Group, Cond_Cat_Item, by = NULL) %>%
mutate(Y = round(rtnorm(nrow(.), 4.5, 2, 1, 9), digits = 0))
## Run initial glmer without any specified starting values
summary(Mod1 <- glmer(Y ~ Category + Condition + Group +
(1 + Group|Subject),
data = Data,
family = gaussian(link = "log")))
plot(powerCurve(extend(Mod1, along = "Subject", n = 50), along = "Subject", nsim = 10))
# Errors suggesting to provide starting values
lastResult()$errors
# "cannot find valid starting values: please specify some"
## Second glmer with starting values
summary(Mod2 <- glmer(Y ~ Category + Condition + Group +
(1 + Group|Subject),
data = Data,
mustart = Data$Y,
family = gaussian(link = "log")))
plot(powerCurve(extend(Mod2, along = "Subject", n = 50), along = "Subject", nsim = 10))
## Errors regarding variable lengths
lastResult()$errors
# "variable lengths differ (found for '(mustart)')"
I am encountering the same problem. I am fairly new to simr
, so pardon me.
I have the following model:
m <- glmer('error ~ type+ (1 | p_id:group) + (1|Q)', data = df, family=gaussian(link="log"))
where type is a categorical variable with 4 levels, _pid is the participant id, group is the group id.
When I run
powerSim(m, test = simr::fixed("type"), nsim=20 )
I get the error "cannot find valid starting values: please specify some" in all 20 simulations, even if the glmer function runs smoothly. Running drop1(m, ~type, test="Chisq")
also does not give any problem and returns a p-value.
Any hint to how to fix it?
It looks like the problem is happening at the simulation stage (it works for simulate.merMod
but simr uses simulate.formula
).
Workaround for now would be a custom simulation function, e.g.:
p <- list(beta=fixef(m), theta=getME(m, "theta"), sigma=sigma(m))
s <- function() simulate(m, newparams=p)[[1]]
powerSim(m, sim=s, etc)
In fitting some gaussian log models, I am getting an error for each simulation in powerSim:
no applicable method for 'test' applied to an object of class "c('glmerMod', 'merMod')"
The model is relatively straightforward:
I then call powerSIm with
powerSim(Mod, test = fixed("A"), nsim = 100)
I get an estimated power of 0.00 with 100 errors (as above)
Has there been a change in how to call powerSim for glmer models?