# Why do linear regression and ANOVA give different $p$-value in case of considering interaction between variable?

I was trying to fit one time-series data (without replicates) using regression model. The data looks like follows:

> xx.2
value time treat
1  8.788269    1     0
2  7.964719    6     0
3  8.204051   12     0
4  9.041368   24     0
5  8.181555   48     0
6  8.041419   96     0
7  7.992336  144     0
8  7.948658    1     1
9  8.090211    6     1
10 8.031459   12     1
11 8.118308   24     1
12 7.699051   48     1
13 7.537120   96     1
14 7.268570  144     1


Because of lack of replicates, I treat the time as continuous variable. Column "treat" shows the case and control data, respectively.

First, I fit the the model "value = time*treat" with "lm" in R:

summary(lm(value~time*treat,data=xx.2))

Call:
lm(formula = value ~ time * treat, data = xx.2)

Residuals:
Min       1Q   Median       3Q      Max
-0.50627 -0.12345  0.00296  0.04124  0.63785

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept)  8.493476   0.156345  54.325 1.08e-13 ***
time        -0.003748   0.002277  -1.646   0.1307
treat       -0.411271   0.221106  -1.860   0.0925 .
time:treat  -0.001938   0.003220  -0.602   0.5606


The pvalue of time and treat is not significant.

While with anova, I got different results:

 summary(aov(value~time*treat,data=xx.2))
Df Sum Sq Mean Sq F value Pr(>F)
time         1 0.7726  0.7726   8.586 0.0150 *
treat        1 0.8852  0.8852   9.837 0.0106 *
time:treat   1 0.0326  0.0326   0.362 0.5606
Residuals   10 0.8998  0.0900


The pvalue for time and treat changed.

With linear regression, if I am right, it means the time and treat has no significant influence on value, but with ANOVA, it means time and treat has significant influence on value.

Could someone explain to me why there is difference in these two methods, and which one to use?

• You may want to look up the different kinds of sums of squares. Specifically, I believe linear regression returns type III sum of squares, while anova returns a different kind. – assumednormal May 22 '12 at 15:50
• If you save the results of lm and aov you can check they produce identical fits; e.g., compare their residuals with the residuals function or examine their coefficients (the $coefficients slot in both cases). – whuber May 22 '12 at 16:13 ## 4 Answers The fit for lm() and aov() are identical but the reporting is different. The t tests are the marginal impact of the variables in question, given the presence of all the other variables. The F tests are sequential - so they test for the importance of time in the presence of nothing but the intercept, of treat in the presence of nothing but the intercept and time, and of the interaction in the presence of all the above. Assuming you are interested in the significance of treat, I suggest you fit two models, one with, and one without, compare the two by putting both models in anova(), and use that F test. This will test treat and the interaction simultaneously. Consider the following: > xx.2 <- as.data.frame(matrix(c(8.788269, 1, 0, + 7.964719, 6, 0, + 8.204051, 12, 0, + 9.041368, 24, 0, + 8.181555, 48, 0, + 8.041419, 96, 0, + 7.992336, 144, 0, + 7.948658, 1, 1, + 8.090211, 6, 1, + 8.031459, 12, 1, + 8.118308, 24, 1, + 7.699051, 48, 1, + 7.537120, 96, 1, + 7.268570, 144, 1), byrow=T, ncol=3)) > names(xx.2) <- c("value", "time", "treat") > > mod1 <- lm(value~time*treat, data=xx.2) > anova(mod1) Analysis of Variance Table Response: value Df Sum Sq Mean Sq F value Pr(>F) time 1 0.77259 0.77259 8.5858 0.01504 * treat 1 0.88520 0.88520 9.8372 0.01057 * time:treat 1 0.03260 0.03260 0.3623 0.56064 Residuals 10 0.89985 0.08998 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > mod2 <- aov(value~time*treat, data=xx.2) > anova(mod2) Analysis of Variance Table Response: value Df Sum Sq Mean Sq F value Pr(>F) time 1 0.77259 0.77259 8.5858 0.01504 * treat 1 0.88520 0.88520 9.8372 0.01057 * time:treat 1 0.03260 0.03260 0.3623 0.56064 Residuals 10 0.89985 0.08998 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > summary(mod2) Df Sum Sq Mean Sq F value Pr(>F) time 1 0.7726 0.7726 8.586 0.0150 * treat 1 0.8852 0.8852 9.837 0.0106 * time:treat 1 0.0326 0.0326 0.362 0.5606 Residuals 10 0.8998 0.0900 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > summary(mod1) Call: lm(formula = value ~ time * treat, data = xx.2) Residuals: Min 1Q Median 3Q Max -0.50627 -0.12345 0.00296 0.04124 0.63785 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 8.493476 0.156345 54.325 1.08e-13 *** time -0.003748 0.002277 -1.646 0.1307 treat -0.411271 0.221106 -1.860 0.0925 . time:treat -0.001938 0.003220 -0.602 0.5606 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 0.3 on 10 degrees of freedom Multiple R-squared: 0.6526, Adjusted R-squared: 0.5484 F-statistic: 6.262 on 3 and 10 DF, p-value: 0.01154  • Thanks for thorough explanation, it reminds me of the ANCOVA (analysis of covariance). The first step of ANCOVA is to test the interaction between categorical factor and covariate to see if they have identical slope for both condition. It is quite similar to what I did here. In ANCOVA, it gives same pvalue for interaction in t-test and F-test since interaction is the last term in aov. – shao May 23 '12 at 7:26 Peter Ellis' answer is excellent, but there is another point to be made. The$t$-test statistic (and its$p$-value) is a test of whether$\beta = 0$. The$F$-test on the anova() printout is whether the added variable significantly reduces the residual sum of squares. The$t$-test is order-independent, while the$F$-test is not. Hence Peter's suggestion that you try the variables in different orders. It is also possible that variables significant in one test may not be significant in the other (and vice-versa). My sense (and other contributors are welcome to correct me) is that when you're trying to predict phenomena (as in a systems application), you are most interested in reducing variance with the fewest predictors, and therefore want the anova() results. If you are trying to establish the marginal effect of$X$on$y$, however, you will be most concerned with significance of your particular$\beta\$ of interest, and all other variables will just control for alternate explanations your peer reviewers will try to find.

