pitakakariki / simr

Power Analysis of Generalised Linear Mixed Models by Simulation
70 stars 19 forks source link

Inconsistent results when using extend() #92

Open smonsays opened 6 years ago

smonsays commented 6 years ago

I am not really shure whether I do / understand something wrong but I get very inconsistent results from powerSim() when using extend(). I do the following:

# Import data set of 5 participants with each 1600 observations
dat <- read.csv("data_N5.csv")

# Calculate mixed effect models 
mixedModel.null   <- lmer(log(rt) ~ group_num + err + freq +
                            (1|vp) + (1|string), data=dat, REML = FALSE)

mixedModel.h1     <- lmer(log(rt) ~ scale(ope) * group_num + err + freq +
                            (1|vp) + (1|string), data=dat, REML = FALSE)

# Specify expected fixed effect size
fixef(mixedModel.h1)["scale(ope)"] <- -0.05

# Perform power calculation with original number of participants
powerTest1 <- powerSim(mixedModel.h1, test = compare(mixedModel.null), nsim = 10)

# Extend number of participants to N = 40 and repeat power calculation
mixedModel.h1Ext <- extend(mixedModel.h1, along="vp", n= "40")
powerTest2 <- powerSim(mixedModel.h1Ext, test = compare(mixedModel.null), nsim = 10)

While powerTest1 returns already a value of around β=0.8, powerTest2 always returns β=0.00. Is there an obvious error in my calculations or should this theoretically work?

Thanks for your help!

pitakakariki commented 6 years ago

Best to check the obvious things first.

1) Compare nrow(getData(mixedModel.h1)) to nrow(getData(mixedModel.h1Ext)) to see if extend really increased your sample size. If I've understood correctly, this should have increased from 8,000 to 64,000?

2) Check if there were any errors during the second set of simulations, e.g. head(powerTest2$errors).

smonsays commented 6 years ago

Thanks for your quick reply. I rechecked both points: (1) The sample size increases correctly to 64,000 and (2) there are no errors nor warnings after both simulations.

pitakakariki commented 6 years ago

I can't see any obvious problems with your code. If you're willing to email me your data (email at this link https://cran.r-project.org/web/packages/simr/) I'd be able to have a closer look.