You can use R stats::anova two ways:
Make up an lm model where you include the following (inclusion order important): covariate mostly correlated with output; next covariate that also strongly correlated with output. You can play with this manually how you want. Fit the model lm_model, and use anova(lm_model).
dat <- data.frame(y = rnorm(100, 0, 1), x1 = rnorm(100, 0, 1), x2 = rnorm(100, 0, 1))
lm_model <- lm(y ~ x1 + x2, data = dat)
summary(lm_model)
anova(lm_model)
> anova(lm_model)
Analysis of Variance Table
Response: y
Df Sum Sq Mean Sq F value Pr(>F)
x1 1 0.026 0.0260 0.0288 0.86550
x2 1 3.828 3.8283 4.2485 0.04196 *
Residuals 97 87.407 0.9011
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
If your model consists of more than 1 coefficient, the anova will make estimation of:
You can just look at Pr(>F) respective to the variable added.
The second way is to feed two models to anova:
lm_model_1 <- lm(y ~ x1, data = dat)
lm_model_2 <- lm(y ~ x1 + x2, data = dat)
anova(lm_model_1, lm_model_2)
> anova(lm_model_1, lm_model_2)
Analysis of Variance Table
Model 1: y ~ x1
Model 2: y ~ x1 + x2
Res.Df RSS Df Sum of Sq F Pr(>F)
1 98 91.235
2 97 87.407 1 3.8283 4.2485 0.04196 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
You get the same p-value for the second term x2 after looking at both approaches to the same thing. If the p-value is lower than precpecified alhpa (significance level) you are fine with inclusion of both predictors.
Also consider this:
lm_model_0 <- lm(y ~ 1, data = dat)
lm_model_1 <- lm(y ~ x1, data = dat)
lm_model_2 <- lm(y ~ x1 + x2, data = dat)
anova(lm_model_0, lm_model_1, lm_model_2)
> anova(lm_model_0, lm_model_1, lm_model_2)
Analysis of Variance Table
Model 1: y ~ 1
Model 2: y ~ x1
Model 3: y ~ x1 + x2
Res.Df RSS Df Sum of Sq F Pr(>F)
1 99 91.261
2 98 91.235 1 0.0260 0.0288 0.86550
3 97 87.407 1 3.8283 4.2485 0.04196 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
- model_0 is based on a constant term only;
- model_1 includes one factor besides the constant;
- model_2 includes one more factor besides what was inside, so it is nested.
If the same result apperas. If you just use model_1 and model_2, then contrasting to a (unseen) model with a constant term is included in anova by default unless you specify the model without intercept.
And this:
> lm_model <- lm(y ~ x1 + x2 - 1, data = dat)
>
> #summary(lm_model)
>
> anova(lm_model)
Analysis of Variance Table
Response: y
Df Sum Sq Mean Sq F value Pr(>F)
x1 1 0.017 0.0166 0.0186 0.89183
x2 1 3.893 3.8926 4.3480 0.03965 *
Residuals 98 87.735 0.8953
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
I removed intercept, and results changed.