library(MASS)
## Warning: package 'MASS' was built under R version 3.5.2
m<-read.csv("Mantel.txt")
We first run lm to obtain the adjusted \(R^2\)
model<-lm(Y~X1+X2+X3,m)
(r123<-summary(model)$adj.r.squared)
## [1] 1
model<-lm(Y~X1+X2,m)
(r12<-summary(model)$adj.r.squared)
## [1] 1
model<-lm(Y~X3,m)
(r3<-summary(model)$adj.r.squared)
## [1] 0.8764847
model<-lm(Y~X2,m)
(r2<-summary(model)$adj.r.squared)
## [1] 0.1641177
model<-lm(Y~X1,m)
(r1<-summary(model)$adj.r.squared)
## [1] 0.1702472
model<-lm(Y~X1+X3,m)
(r13<-summary(model)$adj.r.squared)
## [1] 0.8203869
model<-lm(Y~X2+X3,m)
(r23<-summary(model)$adj.r.squared)
## [1] 0.8205453
Then we calculate all possible AIC
k=2
model<-lm(Y~X1,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(AIC1<-step.model$anova$AIC[1])
## [1] 9.215066
k=2
model<-lm(Y~X1+X2+X3,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(AIC123<-step.model$anova$AIC[1])
## [1] -285.7684
k=2
model<-lm(Y~X1+X2,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(AIC12<-step.model$anova$AIC[1])
## [1] -287.7494
k=2
model<-lm(Y~X2,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(AIC2<-step.model$anova$AIC[1])
## [1] 9.251865
k=2
model<-lm(Y~X2+X3,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(AIC23<-step.model$anova$AIC[1])
## [1] 1.531716
(AIC3<-step.model$anova$AIC[2])
## [1] -0.3087485
k=2
model<-lm(Y~X1+X3,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(AIC13<-step.model$anova$AIC[1])
## [1] 1.536128
Then we calculate all possible BIC. We use the same stepAIC R function, but we use the appropriate penalty
k=log(nrow(m))
model<-lm(Y~X1,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(BIC1<-step.model$anova$AIC[1])
## [1] 8.433942
k=log(nrow(m))
model<-lm(Y~X1+X2+X3,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(BIC123<-step.model$anova$AIC[1])
## [1] -287.3306
k=log(nrow(m))
model<-lm(Y~X1+X2,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(BIC12<-step.model$anova$AIC[1])
## [1] -288.9211
k=log(nrow(m))
model<-lm(Y~X2,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(BIC2<-step.model$anova$AIC[1])
## [1] 8.470741
k=log(nrow(m))
model<-lm(Y~X2+X3,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(BIC23<-step.model$anova$AIC[1])
## [1] 0.3600297
(BIC3<-step.model$anova$AIC[2])
## [1] -1.089873
k=log(nrow(m))
model<-lm(Y~X1+X3,m)
step.model<-stepAIC(model,direction = "backward",trace=FALSE,k=k)
(BIC13<-step.model$anova$AIC[1])
## [1] 0.3644414
$
\[\begin{array}{}
\hline
& R^2 & AIC & BIC \\
X1+X2+X3 & 1 & -285.7683851 & -287.3306334 \\
X1+X2 & 1 & -287.7493779 & -288.9210641 \\
X3 & 0.8764847 & -0.3087485 & -1.0898727 \\
X2 & 0.1641177 & 9.2518654 & 8.4707412 \\
X1 & 0.1702472 & 9.2150658 & 8.4339416 \\
X1+X3 & 0.8203869 & 1.5361276 & 0.3644414 \\
X2+X3 & 0.8205453 & 1.531716 & 0.3600297 \\
\end{array}\]
$
data("mtcars")
mtcars$X=NULL
sum(is.na(mtcars))
## [1] 0
mtcars=na.omit(mtcars)
sum(is.na(mtcars))
## [1] 0
m<-lm(mpg~.,data=mtcars)
AIC(m,k=2)
## [1] 163.7098
stepAIC(m,k=2)
## Start: AIC=70.9
## mpg ~ cyl + disp + hp + drat + wt + qsec + vs + am + gear + carb
##
## Df Sum of Sq RSS AIC
## - cyl 1 0.0799 147.57 68.915
## - vs 1 0.1601 147.66 68.932
## - carb 1 0.4067 147.90 68.986
## - gear 1 1.3531 148.85 69.190
## - drat 1 1.6270 149.12 69.249
## - disp 1 3.9167 151.41 69.736
## - hp 1 6.8399 154.33 70.348
## - qsec 1 8.8641 156.36 70.765
## <none> 147.49 70.898
## - am 1 10.5467 158.04 71.108
## - wt 1 27.0144 174.51 74.280
##
## Step: AIC=68.92
## mpg ~ disp + hp + drat + wt + qsec + vs + am + gear + carb
##
## Df Sum of Sq RSS AIC
## - vs 1 0.2685 147.84 66.973
## - carb 1 0.5201 148.09 67.028
## - gear 1 1.8211 149.40 67.308
## - drat 1 1.9826 149.56 67.342
## - disp 1 3.9009 151.47 67.750
## - hp 1 7.3632 154.94 68.473
## <none> 147.57 68.915
## - qsec 1 10.0933 157.67 69.032
## - am 1 11.8359 159.41 69.384
## - wt 1 27.0280 174.60 72.297
##
## Step: AIC=66.97
## mpg ~ disp + hp + drat + wt + qsec + am + gear + carb
##
## Df Sum of Sq RSS AIC
## - carb 1 0.6855 148.53 65.121
## - gear 1 2.1437 149.99 65.434
## - drat 1 2.2139 150.06 65.449
## - disp 1 3.6467 151.49 65.753
## - hp 1 7.1060 154.95 66.475
## <none> 147.84 66.973
## - am 1 11.5694 159.41 67.384
## - qsec 1 15.6830 163.53 68.200
## - wt 1 27.3799 175.22 70.410
##
## Step: AIC=65.12
## mpg ~ disp + hp + drat + wt + qsec + am + gear
##
## Df Sum of Sq RSS AIC
## - gear 1 1.565 150.09 63.457
## - drat 1 1.932 150.46 63.535
## <none> 148.53 65.121
## - disp 1 10.110 158.64 65.229
## - am 1 12.323 160.85 65.672
## - hp 1 14.826 163.35 66.166
## - qsec 1 26.408 174.94 68.358
## - wt 1 69.127 217.66 75.350
##
## Step: AIC=63.46
## mpg ~ disp + hp + drat + wt + qsec + am
##
## Df Sum of Sq RSS AIC
## - drat 1 3.345 153.44 62.162
## - disp 1 8.545 158.64 63.229
## <none> 150.09 63.457
## - hp 1 13.285 163.38 64.171
## - am 1 20.036 170.13 65.466
## - qsec 1 25.574 175.67 66.491
## - wt 1 67.572 217.66 73.351
##
## Step: AIC=62.16
## mpg ~ disp + hp + wt + qsec + am
##
## Df Sum of Sq RSS AIC
## - disp 1 6.629 160.07 61.515
## <none> 153.44 62.162
## - hp 1 12.572 166.01 62.682
## - qsec 1 26.470 179.91 65.255
## - am 1 32.198 185.63 66.258
## - wt 1 69.043 222.48 72.051
##
## Step: AIC=61.52
## mpg ~ hp + wt + qsec + am
##
## Df Sum of Sq RSS AIC
## - hp 1 9.219 169.29 61.307
## <none> 160.07 61.515
## - qsec 1 20.225 180.29 63.323
## - am 1 25.993 186.06 64.331
## - wt 1 78.494 238.56 72.284
##
## Step: AIC=61.31
## mpg ~ wt + qsec + am
##
## Df Sum of Sq RSS AIC
## <none> 169.29 61.307
## - am 1 26.178 195.46 63.908
## - qsec 1 109.034 278.32 75.217
## - wt 1 183.347 352.63 82.790
##
## Call:
## lm(formula = mpg ~ wt + qsec + am, data = mtcars)
##
## Coefficients:
## (Intercept) wt qsec am
## 9.618 -3.917 1.226 2.936