pitakakariki / simr

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

fcompare doesn't warn when comparing identical models #252

Open KJames5 opened 1 year ago

KJames5 commented 1 year ago

Hi! I'm having issues running a simple GLM in simr(). I want to test the effect of 1 factor (3 groups), with site included as a random factor (3 groups). I don't have any issues when I include 2 independent factors in the model (as a interaction or just including both in the model) so I'm not sure why it is having errors. When I run a powerSim(), the output has a power of 0%, with a warning of "Test did not return a p-value" for each simulation. Here is my code:

tilej <- factor(1:5)
sites_id <- factor(1:3)
tile <- c("control", "A", "B")
##3 types of material
material <- c("standard", "lc concrete", "lc concrete 2")

###Create enough tile values for each tile*material*site(18)
site_full <- rep(sites_id, each=30)
length(site_full)

tile_full <- rep(rep(tile, each=5),6)
length(tile_full)

material_full <- rep(rep(material, each=15),2)
length(material_full)

fact <- data.frame(site=factor(site_full), tile=factor(tile_full), 
                   material=factor(material_full))

##Testing y ~ tile type + (1|site)
fact%>%
  group_by(site,tile)%>%
  summarise(length=length(tile))

fixed <- c(10, -1, 3)
rand <- list(2)
res <- 2

model <- makeLmer(y ~ tile + (1|site), 
                  fixef=fixed, VarCorr=rand, sigma=res, data=fact)
model

sim_treat <- powerSim(model, nsim=1000, test = fcompare(y~tile))
sim_treat

lastResult()$errors

Do you have any idea why this might be? I have no issues when using the same data frame when I include material in the model. Thanks!

pitakakariki commented 1 year ago

Problem is with the test. fcompare keeps random effects the same and modifies the fixed part of the model, which is ~ tile in both cases.

If you were trying to test the fixed effect you probably want test = fcompare(~1).

If you were trying to test the random effect probably test = random().

KJames5 commented 1 year ago

Hi,

Thanks that makes more sense! Thanks for making that clear!