I am looking at the relationship between age and DNA methylation in teeth samples. We collected three tissues from every tooth (dentin, cementum and pulp), but were unable to get measurements from every sample as you can see in the dataset below. The "Sample" variable contains the numbers of the teeth and "Tissue" represents the repeated measures. "ED1" is the response variable.
Sample Tissue Tooth Condition Sex Age ED1
<fctr> <fctr> <fctr> <fctr> <fctr> <dbl> <dbl>
1 1 Dentin Bicuspid Intact F 15.44422 53
2 1 Pulp Bicuspid Intact F 15.44422 53
3 2 Dentin Bicuspid Intact F 15.45243 56
4 2 Pulp Bicuspid Intact F 15.45243 64
5 3 Dentin Molar Intact F 28.04928 50
6 3 Pulp Molar Intact F 28.04928 52
7 4 Cementum Bicuspid Intact M 30.08077 26
8 4 Dentin Bicuspid Intact M 30.08077 64
9 4 Pulp Bicuspid Intact M 30.08077 70
10 5 Cementum Molar Intact M 32.89528 34
# ... with 93 more rows
I wanted to test whether the collected tissue has a significant effect on the methylation measurements. A lot of googling led me to the answer to this question, and constructed my code similarly:
model1 <- lme(ED1 ~ Tissue, data = ANOVA_All, random = ~1|Sample)
But then I realised that I should also include age in this model, so I did. And then I added Tissue in the last argument to make it clear that Tissue is in fact the repeated measure within Sample, and Age is just a fixed effect:
model2 <- lme(ED1 ~ Age + Tissue, data = ANOVA_All, random = ~Tissue|Sample)
> anova(model2)
numDF denDF F-value p-value
(Intercept) 1 51 2110.2166 <.0001
Age 1 48 12.8560 8e-04
Tissue 2 48 87.7027 <.0001
So my first question would be whether this is the correct way to go about that. Changing the "random" argument like that was somewhat of a guess, I am quite bad with constructing syntax so any feedback is very welcome.
We also collected data on sex, type of tooth and condition of the tooth and I would like to include these as well, so I simply expanded the function:
model3 <- lme(ED1 ~ Age + Tissue + Sex + Tooth + Condition, data = ANOVA_All, random = ~Tissue|Sample)
> anova(model3)
numDF denDF F-value p-value
(Intercept) 1 48 2324.9223 <.0001
Age 1 48 12.4313 0.0009
Tissue 2 48 89.3566 <.0001
Sex 1 44 0.0383 0.8458
Tooth 3 44 1.2044 0.3193
Condition 3 44 2.3711 0.0833
But then I noticed that when I switch around the variables, the p-values change, to an extent where an insignificant effect can become significant. In the following model I moved Sex in front of Tissue, and it's p-value became significant:
model4 <- lme(ED1 ~ Age + Sex + Tissue + Tooth + Condition, data = ANOVA_All, random = ~Tissue|Sample)
> anova(model4)
numDF denDF F-value p-value
(Intercept) 1 48 2324.9229 <.0001
Age 1 48 12.4315 0.0009
Sex 1 44 5.0581 0.0296
Tissue 2 48 86.8457 <.0001
Tooth 3 44 1.2044 0.3193
Condition 3 44 2.3710 0.0833
So my second question is if someone could explain to me (in simple terms, preferably) why this happens? Is it because terms are added sequentially and moving Sex to the front allows it to "take credit" for a larger portion of the variance? And considering this, how can I know for sure whether a variable has a significant effect on ED1 or not?
My understanding of statistics and R is very limited and I realize I am trying to do things that are actually going way over my head. But I would greatly appreciate some answers that are not clouded with too many details so even I can understand and implement them.
Thanks in advance!
