What is the difference between the standard errors calculated by predict.lm() and summary.lm()

Question

I am trying to calculate standard errors of group means for a two-way-anova. I found two ways to do this (predict.lm(, se = T) and summary.lm():

set.seed(42234)
exmpl <- data.frame(DV = rnorm(40) + rep(3:6 * 10, each = 10), #  Dependent Variable
                    IV1 = factor(rep(LETTERS[1:2], each = 20)), # Independent Variable (Treatment) 1 
                    IV2 = factor(rep(rep(LETTERS[3:4], each = 10), 2))) #  Independent Variable (Treatment) 2

exmpl.lm <- lm(DV ~ IV1 + IV2, data = exmpl) #  Example data was generated without interactions

summary(exmpl.lm)
as.data.frame(predict(exmpl.lm, data.frame(IV1 = c('A', 'B', 'A', 'B'),
                                           IV2 = c('C', 'C', 'D', 'D')), se = T))

The standard errors of some group means differ. I managed to recalculate the standard errors given by predict() with the explanations from Practical Regression and Anova using R from Faraway (section 3.5). I couldn't find any information about the algorithms used by the summary() function. Any ideas? What I find really confusing is that you can change the output by relabelling one factor:

exmpl2 <- exmpl
exmpl2$IV1 <- factor(exmpl2$IV1, levels = LETTERS[2:1])

exmpl2.lm <- lm(DV ~ IV1 + IV2, data = exmpl2)
summary(exmpl2.lm)
summary(exmpl.lm)

In the first example group A-C has a standard error of 0.2013. In the second example the standard error is given as 0.2324. The data of both examples is the same, only the order of the labels of a categorial (not ordinal) variable were changed. How does this influence the statistical model?

Answer 1

The standard errors are the same for both of your models.

exmpl.lm:

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  30.2917     0.2013  150.48   <2e-16 ***
IV1B         19.6148     0.2324   84.39   <2e-16 ***
IV2D         10.1140     0.2324   43.51   <2e-16 ***

exmpl2.lm:

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  49.9065     0.2013  247.93   <2e-16 ***
IV1A        -19.6148     0.2324  -84.39   <2e-16 ***
IV2D         10.1140     0.2324   43.51   <2e-16 ***

The only difference between these two models are the intercept terms. In the first model, the intercept term is an estimate of the mean when IV1 is 'A' and IV2 is 'C':

> mean(exmpl[exmpl$IV1 == 'A' & exmpl$IV2 == 'C', "DV"])
[1] 30.37096

In the second model, the intercept term is an estimate of the mean when IV1 is 'B' and IV2 is 'C':

> mean(exmpl[exmpl$IV1 == 'B' & exmpl$IV2 == 'C', "DV"])
[1] 49.82731

Your standard errors generated by predict are not different from the intercept standard errors. The predictions correspond to the estimated means (intercepts) for the different combinations of IV1 and IV2.

       fit    se.fit df residual.scale
1 30.29173 0.2012963 37        0.73503
2 49.90654 0.2012963 37        0.73503
3 40.40573 0.2012963 37        0.73503
4 60.02053 0.2012963 37        0.73503