The above two answers are great, but thought I'd add a bit more. Another nugget of information can be gleaned from here.

When you report the lm() results with the interaction term, you're saying something like: "treat 1 is different than treat 0 (beta != 0, p=0.0925), when time is set to the base value of 1". Whereas the anova() results (as previously mentioned) ignore any other variables and concerns itself only with differences in variance.

You can prove this by removing your interaction term and using a simple model with only two main effects (m1):

> m1 = lm(value~time+treat,data=dat)
> summary(m1)

Call:
lm(formula = value ~ time + treat, data = dat)

Residuals:
Min       1Q   Median       3Q      Max
-0.54627 -0.10533 -0.04574  0.11975  0.61528

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept)  8.539293   0.132545  64.426 1.56e-15 ***
time        -0.004717   0.001562  -3.019  0.01168 *
treat       -0.502906   0.155626  -3.232  0.00799 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.2911 on 11 degrees of freedom
Multiple R-squared:   0.64, Adjusted R-squared:  0.5746
F-statistic: 9.778 on 2 and 11 DF,  p-value: 0.003627

> anova(m1)
Analysis of Variance Table

Response: value
Df  Sum Sq Mean Sq F value   Pr(>F)
time       1 0.77259 0.77259  9.1142 0.011677 *
treat      1 0.88520 0.88520 10.4426 0.007994 **
Residuals 11 0.93245 0.08477
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


In this case we see that the reported p-values are the same; that's because in the case of this simpler model,

• This answer unfortunately looks unfinished. Still +1 for the link and for mentioning that the effect is due to different coding schemes. – amoeba Aug 28 '16 at 21:25
• One should also add that summary(lm) and anova(lm) will not always give identical result if there is no interaction term. It just so happens that in these data time and treat are orthogonal and so type I (sequential) and III (marginal) sums of squares yield identical results. – amoeba Jun 6 '17 at 13:29
• The difference has to do with the type pairwise comparisons of cascading models.
• Also, the aov() function has an issue with how it chooses the degrees of freedom. It seems to mix two concepts: 1) the sum of squares from the stepwise comparisons, 2) the degrees of freedom from an overall picture.

PROBLEM REPRODUCTION

> data <- list(value = c (8.788269,7.964719,8.204051,9.041368,8.181555,8.0414149,7.992336,7.948658,8.090211,8.031459,8.118308,7.699051,7.537120,7.268570), time = c(1,6,12,24,48,96,144,1,6,12,24,48,96,144), treat = c(0,0,0,0,0,0,0,1,1,1,1,1,1,1) )
> summary( lm(value ~ treat*time, data=data) )
> summary( aov(value ~ 1 + treat + time + I(treat*time),data=data) )


SOME MODELS USED IN THE EXPLANATION

#all linear models used in the explanation below
> model_0                      <- lm(value ~ 1, data)
> model_time                   <- lm(value ~ 1 + time, data)
> model_treat                  <- lm(value ~ 1 + treat, data)
> model_interaction            <- lm(value ~ 1 + I(treat*time), data)
> model_treat_time             <- lm(value ~ 1 + treat + time, data)
> model_treat_interaction      <- lm(value ~ 1 + treat + I(treat*time), data)
> model_time_interaction       <- lm(value ~ 1 + time + I(treat*time), data)
> model_treat_time_interaction <- lm(value ~ 1 + time + treat + I(treat*time), data)


