If time and Genotype are both categorical predictors as they appear to be, and you are interested in comparing all time/Genotype pairs to each other, then you can just create one interaction variable, and use Tukey contrasts on it:
weights$TimeGeno <- interaction(weigths$Time, weights$Geno)
model <- lme(weight ~ TimeGeno, random = ~1|Animal/time, data=weights)
comp.timegeno <- glht(model, linfct=mcp(TimeGeno="Tukey"))
If you are interested in other contrasts, then you can use the fact that the linfct argument can take a matrix of coefficients for the contrasts - this way you can set up exactly the comparisons you want.
EDIT
There appears some concern in the comments that the model fitted with the TimeGeno predictor is different from the original model fitted with the Time * Genotype predictor. This is not the case, the models are equivalent. The only difference is in the parametrization of the fixed effects, which is set up to make it easier to use the glht function.
I have used one of the built-in datasets (it has Diet instead of Genotype) to demonstrate that the two approaches have the same likelihood, predicted values, etc:
> # extract a subset of a built-in dataset for the example
> data(BodyWeight)
> ex <- as.data.frame(subset(BodyWeight, Time %in% c(1, 22, 44)))
> ex$Time <- factor(ex$Time)
>
> #create interaction variable
> ex$TimeDiet <- interaction(ex$Time, ex$Diet)
>
> model1 <- lme(weight ~ Time * Diet, random = ~1|Rat/Time, data=ex)
> model2 <- lme(weight ~ TimeDiet, random = ~1|Rat/Time, data=ex)
>
> # the degrees of freedom, AIC, BIC, log-likelihood are all the same
> anova(model1, model2)
Model df AIC BIC logLik
model1 1 12 367.4266 387.3893 -171.7133
model2 2 12 367.4266 387.3893 -171.7133
Warning message:
In anova.lme(model1, model2) :
fitted objects with different fixed effects. REML comparisons are not meaningful.
>
> # the second model collapses the main and interaction effects of the first model
> anova(model1)
numDF denDF F-value p-value
(Intercept) 1 26 1719.5059 <.0001
Time 2 26 28.9986 <.0001
Diet 2 13 85.3659 <.0001
Time:Diet 4 26 1.7610 0.1671
> anova(model2)
numDF denDF F-value p-value
(Intercept) 1 24 1719.5059 <.0001
TimeDiet 8 24 29.4716 <.0001
>
> # they give the same predicted values
> newdata <- expand.grid(Time=levels(ex$Time), Diet=levels(ex$Diet))
> newdata$TimeDiet <- interaction(newdata$Time, newdata$Diet)
> newdata$pred1 <- predict(model1, newdata=newdata, level=0)
> newdata$pred2 <- predict(model2, newdata=newdata, level=0)
> newdata
Time Diet TimeDiet pred1 pred2
1 1 1 1.1 250.625 250.625
2 22 1 22.1 261.875 261.875
3 44 1 44.1 267.250 267.250
4 1 2 1.2 453.750 453.750
5 22 2 22.2 475.000 475.000
6 44 2 44.2 488.750 488.750
7 1 3 1.3 508.750 508.750
8 22 3 22.3 518.250 518.250
9 44 3 44.3 530.000 530.000
The only difference is that what hypotheses are easy to test. For example, in the first model it is easy to test whether the two predictors interact, in the second model there is no explicit test for this. On the other hand, the joint effect of the two predictors is easy to test in the second model, but not the first one. The other hypotheses are testable, it is just more work to set those up.
If I understand your question correctly, you are wondering why you got different p-values from your t-tests when they are carried out as post-hoc tests or as separate tests. But did you control the FWER in the second case (because this is what id done with the step down Sidak-Holm method)? Because, in case of simple t-tests, the t-values won't change, unless you use a different pooling method for computing variance at the denominator, but the p-value of the unprotected tests will be lower than the corrected one.
This is easily seen with Bonferroni adjustment, since we multiply the observed p-value by the number of tests. With step-down methods like Holm-Sidak, the idea is rather to sort the null hypothesis tests by increasing p-values and correct the alpha value with Sidak correction factor in a stepwise manner ($\alpha’ = 1 - (1 - \alpha)^k$, with $k$ the number of possible comparisons, updated after each step). Note that, in contrast to Bonferroni-Holm's method, control of the FWER is only guaranteed when comparisons are independent. A more detailed description of the different kind of correction for multiple comparisons is available here: Pairwise Comparisons in SAS and SPSS.
Best Answer
For the
multcomppackage, see the help page forglht; you want to use the"Tukey"option; this does not actually use the Tukey correction, it just sets up all pairwise comparisons. In the example section there's an example that does exactly what you want.This calculates the estimates and se's for each comparison but doesn't do p-values; for that you need
summary.glht; on the help page for that, notice in particular thetestparameter, that allows you to set the actual test that is run. For doing multiply-adjusted p-values, you use theadjustedfunction for this parameter, and to not multiply-adjust, you usetest=adjusted("none")(which isn't mention specifically in the help page, though it does say it takes any of the methods inp.adjust, which is where you'd findnone.)You can also compute the estimates and se's by hand using matrix multiplication, and then get the p-values however you want; this is what the
glhtfunction is doing behind the scenes. To get the matrices you need to start you'd usecoefandvcov.I didn't put complete code as you say it's for a class project (thanks for being honest, by the way!) and the policy here is to provide helpful hints but not solutions.