Test for differences in Mesor between levels of grouping variables???

[MilaSMayor]

Dear GLMMcosinor developer,

I would like to know if there is some function where I can find the Test for differences in Mesor between levels of grouping variable. In the function test_cosinor_levels the cosinor parameter to compare are acrophase and amplitude, not for Mesor.

Best,

Mila

[RWParsons]

Hi Mila,

Great question! I have intentionally not implemented the testing of a difference for the MESOR. By creating a rhythmic component to the model using amp_acro(), it doesn't necessarily add a parameter for a difference in MESOR. If you want to add a parameter to the model for a difference in MESOR, this can be done by adding a parameter for the grouping variable.

# this model has a parameter for the grouping variable ("X")
cglmm(
    vit_d ~ X + amp_acro(time, group = "X", period = 12),
    data = vitamind
  )
autoplot(object)

# but this one doesn't, meaning that both groups will be assumed to share the same MESOR
cglmm(
    vit_d ~ amp_acro(time, group = "X", period = 12),
    data = vitamind
  )
autoplot(object)

If you've added a parameter for the group to the model, then you can test to see if this is different from the other in the same way that you would if it were a normal regression model.

In this case, since we have two groups, we can check to see whether the X parameter estimate is different from zero and interpret the confidence interval (which represents the difference in the X=1 versus X=0... The MESOR for the X=0 group is estimated by the intercept term).

object <- cglmm(
  vit_d ~ X + amp_acro(time, group = "X", period = 12),
  data = vitamind
)
summary(object)

Hopefully this clarifies why I haven't incorporated it into test_cosinor() and helps you test differences in your model. If you have multiple levels to your grouping variable and things are a bit more complex than this example, then you can access the underlying glmmTMB model from cglmm with object$fit and use this on functions from other packages (maybe {marginaleffects} or {emmeans}).

Thanks!

Rex

[MilaSMayor]

Many thanks Rex. I did your suggestion, I applied emmeans to my model$fit, and I got the comparisons. Thanks again. :)