[R] R 정리 11 - 다중 선형 회귀

11. 다중 회귀 분석 :  독립 변수가 복수 개

 

 - data set

state.x77                          # 미국 50개 주에 대한 dataset
dim(state.x77)                     # 50  8
str(state.x77)
class(state.x77)                   # matrix array
colnames(state.x77)
states <- as.data.frame(state.x77[, c('Murder','Population','Income','Illiteracy', 'Frost')])
str(states)
head(states)

#             Murder Population Income Illiteracy Frost
# Alabama      15.1       3615   3624        2.1    20
# Alaska       11.3        365   6315        1.5   152
# Arizona       7.8       2212   4530        1.8    15
# Arkansas     10.1       2110   3378        1.9    65
# California   10.3      21198   5114        1.1    20
# Colorado      6.8       2541   4884        0.7   166

 

 - 상관계수

# 상관계수
cor(states)
#                Murder Population     Income Illiteracy      Frost
# Murder      1.0000000  0.3436428 -0.2300776  0.7029752 -0.5388834

 

- 다중회귀모델

fit <- lm(Murder ~ Population + Income + Illiteracy + Frost, data=states)
fit
summary(fit) # Multiple R-squared:  0.567

# 값이 5개 이하일 경우 Multiple 사용
# Multiple R-squared과 Adjusted R-squared 차이가 클 경우 data 조정 필요.

par(mfrow=c(2,2))
plot(fit)

---

- 새로운 모델로 예측값 산출

head(states,3)
part_states <- states[1:3, ]
part_states
pre_statues <- edit(part_states) # 새로운 데이터로 값으로 예측
pre_statues
# Murder Population Income Illiteracy Frost
# Alabama   15.1       6000   6000        0.5     5
# Alaska    11.3       1000   1000        8.9   200
# Arizona    7.8       2212   4530        1.8    15

predict(fit, pre_statues)
# Alabama    Alaska   Arizona 
# 5.037489 38.510170  9.487003 

 

 - 예측 변수 선택 : 다양한 고려가 필요

fit1 <- lm(Murder ~., data=states)
summary(fit1) # Adjusted R-squared:  0.5285

fit2 <- lm(Murder ~ Population+Illiteracy, data=states)
summary(fit2) # Adjusted R-squared:  0.5484

 

 - AIC 통계량 : 최적의 모형의 적합도를 위한 독립변수를 추천
                  : 다수의 모델을 비교할 때 사용

                  : 숫자가 작은 모델이 우수한 모델.

AIC(fit1, fit2)
#     df      AIC
# fit1  6 241.6429
# fit2  4 237.6565

---

 - 위의 작업을 자동으로 처리하기 - stepwise regression

  1) backward : 모든 변수를 참여시킨 후 기여도가 낮은 것부터 제거하는 방법

full_model <- lm(Murder ~., data=states)
reduced_model <- step(full_model, direction = "backward")
# Start:  AIC=97.75
# Murder ~ Population + Income + Illiteracy + Frost
# 
# Df Sum of Sq    RSS     AIC
# - Frost       1     0.021 289.19  95.753
# - Income      1     0.057 289.22  95.759
# <none>                    289.17  97.749
# - Population  1    39.238 328.41 102.111
# - Illiteracy  1   144.264 433.43 115.986
# => Frost 제거
# Step:  AIC=95.75
# Murder ~ Population + Income + Illiteracy
# 
# Df Sum of Sq    RSS     AIC
# - Income      1     0.057 289.25  93.763
# <none>                    289.19  95.753
# - Population  1    43.658 332.85 100.783
# - Illiteracy  1   236.196 525.38 123.605
# => Income 제거
# Step:  AIC=93.76
# Murder ~ Population + Illiteracy
# 
# Df Sum of Sq    RSS     AIC
# <none>                    289.25  93.763
# - Population  1    48.517 337.76  99.516
# - Illiteracy  1   299.646 588.89 127.311
summary(reduced_model)

 

  2) forward : 유익한 변수부터 하나씩 추가하는 방법