HOW LM T_TEST WORKS AND RELATES TO F-TEST

# the t-test with the estimator and it's variance, mean square error, is
# related to the F test of pairwise comparison of models by dropping 1
# model parameter

> anova(model_treat_time_interaction, model_time_interaction)

Analysis of Variance Table

Model 1: value ~ 1 + time + treat + I(treat * time)
Model 2: value ~ 1 + time + I(treat * time)
Res.Df     RSS Df Sum of Sq      F  Pr(>F)
1     10 0.89985
2     11 1.21118 -1  -0.31133 3.4598 0.09251 .
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

> anova(model_treat_time_interaction, model_treat_interaction)

Analysis of Variance Table

Model 1: value ~ 1 + time + treat + I(treat * time)
Model 2: value ~ 1 + treat + I(treat * time)
Res.Df     RSS Df Sum of Sq      F Pr(>F)
1     10 0.89985
2     11 1.14374 -1   -0.2439 2.7104 0.1307

> anova(model_treat_time_interaction, model_treat_time)

Analysis of Variance Table

Model 1: value ~ 1 + time + treat + I(treat * time)
Model 2: value ~ 1 + treat + time
Res.Df     RSS Df Sum of Sq      F Pr(>F)
1     10 0.89985
2     11 0.93245 -1 -0.032599 0.3623 0.5606

> # which is the same as
> drop1(model_treat_time_interaction, scope  = ~time+treat+I(treat*time), test="F")

Single term deletions

Model:
value ~ 1 + time + treat + I(treat * time)
Df Sum of Sq     RSS     AIC F value  Pr(>F)
<none>                       0.89985 -30.424
time             1  0.243896 1.14374 -29.067  2.7104 0.13072
treat            1  0.311333 1.21118 -28.264  3.4598 0.09251 .
I(treat * time)  1  0.032599 0.93245 -31.926  0.3623 0.56064
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


HOW AOV WORKS AND CHOOSES DF IN F-TESTS

> #the aov function makes stepwise additions/drops
>
> #first the time, then treat, then the interaction
> anova(model_0, model_time)

Analysis of Variance Table

Model 1: value ~ 1
Model 2: value ~ 1 + time
Res.Df    RSS Df Sum of Sq      F  Pr(>F)
1     13 2.5902
2     12 1.8176  1    0.7726 5.1006 0.04333 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

> anova(model_time, model_treat_time)

Analysis of Variance Table

Model 1: value ~ 1 + time
Model 2: value ~ 1 + treat + time
Res.Df     RSS Df Sum of Sq      F   Pr(>F)
1     12 1.81764
2     11 0.93245  1    0.8852 10.443 0.007994 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

> anova(model_treat_time, model_treat_time_interaction)

Analysis of Variance Table

Model 1: value ~ 1 + treat + time
Model 2: value ~ 1 + time + treat + I(treat * time)
Res.Df     RSS Df Sum of Sq      F Pr(>F)
1     11 0.93245
2     10 0.89985  1  0.032599 0.3623 0.5606

>
> # note that the sum of squares for within model variation is the same
> # but the F values and p-values are not the same because the aov
> # function somehow chooses to use the degrees of freedom in the
> # complete model in all stepwise changes
>


IMPORTANT NOTE

> # Although the p and F values do not exactly match, it is this effect
> # of order and selection of cascading or not in model comparisons.
> # An important note to make is that the comparisons are made by
> # stepwise additions and changing the order of variables has an
> # influence on the outcome!
>
> # Additional note changing the order of 'treat' and 'time' has no
> # effect because they are not correlated

> summary( aov(value ~ 1 + treat + time +I(treat*time), data=data) )

Df Sum Sq Mean Sq F value Pr(>F)
treat            1 0.8852  0.8852   9.837 0.0106 *
time             1 0.7726  0.7726   8.586 0.0150 *
I(treat * time)  1 0.0326  0.0326   0.362 0.5606
Residuals       10 0.8998  0.0900
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

> summary( aov(value ~ 1 + I(treat*time) + treat + time, data=data) )

Df Sum Sq Mean Sq F value  Pr(>F)
I(treat * time)  1 1.3144  1.3144  14.606 0.00336 **
treat            1 0.1321  0.1321   1.469 0.25343
time             1 0.2439  0.2439   2.710 0.13072
Residuals       10 0.8998  0.0900
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

> # This is an often forgotten quirck
> # best is to use manual comparisons such that you know
> # and understand your hypotheses
> # (which is often forgotten in the click and
> #     point anova modelling tools)
> #
> # anova(model1, model2)
> #     or use
> # stepAIC from the MASS library