http://www-bcf.usc.edu/~gareth/ISL/
Linear Model Selection and Regularization은 교재 6장의 내용입니다.
여기서는 Hitters라는 데이터셋을 이용하고 있습니다.
library(ISLR)
Hitters
names(Hitters)
> names(Hitters)
[1] "AtBat" "Hits" "HmRun" "Runs" "RBI" "Walks" "Years"
[8] "CAtBat" "CHits" "CHmRun" "CRuns" "CRBI" "CWalks" "League"
[15] "Division" "PutOuts" "Assists" "Errors" "Salary" "NewLeague"Hitters는 MLB 1986-1987 시즌 데이터입니다.
Hitters=na.omit(Hitters)여기서는 결측치를 제외하고 분석을 진행합니다.
Linear Model Selection에서는 다음의 3가지 모델 선택 방법에 대해 다룹니다.
1. Best Subset Selection
2. Forward Selection
3. Backward Elimination
먼저 Best Subset Selection은 가능한 모든 모델을 고려하여 그 중에 가장 좋은 것을 찾는 방법입니다. 3가지 방법 중 가장 이상적이지만 p(독립변수의 개수)가 커지면 계산 자체가 불가능하다는 단점이 있습니다.
그래서 나온 방법이 Forward Selection과 Backward Elimination입니다. Forward Selection(FS)은 절편만 포함된 모델에서 출발하여 predictor(독립변수)를 하나씩 추가하면서 가장 좋은 모델을 선택하는 방법입니다.
이와 반대로 Backward Elimination(BE)은 모든 변수가 포함된 모델(Full Model)에서 출발하여 predictor를 하나씩 제거하면서 가장 좋은 모델을 선택하는 방법입니다.
FS의 경우 한 번 추가된 변수는 다시 제거되지 않고, BE의 경우 한 번 제거된 변수는 다시 추가되지 않기 때문에 이 둘을 greedy algorithm이라고 합니다.
R에서 Best Subset Selection을 사용하려면 leaps 라이브러리를 불러와야 합니다.
library(leaps)
regfit.full=regsubsets(Salary~.,Hitters)
summary(regfit.full)leaps 라이브러리에 있는 regsubsets함수를 이용하면 Best Subset Selection을 구현할 수 있습니다.
> summary(regfit.full)
Subset selection object
Call: regsubsets.formula(Salary ~ ., Hitters)
19 Variables (and intercept)
Forced in Forced out
AtBat FALSE FALSE
Hits FALSE FALSE
HmRun FALSE FALSE
Runs FALSE FALSE
RBI FALSE FALSE
Walks FALSE FALSE
Years FALSE FALSE
CAtBat FALSE FALSE
CHits FALSE FALSE
CHmRun FALSE FALSE
CRuns FALSE FALSE
CRBI FALSE FALSE
CWalks FALSE FALSE
LeagueN FALSE FALSE
DivisionW FALSE FALSE
PutOuts FALSE FALSE
Assists FALSE FALSE
Errors FALSE FALSE
NewLeagueN FALSE FALSE
1 subsets of each size up to 8
Selection Algorithm: exhaustive
AtBat Hits HmRun Runs RBI Walks Years CAtBat CHits CHmRun CRuns CRBI CWalks
1 ( 1 ) " " " " " " " " " " " " " " " " " " " " " " "*" " "
2 ( 1 ) " " "*" " " " " " " " " " " " " " " " " " " "*" " "
3 ( 1 ) " " "*" " " " " " " " " " " " " " " " " " " "*" " "
4 ( 1 ) " " "*" " " " " " " " " " " " " " " " " " " "*" " "
5 ( 1 ) "*" "*" " " " " " " " " " " " " " " " " " " "*" " "
6 ( 1 ) "*" "*" " " " " " " "*" " " " " " " " " " " "*" " "
7 ( 1 ) " " "*" " " " " " " "*" " " "*" "*" "*" " " " " " "
8 ( 1 ) "*" "*" " " " " " " "*" " " " " " " "*" "*" " " "*"
LeagueN DivisionW PutOuts Assists Errors NewLeagueN
1 ( 1 ) " " " " " " " " " " " "
2 ( 1 ) " " " " " " " " " " " "
3 ( 1 ) " " " " "*" " " " " " "
4 ( 1 ) " " "*" "*" " " " " " "
5 ( 1 ) " " "*" "*" " " " " " "
6 ( 1 ) " " "*" "*" " " " " " "
7 ( 1 ) " " "*" "*" " " " " " "
8 ( 1 ) " " "*" "*" " " " " " "맨 왼쪽에 있는 숫자는 모델에 포함되는 predictor의 개수를 의미합니다. 별표는 predictor 개수에 따른 모델 집단에서 가장 좋은 모델에 포함된 변수를 표시해주고 있습니다. 예를 들어, predictor 개수가 1인 모델의 경우 CRBI 변수가 포함된 모델이 가장 좋은 모델이라는 것을 알 수 있습니다.
regfit.full=regsubsets(Salary~.,data=Hitters,nvmax=19)
reg.summary=summary(regfit.full)
names(reg.summary)nvmax옵션을 이용하면 모델을 최대 몇개까지 만들지 결정할 수 있습니다. 디폴트값이 8이기 때문에 아무것도 지정하지 않은 경우에 위에 있는 결과에서 알 수 있듯이 predictor가 8개 포함된 모델까지 표시됩니다.
> names(reg.summary) # rsq = R^2 / rss = SSE
[1] "which" "rsq" "rss" "adjr2" "cp" "bic" "outmat" "obj"summary에 어떤 인자들이 있나 확인해보니 모델이 좋은지 아닌지를 판단할 수 있는 기준들이 보입니다.
- rsq = R^2 = 결정계수
- rss = residual sum of square = 잔차제곱합 = SSE
- adjr2 = 수정된 결정계수
- cp = Mallow's Cp
- bic = Bayesian Information Criterion
reg.summary$rsq
> reg.summary$rsq
[1] 0.3214501 0.4252237 0.4514294 0.4754067 0.4908036 0.5087146 0.5141227 0.5285569
[9] 0.5346124 0.5404950 0.5426153 0.5436302 0.5444570 0.5452164 0.5454692 0.5457656
[17] 0.5459518 0.5460945 0.5461159predictor의 수가 증가할수록 R^2도 단조 증가 한다는 것을 알 수 있습니다. 또한 모든 변수를 다 넣어도 R^2 값이 약 55%정도 밖에 되지 않는다는 것도 알 수 있습니다.
모든 모델에 대한 RSS, adjR^2, Cp, BIC를 한꺼번에 그려보면 다음과 같습니다.
par(mfrow=c(2,2))
plot(reg.summary$rss,xlab="Number of Variables",ylab="RSS",type="l")
plot(reg.summary$adjr2,xlab="Number of Variables",ylab="Adjusted RSq",type="l") # Adj R^2
which.max(reg.summary$adjr2)
points(11,reg.summary$adjr2[11], col="red",cex=2,pch=20)
plot(reg.summary$cp,xlab="Number of Variables",ylab="Cp",type='l')
which.min(reg.summary$cp)
points(10,reg.summary$cp[10],col="red",cex=2,pch=20)
plot(reg.summary$bic,xlab="Number of Variables",ylab="BIC",type='l')
which.min(reg.summary$bic)
points(6,reg.summary$bic[6],col="red",cex=2,pch=20)
첫 번째 그래프를 보면 RSS는 더 복잡한 모델일수록 줄어든다는 것을 알 수 있습니다. 따라서 RSS는 predictor 개수가 다른 모델 비교에 사용할 수 없습니다. (반대로 R^2는 더 복잡한 모델일수록 증가합니다. 이 또한 predictor 개수가 다른 모델 비교에 사용할 수 없습니다.)
두 번째 그래프를 보면 adj R^2는 모델의 복잡도가 증가할수록 커졌다가 일정시점이 지나면 다시 작아지는 것을 볼 수 있습니다. 따라서 predictor 개수가 다른 모델 간 비교에 사용할 수 있습니다.
세 번째 그래프는 Mellow's Cp를 기준으로 그린 것입니다. Cp 값이 작을수록 좋은 모델입니다.
네 번째 그래프는 BIC를 기준으로 그린 것입니다. BIC는 predictor 개수가 적은 것을 선택하는 경향이 있습니다. 따라서 과적합(overfitting)의 문제를 피할 수 있지만 중요한 변수를 빠뜨리는 경우가 발생할 수 있습니다.
보통 prediction하는데는 AIC나 Mellow Cp를 사용한다고 합니다. AIC를 사용해서 골랐을 때는 남아있는 factor들이 중요한 것이라고 말할 때 조심해야 합니다. 하지만 BIC를 사용한 경우 남아있는 factor들이 어느 정도는 중요한 것이라고 말할 수 있습니다.
FS와 BE는 다음 코드를 사용해서 구현할 수 있습니다.
regfit.fwd=regsubsets(Salary~., data=Hitters, nvmax=19, method="forward")
summary(regfit.fwd)
regfit.bwd=regsubsets(Salary~., data=Hitters, nvmax=19, method="backward")
summary(regfit.bwd)
coef(regfit.full,7)
coef(regfit.fwd,7)
coef(regfit.bwd,7) # 7은 predictor 개수입니다.regsubsets 함수 말고도 'step'이라는 함수를 사용하는 방법도 있습니다.