0

I am examining a certain fixed effect (effect3 here) by fitting two negative binomial mixed models using the glmmTMB function: one with and one without the fixed effect. Next, I perform a likelihood ratio test to test the significance of effect3 (example data):

eff1 <- c(0.026, 0.003, -0.008, -0.057, -0.022)

eff2 <- c(-0.002, -0.013, 0.036, 0.005, 0.074)

eff3 <- c(0.027, 0.021, -0.015, 0.008, 0)

counts <- c(18317, 19899, 11048, 23920, 20656)

data.eff <- data.frame(effect1 = eff1, effect2 = eff2, effect3 = eff3, Val = counts)

modeling <- function(m.data, vars){
  # This is how the model should look like
  f <- reformulate(termlabels = vars, response = outcome)
  model <- glmmTMB(formula = f, data = m.data, family = nbinom2())
}

outcome <- "Val"
variables <- colnames(data.eff[,-which(colnames(data.eff) == "Val")])
variableswo <- colnames(data.eff[,-which(colnames(data.eff) == "Val")])[-3]

model1 <- modeling(m.data = data.eff, vars = variables)
model2 <- modeling(m.data = data.eff, vars = variableswo)

res <- lrtest(model1, model2)

Sadly, I cannot find a way to evaluate the effect size; or rather none of the numerous potential approaches seem to work. One of the ideas that failed was calculating Cohen's f-squared:

r22 <- r2(model1)
r21 <- r2(model2)
f2 <- (r22-r21)/(1-r22)

However, r2(model) always results in NA values when using the above code. I think the bug is in the function, because the commands within the function work when run one-by-one but calling the arguments seems to crash it. Maybe it has got to do with the environment?

Therefore, my questions are:

  1. Can I use glmmTMB() for a fixed effects model so a model without random effects?
  2. How do I calculate such an effect size? Is Cohen's f2 the right choice?
  3. How do I implement that in R given the code above? How can I fix that function?

Thank you!

Friedebert
  • 11
  • 2
  • can you be more specific what you mean by "effect size"? There are a **lot** of definitions ... – Ben Bolker Oct 05 '20 at 23:58
  • Hmmm, I was thinking about something like Cohen's f^2 (like shown above). Maybe that approach didn't work because of wrong R-squared calculations but I didn't find any alternative for calculating R^2 with such a glmmTMB model so far... – Friedebert Oct 06 '20 at 15:33
  • can you give us a [mcve]? the couple of examples I tried with `r2_efron()` and `r2()$R2_marginal()` didn't have the pathology you described ... – Ben Bolker Oct 06 '20 at 17:59
  • The examples I can generate here all work with `r2()`. My data, however, produce NAs as a result from `r2()`. `r2_efron()` works in all cases, but the effect sizes I get using the code above are as high as 14 which cannot be true. I have got 18 fixed effects with numeric values between -1 and 1. I get a result for `r2()` when using the same "Counts" but only a fraction of the fixed effects. I also get a result when using less sparse and generally higher "Counts" data and all fixed effects. Is this an issue with the variance? Using glm instead of glmmTMB also produces `r2()` results. – Friedebert Oct 07 '20 at 11:12
  • Sorry, to clarify the last sentence: this `model1 <- glm(formula = Counts ~ effect1 + effect2 + effect3, data = data, family = nbinom2())` gives a result for `r2()` while this `model1 <- glmmTMB(formula = Counts ~ effect1 + effect2 + effect3, data = data, family = nbinom2())` does not. Both models look very different, though. Which one is the better option, though? – Friedebert Oct 07 '20 at 11:26
  • That 'r2()' cannot be calculated from some of the models must be an issue with the data frames in conjunction with the function I am using to compute glmmTMB (see question above). Two identical data frames, contianing the exact same numbers in the same format and order etc., perform differently. One is created in a long process involving `merge()`, one is directly made from vectors and `data.frame()`. Both result in the exact same model, but for the newly created, 'r2()' gives a result while for the original one, `r2()`gives NA. – Friedebert Oct 09 '20 at 16:37
  • again, can you post a [mcve] (see also [here](https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example))? Without the actual data it's almost impossible to answer this. – Ben Bolker Oct 09 '20 at 18:31
  • Sorry for the mess... I changed the question above incl. an example, this should sum up the issue? Thank you! – Friedebert Oct 14 '20 at 09:40
  • Update: The problem seems to be the function: the arguments weren't read appropriately or rather lost in the environment somewhere. When assigning the data do the global environment, it suddenly worked. This is a very dumb workaround, but I couldn't get any further than this. And yes, cohen's f2 seems to be the correct option to calculate effect size here. – Friedebert Oct 29 '20 at 15:07

0 Answers0