min_model <- lm(Murder ~1, data=states)
fwd_model <- step(min_model, direction = "forward",
                  scop=Murder ~ Population+Income+Illiteracy+Frost, trace = 1) # trace : log
# Start:  AIC=131.59
# Murder ~ 1
# 
# Df Sum of Sq    RSS     AIC
# + Illiteracy  1    329.98 337.76  99.516
# + Frost       1    193.91 473.84 116.442
# + Population  1     78.85 588.89 127.311
# + Income      1     35.35 632.40 130.875
# <none>                    667.75 131.594
# 
# Step:  AIC=99.52
# Murder ~ Illiteracy
# 
# Df Sum of Sq    RSS     AIC
# + Population  1    48.517 289.25  93.763
# <none>                    337.76  99.516
# + Frost       1     5.387 332.38 100.712
# + Income      1     4.916 332.85 100.783
# 
# Step:  AIC=93.76
# Murder ~ Illiteracy + Population
# 
# Df Sum of Sq    RSS    AIC
# <none>                289.25 93.763
# + Income  1  0.057022 289.19 95.753
# + Frost   1  0.021447 289.22 95.759
summary(fwd_model)

---

- 다중 선형회귀 모델 : iris dataset

head(iris)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species
# 1          5.1         3.5          1.4         0.2  setosa
# 2          4.9         3.0          1.4         0.2  setosa
# 3          4.7         3.2          1.3         0.2  setosa
# 4          4.6         3.1          1.5         0.2  setosa
# 5          5.0         3.6          1.4         0.2  setosa
# 6          5.4         3.9          1.7         0.4  setosa
str(iris)
cor(iris[, -5]) # 범주형 데이터 제외

 

 - train / test로 자료를 분리. 7 : 3
  : 표본의 70% 훈련모델(train), 30% 검증 모델(test)

nrow(iris)   # 150
set.seed(12) # 난수의 시작값 고정
sam_tt <- sample(1:nrow(iris), nrow(iris) * 0.7, replace = F)
sam_tt
NROW(sam_tt)          # 105, NROW() : 배열일 경우 
train <- iris[sam_tt, ]
test <- iris[-sam_tt, ]
dim(train); dim(test)  # 105 5 / 45 5
train
test

 

- 모델 작성

model <- lm(Sepal.Length ~ Sepal.Width+Petal.Length+Petal.Width, data = train) # 학습용 데이터로 모델 작성
options(scipen = 999) # 소수점을 표현
options(scipen = 0)   # E-1
summary(model) # p-value: < 0.00000000000000022

 

 - 선형모델의 만족 조건

# 선형성
aa <- residuals(model)
shapiro.test(aa)              # p-value = 0.8316 > 0.05 이므로 선형성 만족

# 정규성
shapiro.test(model$residuals) # p-value = 0.8316 > 0.05 이므로 정규성 만족

# 독립성(잔차)
install.packages("lmtest")
library(lmtest)
dwtest(model)                 # DW = 2.4964  => 2에 가까울 수록 좋음. 만족.

# 등분산성
install.packages("car")
library(car)
ncvTest(model)                # p = 0.12729 > 0.05 이므로 등분산성 만족

# 다중공선성
model2 <- lm(Sepal.Length ~ Sepal.Width+Petal.Length, data = train)
vif(model2)
vif(model)                    # Petal.Width 14.014160 > 10을 넘으면 다중공선성을 의심해야함.

 

 - 모델평가

pred <- predict(model, test) # test data로 검증
pred[1:3]
# 3        8        9 
# 4.770099 5.048755 4.646646

head(test[,1], 2) # 실제 값              - 4.7 5.0
head(pred, 2)     # 모델에 의한 예측 값  - 4.770099 5.048755