1
$\begingroup$

I'm trying to make sure that I'm using mixed models in the correct way. I measured enzimatic activity in 10 specimens (5 female and 5 male), each one measured in three technical replicates. These measurements were made during the four seasons in two different years. I would like to know if Season, Sex and the interaction between Season and Sex have effect on enzimatic activity. So I considered as fixed effects Season (4 levels), Sex (2 levels) and Season*Sex, and as random effects Year.

head(data)

Specimen Replicates Activity Sex Season Year

  1          A   184.78   F    AUT    1
  1          A   179.18   F    AUT    1
  1          A   183.48   F    AUT    1
  2          B   162.77   F    AUT    1
  2          B   161.01   F    AUT    1
  2          B   154.53   F    AUT    1

Iam using these:

model1=lmer(Activity~Season+Sex+Season*Sex+(1|Year))

but interaction wasn't significant, so I simplified to:

model2=lmer(Activity~Season+Sex+(1|Year))

Is it right?

I know that R typically treat all observations having the same set of predictor variables as technical replicates, but Iam not sure if I should to specify the "Specimen" in my model, like this...

model3=lmer(Activity~Season+Sex+(1|Specimen)+(1|Year))

...to R recognize my technical replicates as technical replicates instead outcomes from different subjects.

I know that I could use the mean of my technical replicates as my dependent variable (Activity) but I would like to use all my data in this analysis.

Someone can help me with this issue?

Thanks in advance.

Best regards

D.

$\endgroup$
4
  • 1
    $\begingroup$ (i) 2 levels for Year is far too low for a random effect, change that to a fixed effect. $\endgroup$
    – mkt
    Sep 13 at 18:56
  • 1
    $\begingroup$ (ii) Season+Sex+Season*Sex is redundant, all you need is Season*Sex because that includes both main effects and interactions. R might deal with that intelligently, though, so it may not make a difference. $\endgroup$
    – mkt
    Sep 13 at 18:57
  • $\begingroup$ (iii) What does Specimen mean exactly, and how is it different from Replicate? Can you post your data using dput() or at least a reproducible example illustrating the structure? $\endgroup$
    – mkt
    Sep 13 at 18:58
  • $\begingroup$ Hi mkt, thanks. Specimen is my Subject variable. I will rewrite using Year as fixed effects and post it again using dput(). I will also try to explain in more details my issues on this analysis. Thank you very much. $\endgroup$ Sep 14 at 14:16

1 Answer 1

2
$\begingroup$

In this design, one expects that the functional form of the random effect, using lme4's formula syntax, would be ~ ... + (1|Specimen).

The objective of adjusting for fixed and random effects is that the errors are conditionally independent, and consequently the 95% profile CIs are accurate, the parameter estimates are precise, and the likelihood ratio test is well calibrated.

Adjusting for a fixed effect "costs" one degree of freedom for each non-referent level, when the number of levels becomes too large, the estimation is imprecise; in this case a random effect does better. As the comments point out, Year has two levels, so the fixed effect "costs" barely anything while having much greater predictive power than if it were handled as a random effect.

Conversely, Specimen has 10 levels, which would require 9 covariates to handle using dummy encoding. In an analysis set comprised of only 30 observations, adjusting for these values is not at all tenable. The specimen-effect (which isn't directly handled at all in the proposed models) is a candidate for adjustment as a random effect. One possible approach is to specify specimen as a random intercept, i.e. ~ ... + (1|Specimen). This is a special case of a repeated measures ANOVA. Each cluster of observations within an instance of specimen has a sort of randomly assigned intercept according to a normal $N(0, \sigma_s)$ distribution where $\sigma_s$ is the variance of the random intercept.

The details of fitting mixed models are complicated. The actual value of the "random intercept" is sort of post-estimated as a nuisance parameter, and the level and type of inferences you can make about random effects are completely unlike those you can make about fixed effects. You wouldn't for instance report what the estimated random intercept values are, but it's useful to inspect them as an analyst. Similarly, the method of constructing usual summaries for the fixed effects are completely unlike standard linear models. Small perturbations in fixed effects can completely change the structure and estimate of the random effects, which is why profiling or bootstrapping is often necessary to get p-values and 95% CIs.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.