Assignment #2
Data Analytics
Seongeun Son
29 8월, 2026
1. 대학 신입생 경쟁률 예측
kor_univ_stat.csv는 교육부에서 발표하는 대학
공시정보로 국내 4년제 대학의 주요 정보를 포함한다.
•
년도 : 공시자료 발표년도 (2022, 2023)
•
학교명 : 학교명
• 설립유형 :
국립/사립/특별법법인
• 지역 : 서울/경기/인천…
(광역자치단체)
• 입학정원 : 단위 (명)
•
졸업생수 : 단위 (명)
• 전임교원수
: 단위 (명)
• 재학생수 : 단위 (명)
•
신입생경쟁률 : 단위 (몇 대 1)
•
신입생충원률 : 단위 (%)
• 취업률
: 단위 (%)
• 외국인학생수 : 단위 (명)
•
전임교원확보율 : 단위 (%)
•
전임교원강의담당비율 : 단위 (%)
•
인당연간장학금 : 단위 (원)
•
연평균등록금 : 단위 (천원)
•
학생1인당교육비 : 단위 (천원)
•
기숙사수용율 : 단위 (%)
•
학생1인당도서자료수 : 단위 (권)
데이터셋의 나머지 정보를 활용하여 신입생경쟁률을 예측하는 linear regression 모델을 수립하고자 한다.
1. 다음 과정을 통해 분석에 필요한 데이터프레임을 준비한다.
- 범주형 변수는 factor로 변환
- 결측치 포함한 행 제거
- 재학생 100명 미만인 대학의 데이터 제거
- 표본 수가 적은 제주, 세종, 울산 지역 대학의 데이터 제거
- “학교명”과 “신입생충원율”은 feature에서 제외
library(car)
library(pls)
library(ggplot2)
library(dplyr)
library(GGally)
library(psych)
library(caret)
library(tidyr)
library(rsample)
library(leaps)
library(glmnet)
library(tidyverse)
library(ggcorrplot)
library(viridis)#차트 테마 설정
my_theme <- function() {
theme_bw() +
theme(
panel.border = element_blank(),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "grey90"),
axis.line = element_line(color = "grey40"),
plot.title = element_text(size = 16, face = "bold"),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
legend.title = element_text(size = 12, face = "bold"),
legend.text = element_text(size = 11)
)
}# grok3의 도움을 받아 작성한 코드
#데이터 불러오기
kor_univ_stat <- read.csv("data/kor_univ_stat.csv")
str(kor_univ_stat)## 'data.frame': 443 obs. of 20 variables:
## $ 년도 : int 2022 2022 2022 2022 2022 2022 2022 2022 2022 2022 ...
## $ 학교명 : chr "가야대학교(고령)-제2캠퍼스" "가야대학교(김해)" "가천대학교" "가톨릭관동대학교" ...
## $ 설립유형 : chr "사립" "사립" "사립" "사립" ...
## $ 지역명 : chr "경북" "경남" "경기" "강원" ...
## $ 입학정원 : int 0 460 4081 1824 108 1574 173 50 200 1491 ...
## $ 졸업생수 : int 0 554 4437 1685 123 1626 179 39 187 1611 ...
## $ 전임교원수 : int 0 81 937 410 23 269 931 15 26 238 ...
## $ 재학생 : int 0 1837 19370 6856 401 6799 913 216 635 6325 ...
## $ 신입생경쟁률 : num 0 4.6 15.4 3.8 4.8 12.7 45.5 1.1 0.7 7.4 ...
## $ 신입생충원율 : num 0 77.1 100 79.8 95.4 99.8 100 98 44.2 99.9 ...
## $ 취업률 : num 0 74.6 67.4 60.8 73.7 60.2 90.1 0 0 60.8 ...
## $ 외국인학생수 : int 0 21 2131 49 0 499 1 6 2 79 ...
## $ 전임교원1인당학생수 : num 0 30 19.6 22.3 25.3 ...
## $ 전임교원확보율 : num 0 86.2 87.1 113 100 ...
## $ 전임교원강의담당비율: num 0 73 63.8 79.2 66.2 62.5 83 66.2 58.2 69.7 ...
## $ 학생1인당연간장학금 : num 0 3710933 3374647 3198333 4926988 ...
## $ 연평균등록금 : num 0 6810 8263 7200 7294 ...
## $ 학생1인당교육비 : num 0 11519 14496 12724 16921 ...
## $ 기숙사수용율 : num 0 34.3 14.9 42.5 45.4 16.7 7.5 102 26 14.5 ...
## $ 학생1인당도서자료수 : num 0 109.7 53.3 95.7 134.2 ...# a) 범주형 변수 factor로 변환
# 설립유형과 지역명을 범주형으로 변환 (숫자가 아닌 변수들)
kor_univ_stat$설립유형 <- as.factor(kor_univ_stat$설립유형)
kor_univ_stat$지역명 <- as.factor(kor_univ_stat$지역명)
# b) 결측치 포함한 행 제거
kor_univ_stat <- na.omit(kor_univ_stat)
# c) 재학생 100명 미만인 대학 데이터 제거
kor_univ_stat <- kor_univ_stat[kor_univ_stat$재학생 >= 100, ]
# d) 제주, 세종, 울산 지역 대학 데이터 제거
kor_univ_stat <- kor_univ_stat[!kor_univ_stat$지역명 %in% c("제주", "세종", "울산"), ]
# e) "학교명"과 "신입생충원율" 열 제거
kor_univ_stat <- kor_univ_stat[, !names(kor_univ_stat) %in% c("학교명", "신입생충원율")]
# 결과 확인
head(kor_univ_stat)## 년도 설립유형 지역명 입학정원 졸업생수 전임교원수 재학생 신입생경쟁률 취업률
## 2 2022 사립 경남 460 554 81 1837 4.6 74.6
## 3 2022 사립 경기 4081 4437 937 19370 15.4 67.4
## 4 2022 사립 강원 1824 1685 410 6856 3.8 60.8
## 5 2022 사립 충북 108 123 23 401 4.8 73.7
## 6 2022 사립 경기 1574 1626 269 6799 12.7 60.2
## 7 2022 사립 서울 173 179 931 913 45.5 90.1
## 외국인학생수 전임교원1인당학생수 전임교원확보율 전임교원강의담당비율
## 2 21 29.99 86.17 73.0
## 3 2131 19.57 87.08 63.8
## 4 49 22.28 112.95 79.2
## 5 0 25.30 100.00 66.2
## 6 499 29.76 70.42 62.5
## 7 1 2.96 479.90 83.0
## 학생1인당연간장학금 연평균등록금 학생1인당교육비 기숙사수용율
## 2 3710933 6809.8 11518.9 34.3
## 3 3374647 8263.2 14496.4 14.9
## 4 3198333 7199.7 12723.8 42.5
## 5 4926988 7293.6 16921.4 45.4
## 6 3287369 6997.6 26257.9 16.7
## 7 2662821 9263.0 0.0 7.5
## 학생1인당도서자료수
## 2 109.7
## 3 53.3
## 4 95.7
## 5 134.2
## 6 110.5
## 7 0.0
2. 변수들 간의 상관 관계를 다양한 그래프를 활용하여
시각화해보고, 이로부터 데이터의 특성을 분석해보자.
# 신입생경쟁률(target) 히스토그램 생성
ggplot(kor_univ_stat, aes(x=신입생경쟁률)) + geom_histogram(color="white", fill = "cornflowerblue") +
labs(title = "신입생 경쟁률 분포", x = "신입생경쟁률", y = "빈도") +
my_theme()
[분석]
● 신입생 경쟁률(target)
히스토그램
경쟁률 분포를 확인해보면, 경쟁률이 5를 전후로
가장 빈도가 높으며, 대부분의 값들이 0에서 20 사이에 집중되어 있는 것을
확인할 수 있다.
# 신입생경쟁률과 취업률, 입학정원 사이의 상관관계 시각화
# Grok3의 도움을 받아 작성한 이상치 제거 dataframe 생성 코드
# 이상치 제거 함수 (IQR 방식)
remove_outliers <- function(x) {
Q1 <- quantile(x, 0.25, na.rm = TRUE)
Q3 <- quantile(x, 0.75, na.rm = TRUE)
IQR <- Q3 - Q1
lower_bound <- Q1 - 1.5 * IQR
upper_bound <- Q3 + 1.5 * IQR
return(x >= lower_bound & x <= upper_bound)
}
# 이상치 제거된 데이터프레임 생성
kor_univ_stat_clean <- kor_univ_stat %>%
filter(remove_outliers(신입생경쟁률) & remove_outliers(취업률) & remove_outliers(입학정원))
# 신입생경쟁률 vs 취업률 산점도 생성
ggplot(kor_univ_stat_clean, aes(x = 신입생경쟁률, y = 취업률)) +
geom_point(color = "cornflowerblue", size = 2) +
geom_smooth(method = "lm", color = "indianred", se = FALSE) +
labs(title = "신입생경쟁률 vs 취업률", x = "신입생경쟁률", y = "취업률") +
my_theme()
# 신입생경쟁률 vs 입학정원 산점도 생성
ggplot(kor_univ_stat_clean, aes(x = 신입생경쟁률, y = 입학정원)) +
geom_point(color = "cornflowerblue", size = 2) +
geom_smooth(method = "lm", color = "indianred", se = FALSE) +
labs(title = "신입생경쟁률 vs 입학정원", x = "신입생경쟁률", y = "입학정원") +
my_theme()
[분석]
● 신입생 경쟁률 vs 취업률
일반적으로 경쟁률이 높은 학교일수록 더 좋은 학교일 가능성이
높다. 그렇다면 좋은 학교를 판단할 수 있는 기준은
무엇일까?
데이터에 포함된 변수 중에서는 취업률이
그 기준이 될 수 있다고 생각되어, 두 변수 간의 관계를 시각화해 보았다.
경쟁률과 취업률 사이에 양의 상관관계가 있을 것이라고 기대했지만,
그래프를 보면 두 변수 사이의 뚜렷한 관계는 보이지 않는다. 회귀선의
기울기가 약간 있는 듯하지만, 회귀선이 없다고 생각하면 상관관계가 거의
0에 가까워 보인다. 따라서 신입생 경쟁률과 취업률 사이에는 의미
있는 상관관계가 없다고 판단하였다.
● 신입생경쟁률
vs 입학정원
두 번째로 경쟁률과 연관이 있어 보이는 변수는
입학정원이었다. 일반적으로 정원이 많을수록 경쟁률은
낮아질 것으로 예상되므로, 음의 상관관계가 있을 것이라 생각했다.
하지만 실제로 시각화된 그래프를 보면, 오히려 양의
상관관계가 나타나는 것으로 보인다. 이유를 생각해보면, 좋은
학교일수록 교육 및 시설 기반이 잘 갖춰져 있어 더 많은 학생을 수용할 수
있는 능력을 갖추고 있고, 그로 인해 입학정원이 많으면서도 경쟁률 역시
높을 수 있다는 해석이 가능하다. 즉, 예상과는 다르게 경쟁률과 입학정원
사이에는 약한 양의 상관관계가 존재한다고 판단할 수
있다.
[분석]
● 입학정원, 졸업생수, 재학생,
외국인학생수 관계
학생 수와 관련된 변수들 간의 상관관계를
파악하고자, 입학정원, 졸업생수, 재학생, 외국인학생수에
대한 상관관계 그래프를 그려보았다. 예상대로 입학정원, 졸업생수,
재학생 간에는 약 0.98의 매우 높은 상관계수가
나타났으며, 외국인학생수 역시 다른 변수들과 약
0.7 수준의 높은 상관관계를 보였다.
이러한 결과를
바탕으로, 해당 변수들은 후속 분석 과정에서 하나의 주성분으로 통합될
가능성이 높다고 판단하였다.
[분석]
● 재학생, 전임교원수,
전임교원1인당학생수, 전임교원확보율 관계
이번에는
전임교원과 관련된 변수들 간의 상관관계를 파악하고자,
재학생, 전임교원수, 전임교원 1인당 학생 수, 전임교원
확보율을 시각화하였다. 예상대로 재학생 수와 전임교원
수 간에는 뚜렷한 상관관계가 나타났으나,
나머지 변수들 간의 상관관계는 크지 않았다.
그
이유는 재학생 수가 증가하면 전임교원 수도 함께 증가하는 경향이 있어,
전임교원 1인당 학생 수는 상대적으로 일정하게 유지되었기 때문으로
해석된다.
또한, 전임교원 확보율 역시 재학생 수 증가율이 일정할
경우 큰 변동이 없기 때문에, 상관관계가 뚜렷하게 나타나지 않은 것으로
판단된다.
- Linear regression 실행 결과를 분석해보자.
- 위 모델에 대해 CV를 수행하고, Training Error, CV Error, Test Error의 크기를 비교해보자.
- Random seed를 바꾸어 Training/Test set의 분할을 다르게 하여 a), b)를 10회 반복해보자. 특히 Training error와 CV Error의 크기에 대한 어떠한 경향을 관찰할 수 있는지 분석해보자.
- 10회 반복한 결과 중 Training Error와 CV Error의 차이가 가장 큰 데이터셋 분할을 아래 문제에서 사용하자.
a) Linear regression 실행 결과를 분석해보자.
# 시드 설정
set.seed(30)
# training set, test set split
split <- initial_split(kor_univ_stat, prop = 0.7, strata = "신입생경쟁률")
kor_univ_stat_train <- training(split)
kor_univ_stat_test <- testing(split)
# pdf 추정 그래프
ggplot() + geom_density(data = kor_univ_stat_train, aes(x = 신입생경쟁률, color = "Train"), linewidth = 1.0) + geom_density(data = kor_univ_stat_test, aes(x = 신입생경쟁률, color = "Test"), linewidth = 1.0) + scale_color_manual(values = c("Train" = "cornflowerblue", "Test" = "indianred")) + labs(title = "Train/Test 신입생 경쟁률 분포 비교", x = "신입생경쟁률", y = "밀도", color = "구분") + my_theme()
# Linear Regression 적용
kor_univ_stat_lm <- lm(신입생경쟁률 ~ ., data = kor_univ_stat_train)
summary(kor_univ_stat_lm)##
## Call:
## lm(formula = 신입생경쟁률 ~ ., data = kor_univ_stat_train)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.8816 -1.4518 -0.2612 1.3575 27.3498
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.995e+02 9.099e+02 0.219 0.826663
## 년도 -9.818e-02 4.499e-01 -0.218 0.827445
## 설립유형사립 -1.069e+00 1.094e+00 -0.977 0.329575
## 설립유형특별법법인 1.217e+01 2.043e+00 5.956 9.13e-09 ***
## 지역명경기 2.740e+00 1.107e+00 2.475 0.014013 *
## 지역명경남 -1.229e+00 1.323e+00 -0.929 0.353671
## 지역명경북 -2.172e-01 1.165e+00 -0.186 0.852267
## 지역명광주 -8.880e-02 1.331e+00 -0.067 0.946867
## 지역명대구 1.358e-01 2.272e+00 0.060 0.952375
## 지역명대전 -1.950e+00 1.361e+00 -1.432 0.153351
## 지역명부산 6.837e-02 1.237e+00 0.055 0.955949
## 지역명서울 5.168e+00 1.189e+00 4.345 2.05e-05 ***
## 지역명인천 1.322e+00 1.667e+00 0.793 0.428559
## 지역명전남 2.521e-01 1.353e+00 0.186 0.852384
## 지역명전북 -9.520e-01 1.330e+00 -0.716 0.474749
## 지역명충남 8.058e-01 1.266e+00 0.636 0.525089
## 지역명충북 1.016e+00 1.217e+00 0.834 0.404846
## 입학정원 -1.093e-03 1.489e-03 -0.734 0.463722
## 졸업생수 -1.603e-04 1.316e-03 -0.122 0.903126
## 전임교원수 -2.220e-03 1.445e-03 -1.537 0.125663
## 재학생 5.583e-04 4.919e-04 1.135 0.257514
## 취업률 6.368e-02 1.866e-02 3.413 0.000754 ***
## 외국인학생수 3.380e-04 5.857e-04 0.577 0.564388
## 전임교원1인당학생수 -3.548e-02 4.136e-02 -0.858 0.391850
## 전임교원확보율 1.442e-02 3.051e-03 4.725 3.92e-06 ***
## 전임교원강의담당비율 8.435e-03 2.263e-02 0.373 0.709690
## 학생1인당연간장학금 -7.301e-07 3.230e-07 -2.260 0.024693 *
## 연평균등록금 3.918e-04 2.660e-04 1.473 0.142049
## 학생1인당교육비 1.540e-05 1.295e-05 1.189 0.235599
## 기숙사수용율 -8.482e-03 9.802e-03 -0.865 0.387754
## 학생1인당도서자료수 -3.713e-03 4.861e-03 -0.764 0.445661
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.329 on 240 degrees of freedom
## Multiple R-squared: 0.6229, Adjusted R-squared: 0.5757
## F-statistic: 13.21 on 30 and 240 DF, p-value: < 2.2e-16
[분석]
통계적으로 유의미한
feature들을 살펴보면, 설립유형특별법법인, 지역(서울, 경기), 취업률,
전임교원확보율, 학생1인당연간장학금이 있다.
이러한 결과는
설립유형특별법법인(KAIST, POSTECH), 수도권에 위치한 대학, 취업률과
전임교원 확보율이 높은 대학일수록 신입생 경쟁률에 긍정적인 영향을
미친다는 점을 의미하며, 충분히 납득이 가는 결과이다.
학생 1인당
연간 장학금이 증가할수록 신입생 경쟁률이 줄어드는 이유를 생각해보면, 잘
되는 학교일수록 굳이 학생들에게 많은 장학금을 지급할 필요가 없기 때문일
수 있다. 이미 학교의 명성이나 다른 경쟁력 있는 요소들로 학생들이
자연스럽게 지원하고, 장학금이 크게 중요한 요소로 작용하지 않는다고 볼 수
있다.
b) 위 모델에 대해 CV를 수행하고, Training Error, CV Error, Test Error의 크기를 비교해보자.
# cv 설정
set.seed(30)
cv <- trainControl(method = "repeatedcv", number = 5)
# Linear Regression 적용
lm_fit <- train(data = kor_univ_stat_train, 신입생경쟁률~., method = "lm", trControl = cv, preProcess = c("center", "scale"))
# training set, test set 예측값 생성
train_pred_lm <- predict(kor_univ_stat_lm, kor_univ_stat_train)
test_pred_lm <- predict(lm_fit, kor_univ_stat_test)
# 결과 출력
cat("Training Error:", RMSE(train_pred_lm, kor_univ_stat_train$신입생경쟁률), "\n")## Training Error: 3.132973## CV Error: 4.129322## Test Error: 3.064756[분석]
각각의 크기를 비교해보면 CV Error
> Training Error > Test Error 이다.
CV, Training Error만
보면 CV가 살짝 더 높아 overfitting 경향이 있나? 라는 생각이 들지만 Test
Error가 제일 작으므로 성능이 괜찮다고 판단된다.
c) Random seed를 바꾸어 Training/Test set의 분할을 다르게
하여 a), b)를 10회 반복해보자.
특히 Training
error와 CV Error의 크기에 대한 어떠한 경향을 관찰할 수 있는지
분석해보자.
# 결과 저장용 column 설정
train_errors <- c()
cv_errors <- c()
test_errors <- c()
# 10번 반복
for (i in 1:10) {
set.seed(i)
# 데이터 분할
split <- initial_split(kor_univ_stat, prop = 0.7, strata = "신입생경쟁률")
kor_univ_stat_train_repeat <- training(split)
kor_univ_stat_test_repeat <- testing(split)
# lm 생성
kor_univ_stat_lm_repeat <- lm(신입생경쟁률 ~ ., data = kor_univ_stat_train_repeat)
# Training Error 생성
train_pred_lm_repeat <- predict(kor_univ_stat_lm_repeat, kor_univ_stat_train_repeat)
rmse_train_repeat <- RMSE(train_pred_lm_repeat, kor_univ_stat_train_repeat$신입생경쟁률)
train_errors <- c(train_errors, rmse_train_repeat)
# CV 설정
set.seed(30)
cv <- trainControl(method = "repeatedcv", number = 5)
# Linear Regression 적용
lm_fit_repeat <- train(data = kor_univ_stat_train_repeat, 신입생경쟁률~., method = "lm", trControl = cv, preProcess = c("center", "scale"))
rmse_cv_repeat <- lm_fit_repeat$results$RMSE
cv_errors <- c(cv_errors, rmse_cv_repeat)
# Test Error 생성
test_pred_lm_repeat <- predict(lm_fit_repeat, kor_univ_stat_test_repeat)
rmse_test_repeat <- RMSE(test_pred_lm_repeat, kor_univ_stat_test_repeat$신입생경쟁률)
test_errors <- c(test_errors, rmse_test_repeat)
}
# 결과 데이터 프레임 생성
results_df <- data.frame(Repeat = 1:10, Train_RMSE = train_errors,
CV_RMSE = cv_errors,Test_RMSE = test_errors)
# 결과 출력
print(results_df)## Repeat Train_RMSE CV_RMSE Test_RMSE
## 1 1 3.170488 3.831948 3.825524
## 2 2 3.206488 6.031362 3.945041
## 3 3 3.110231 4.945344 3.339783
## 4 4 2.378263 2.910006 4.853316
## 5 5 2.378042 2.701075 4.626961
## 6 6 3.179443 3.727167 3.054395
## 7 7 2.950360 3.807972 4.523557
## 8 8 2.542711 3.161206 4.255613
## 9 9 3.089400 4.785383 3.360000
## 10 10 2.622937 3.497443 4.259724# Chat Gpt의 도움을 받아 작성한 그래프 코드
# 길게 변형 (long format)
results_long <- pivot_longer(results_df, cols = c("Train_RMSE", "CV_RMSE"),
names_to = "Error_Type", values_to = "RMSE")
# 선 그래프
ggplot(results_long, aes(x = Repeat, y = RMSE, color = Error_Type)) +
geom_line(linewidth = 1.2) +
geom_point(size = 2) +
scale_x_continuous(breaks = seq(1, 10, 1), labels = as.integer) + # x축을 정수형으로 처리
scale_color_manual(values = c("Train_RMSE" = "cornflowerblue", "CV_RMSE" = "indianred"),
labels = c("Train_RMSE" = "Train RMSE", "CV_RMSE" = "CV RMSE")) +
labs(title = "반복횟수에 따른 Training vs CV RMSE",
x = "반복횟수",
y = "RMSE",
color = "구분") +
my_theme()
[분석]
그래프를 분석한 결과, Training
Error와 CV Error는 대체로 비슷한 경향을 보이며, 반복 횟수에 따라 함께
상승하거나 하락하는 패턴을 보인다. 하지만 모든 반복 횟수에서 Training
Error가 CV Error보다 낮게 유지되며, CV Error가 몇 차례 급격히 튀는
현상을 보이므로 이는 모델이 overfitting되어 검증 데이터에 대한
variance가 커졌다고 판단할 수 있다.
d) 10회 반복한 결과 중 Training Error와 CV Error의 차이가
가장 큰 데이터셋 분할을 아래 문제에서 사용하자.
seed가 2일 때 두 값의 차이가 가장 크므로 이 값을 아래 문제에
적용한다.
4. 강의와 실습에서 다룬 여러 기법을 활용하여 3번의 linear
regression 모델을 개선시켜 본다. 아래 세 가지 기법에
대해 각각 Parameter tuning을 수행하여 best parameter를 선택해
본다. 어떠한 모델이 만들어지는가? 그리고 만들어진
모델의 Training Error와 CV Error를 비교해 보자. 어떠한
경향을 관찰할 수 있는가?
a) Stepwise selection
b) Regularization
c) Principal components regression
최종적으로 하나의 best 모델을 선택해보자.
마지막으로 이 모델에 Test set을 적용해본다.
a)
Stepwise selection
# rmse 차이가 가장 많이 나는 seed 설정
set.seed(2)
# data split
split <- initial_split(kor_univ_stat, prop = 0.7, strata = "신입생경쟁률")
train_data_4 <- training(split)
test_data_4 <- testing(split)
# cv 설정
set.seed(30)
cv <- trainControl(method = "repeatedcv", number = 5, repeats = 5)
# Forward Stepwise Selection 적용
invisible(capture.output(
fwd_model <- train(신입생경쟁률~., data = train_data_4, method = "leapForward",
tuneGrid = data.frame(nvmax = 1:17), trControl = cv)
))
fwd_model## Linear Regression with Forward Selection
##
## 271 samples
## 17 predictor
##
## No pre-processing
## Resampling: Cross-Validated (5 fold, repeated 5 times)
## Summary of sample sizes: 217, 219, 216, 216, 216, 217, ...
## Resampling results across tuning parameters:
##
## nvmax RMSE Rsquared MAE
## 1 5.220872 0.06332852 3.736752
## 2 5.157034 0.10688561 3.633515
## 3 5.067317 0.13300081 3.533805
## 4 5.036387 0.15596232 3.497557
## 5 4.645440 0.26005987 3.227181
## 6 4.471195 0.31639120 3.084548
## 7 4.600701 0.31751161 3.021217
## 8 4.526348 0.34487602 2.919817
## 9 4.521498 0.36340603 2.896437
## 10 4.436288 0.38405834 2.868220
## 11 4.400556 0.38837897 2.850808
## 12 4.371787 0.40926369 2.829053
## 13 4.305021 0.42638516 2.801736
## 14 4.344729 0.44339096 2.765356
## 15 4.324152 0.45088277 2.749093
## 16 4.350850 0.44830782 2.744644
## 17 4.313407 0.45523389 2.736723
##
## RMSE was used to select the optimal model using the smallest value.
## The final value used for the model was nvmax = 13.## Best feature 수: 13
[분석]
RMSE값이 가장 낮은
feature 수 = 13 을 best parameter로 설정
b) Regularization
set.seed(30)
# target 설정
Y <- train_data_4$신입생경쟁률
# feature 설정
X <- model.matrix(신입생경쟁률 ~ ., data = train_data_4)[, -1]
# Ridge 실행
ridge_model <- cv.glmnet(x = X, y = Y, alpha = 0, nfolds = 5)
# Lasso 실행
lasso_model <- cv.glmnet(x = X, y = Y, alpha = 1, nfolds = 5)
# 그래프 그리기 (cv.glmnet 결과를 ggplot으로 재구성)
plot_cv_glmnet <- function(model, title) {
cv_df <- data.frame(
log_lambda = log(model$lambda),
cvm = model$cvm,
cvup = model$cvup,
cvlo = model$cvlo
)
ggplot(cv_df, aes(x = log_lambda, y = cvm)) +
geom_errorbar(aes(ymin = cvlo, ymax = cvup), color = "grey75", width = 0) +
geom_line(color = "cornflowerblue", linewidth = 1.0) +
geom_point(color = "cornflowerblue", size = 2) +
geom_vline(xintercept = log(model$lambda.min), linetype = "dashed", color = "indianred") +
geom_vline(xintercept = log(model$lambda.1se), linetype = "dashed", color = "indianred") +
labs(title = title, x = "log(λ)", y = "CV Mean-Squared Error") +
my_theme()
}
plot_cv_glmnet(ridge_model, "Lambda 값에 따른 CV 오차 (Ridge)")
# best lambda 설정(1se)
best_prm_ridge = ridge_model$lambda.1se
best_prm_lasso = lasso_model$lambda.1se
cat("Best Lambda Ridge:", best_prm_ridge, "\n")## Best Lambda Ridge: 16.0815## Best Lambda lasso: 0.6644898
[분석]
1se rule을 적용하여 best
lambda값 ridge : 16.0815, lasso : 0.6645
c) Principal components regression
# seed 설정
set.seed(30)
# PCR 적용
pcr_model <- train(신입생경쟁률~., data = train_data_4, method = "pcr", trControl = cv,
preProcess = c("center", "scale"), tuneGrid = data.frame(ncomp=1:17))
pcr_model## Principal Component Analysis
##
## 271 samples
## 17 predictor
##
## Pre-processing: centered (33), scaled (33)
## Resampling: Cross-Validated (5 fold, repeated 5 times)
## Summary of sample sizes: 217, 219, 216, 216, 216, 217, ...
## Resampling results across tuning parameters:
##
## ncomp RMSE Rsquared MAE
## 1 4.670053 0.2156703 3.248658
## 2 4.788612 0.1834180 3.339346
## 3 4.763120 0.1961283 3.367437
## 4 4.578451 0.2530715 3.191470
## 5 3.974011 0.4505780 2.605465
## 6 4.038411 0.4418390 2.656590
## 7 3.999901 0.4492088 2.627206
## 8 3.962470 0.4657309 2.612442
## 9 3.941072 0.4720183 2.581830
## 10 3.922276 0.4763013 2.580970
## 11 3.928585 0.4756561 2.588586
## 12 3.932455 0.4769784 2.592760
## 13 3.937230 0.4771692 2.608382
## 14 3.952358 0.4732386 2.620583
## 15 3.936737 0.4773638 2.635582
## 16 3.956618 0.4708150 2.658600
## 17 3.974106 0.4666846 2.665000
##
## RMSE was used to select the optimal model using the smallest value.
## The final value used for the model was ncomp = 10.
## Best 주성분 수: 10
[분석]
분석 결과 Best
주성분 수 = 10 으로 설정
Training Error, CV
Error 비교
# Training Error, CV Error 구하기
# fwd model
fwd_train_pred <- predict(fwd_model, train_data_4)
fwd_train_rmse <- RMSE(fwd_train_pred, train_data_4$신입생경쟁률)
fwd_cv_rmse <- mean(fwd_model$results$RMSE)
# lasso model
lasso_train_pred <- predict(lasso_model, X, s = best_prm_lasso)
lasso_train_rmse <- RMSE(lasso_train_pred, Y)
lasso_cv_rmse <- lasso_model$cvm[lasso_model$lambda == best_prm_lasso]
# ridge model
ridge_train_pred <- predict(ridge_model, X, s = best_prm_ridge)
ridge_train_rmse <- RMSE(ridge_train_pred, Y)
ridge_cv_rmse <- ridge_model$cvm[ridge_model$lambda == best_prm_ridge]
# pcr model
pcr_train_pred <- predict(pcr_model, train_data_4)
pcr_train_rmse <- RMSE(pcr_train_pred, train_data_4$신입생경쟁률)
pcr_cv_rmse <- mean(pcr_model$results$RMSE)
# Chat Gpt의 도움을 받아 작성한 그래프 코드
# 결과 데이터프레임 생성
rmse_df <- data.frame(
Model = rep(c("Forward", "Lasso", "Ridge", "PCR"), each = 2),
Error_Type = rep(c("Train RMSE", "CV RMSE"), times = 4),
RMSE = c(fwd_train_rmse, fwd_cv_rmse,
lasso_train_rmse, lasso_cv_rmse,
ridge_train_rmse, ridge_cv_rmse,
pcr_train_rmse, pcr_cv_rmse)
)
rmse_df## Model Error_Type RMSE
## 1 Forward Train RMSE 3.483363
## 2 Forward CV RMSE 4.593740
## 3 Lasso Train RMSE 3.787021
## 4 Lasso CV RMSE 17.269959
## 5 Ridge Train RMSE 4.165174
## 6 Ridge CV RMSE 18.818810
## 7 PCR Train RMSE 3.600002
## 8 PCR CV RMSE 4.132733# 그래프 그리기
ggplot(rmse_df, aes(x = Model, y = RMSE, fill = Error_Type)) +
geom_bar(stat = "identity", position = "dodge", width = 0.6) +
scale_fill_manual(values = c("Train RMSE" = "cornflowerblue", "CV RMSE" = "indianred")) +
labs(title = "모델별 Training vs CV RMSE 비교",
x = "모델", y = "RMSE", fill = "구분") +
my_theme()
[분석]
Forward와 PCR 모델은 CV
Error가 Train Error보다 약간 크지만 차이가 크지 않아 안정적인 성능을
보이는 것으로 판단된다. 반면, Regularization 기법인 Lasso와 Ridge는 CV
Error와 Train Error의 차이가 크게 나타나 과적합된 것으로 보인다.
최종 모델 선택 및 Test set 적용
네 모델 중
Forward selection과 PCR이 CV Error와 Training Error의 차이가 작아 가장
안정적인 성능을 보였다. 이 중 CV RMSE가 더 낮은 PCR을 최종 모델로
선택하고, Test set에 적용해 실제 일반화 성능을 확인해보자.
# 최종 선택 모델(PCR)을 Test set에 적용
pcr_test_pred <- predict(pcr_model, test_data_4)
pcr_test_rmse <- RMSE(pcr_test_pred, test_data_4$신입생경쟁률)
cat("PCR Training RMSE:", pcr_train_rmse, "\n")## PCR Training RMSE: 3.600002## PCR CV RMSE: 4.132733## PCR Test RMSE: 3.501565# Train / CV / Test RMSE 비교 그래프
pcr_final_df <- data.frame(
Error_Type = factor(c("Train RMSE", "CV RMSE", "Test RMSE"),
levels = c("Train RMSE", "CV RMSE", "Test RMSE")),
RMSE = c(pcr_train_rmse, pcr_cv_rmse, pcr_test_rmse)
)
ggplot(pcr_final_df, aes(x = Error_Type, y = RMSE, fill = Error_Type)) +
geom_bar(stat = "identity", width = 0.5) +
geom_text(aes(label = round(RMSE, 3)), vjust = -0.8, size = 4) +
scale_fill_manual(values = c("Train RMSE" = "cornflowerblue", "CV RMSE" = "indianred", "Test RMSE" = "grey40")) +
scale_y_continuous(expand = expansion(mult = c(0, 0.15))) +
labs(title = "최종 모델(PCR)의 Train / CV / Test RMSE 비교",
x = NULL, y = "RMSE") +
guides(fill = "none") +
my_theme()
[분석]
Test RMSE가 3.502로 Train RMSE, CV
RMSE보다 오히려 낮게 나타났다. 이는 PCR 모델이 과적합 없이 새로운
데이터에도 안정적인 성능을 보인다는 것을 의미한다.
###
2. Lasso Regression의 효과
먼저 아래와 같이 랜덤으로 데이터를 생성하자.
(i) rnorm() 함수를
활용해서 평균이 0, 표준편차가 1인 표준정규분포로부터 크기가 100인 vector
를 생성하고, 평 균이 0, 표준편차가 4인 정규분포로부터 크기가 100인 오차
vector 을 생성한다. 와 을 생성하기 위한 rnorm() 함수에 대해서 동일한
random seed 값을 사용하지 않도록 주의하자.
(ii) 크기가 인
target vector 를 다음 식을 사용하여 생성한다.
Y = 1 − 2X + 3X2 −
4X3 + ϵ
# Chat Gpt의 도움을 받아 작성한 코드
# seed 설정
set.seed(30)
# X:평균이 0, 표준편차가 1인 표준정규분포에서 크기 100만큼 생성
X <- rnorm(100, mean = 0, sd = 1)
# seed 설정
set.seed(23)
# ϵ:평균이 0, 표준편차가 4인 정규분포에서 크기 100만큼 생성
epsilon <- rnorm(100, mean = 0, sd = 4)
# 주어진 식을 사용하여 target Y 생성
Y <- 1 - 2 * X + 3 * X^2 - 4 * X^3 + epsilon
# 결과 확인
head(X) # X 벡터의 처음 6개 값## [1] -1.2885182 -0.3476894 -0.5216288 1.2734732 1.8245206 -1.5113079## [1] 17.8879228 0.4874399 7.0803499 2.2308693 -12.9704017 29.1123548
즉, i번째 관측치 값은 세 가지 feature X, X^2, X^3에 대한
선형식에 오차 를 더한 것과 같다. 위의 선형 관계 식을 모른 채 개의
관측치만 주어졌을 때 이를 추정하기 위한 regression model을 아래의
순서대로 만들어 보자. 즉, 실제 regression coefficient B_0 = 1, B_1 = -2,
B_2 = 3, B_3 = -4를 데이터로부터 추정해야 한다.
X,
X^2, X^3…, X^10의 10개 변수를 feature로, Y를 target으로 설정하자. 이때
feature 변수들과 target변수 사이의 상관관계를 시각화해보자.
# Chat Gpt의 도움을 받아 작성한 코드
# X, X^2, ..., X^10을 생성
X_features <- data.frame(
X = X,
X2 = X^2,
X3 = X^3,
X4 = X^4,
X5 = X^5,
X6 = X^6,
X7 = X^7,
X8 = X^8,
X9 = X^9,
X10 = X^10
)
# 전체 데이터 프레임
data_2 <- cbind(X_features, Y)
# ggcorrplot을 통한 시각화
cor_matrix <- cor(data_2)
ggcorrplot(cor_matrix,
lab = FALSE,
colors = c("indianred", "white", "cornflowerblue")) +
labs(title = "Feature 간 상관관계") +
my_theme() +
theme(panel.grid.major = element_blank(), axis.line = element_blank())
[분석]
상관계수 heat map을 살펴보면
흥미로운 패턴을 확인할 수 있다. X의 지수가 홀수인 경우, 다른 홀수 지수
변수들과 높은 상관관계를 보이며, 짝수인 경우에도 마찬가지로 짝수 지수
변수들과의 상관관계가 높다. 이는 예를 들어, X의 그래프가 와 X^2의
그래프보단 X^3의 그래프와 형태가 유사하다는 점에서 유추할 수 있다.
2. 10개의 feature를 모두 포함하는 linear regression
model을 만들어보자. 통계적으로 유의한 변수가 있는가? regression
coefficient 값을 실제 값과 비교해보자.
##
## Call:
## lm(formula = Y ~ ., data = data_2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.5353 -2.4504 -0.2467 2.6547 7.8305
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.47260 0.82989 2.979 0.00372 **
## X -2.55100 2.03645 -1.253 0.21361
## X2 -3.61920 4.68766 -0.772 0.44212
## X3 -3.95602 4.39356 -0.900 0.37033
## X4 6.48207 6.47636 1.001 0.31960
## X5 -0.35672 2.71802 -0.131 0.89588
## X6 -2.22006 3.20699 -0.692 0.49058
## X7 0.20837 0.60812 0.343 0.73268
## X8 0.29501 0.62799 0.470 0.63967
## X9 -0.02104 0.04243 -0.496 0.62127
## X10 -0.01385 0.04032 -0.344 0.73197
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.77 on 89 degrees of freedom
## Multiple R-squared: 0.9655, Adjusted R-squared: 0.9616
## F-statistic: 248.7 on 10 and 89 DF, p-value: < 2.2e-16# Chat Gpt의 도움을 받아 작성한 regression coefficient 값 비교 코드
# 실제 계수
true_coef <- c(Intercept = 1, X1 = -2, X2 = 3, X3 = -4,
X4 = 0, X5 = 0, X6 = 0, X7 = 0, X8 = 0, X9 = 0, X10 = 0)
# 추정된 계수
est_coef <- coef(lm_model_2)
# 비교 dataframe
coef_compare <- data.frame(
Variable = names(est_coef),
Estimated = round(est_coef, 2),
True = true_coef,
Difference = round(abs(est_coef - true_coef), 3)
)
print(coef_compare)## Variable Estimated True Difference
## (Intercept) (Intercept) 2.47 1 1.473
## X X -2.55 -2 0.551
## X2 X2 -3.62 3 6.619
## X3 X3 -3.96 -4 0.044
## X4 X4 6.48 0 6.482
## X5 X5 -0.36 0 0.357
## X6 X6 -2.22 0 2.220
## X7 X7 0.21 0 0.208
## X8 X8 0.30 0 0.295
## X9 X9 -0.02 0 0.021
## X10 X10 -0.01 0 0.014
[분석]
결과를 보면 통계적으로
유의한 feature는 없다. 실제값과 추정값을 비교했을 때, 특히 X2와 X4에서의
차이가 컸고, 나머지 변수들도 큰 차이는 아니지만 그렇다고 무시할 수 있을
만큼 작지도 않았다.
3. X, X^2, X^3의 3개 변수를
feature로, Y를 target으로 linear regression model을 만들어보자.모든
feature들이 통계적으로 유의한가? regression coefficient 값을 실제 값과
비교해보자.
# X, X^2, X^3 생성
X_features_3 <- data.frame(
X = X,
X2 = X^2,
X3 = X^3)
# feature 3개와 Y를 넣은 데이터 프레임
data_3 <- cbind(X_features_3, Y)
# Linear Regression model
lm_model_3 <- lm(Y~., data = data_3)
summary(lm_model_3)##
## Call:
## lm(formula = Y ~ ., data = data_3)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.8381 -2.4495 -0.2242 2.9823 7.4763
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.4026 0.4787 2.930 0.004238 **
## X -2.1839 0.6180 -3.534 0.000632 ***
## X2 2.8993 0.2555 11.349 < 2e-16 ***
## X3 -3.9342 0.1593 -24.696 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.808 on 96 degrees of freedom
## Multiple R-squared: 0.962, Adjusted R-squared: 0.9608
## F-statistic: 809.9 on 3 and 96 DF, p-value: < 2.2e-16# 실제 계수
true_coef_3 <- c(Intercept = 1, X1 = -2, X2 = 3, X3 = -4)
# 추정된 계수
est_coef_3 <- coef(lm_model_3)
# 비교 dataframe
coef_compare_3 <- data.frame(
Variable = names(est_coef_3),
Estimated = round(est_coef_3, 2),
True = true_coef_3,
Difference = round(abs(est_coef_3 - true_coef_3), 3)
)
print(coef_compare_3)## Variable Estimated True Difference
## (Intercept) (Intercept) 1.40 1 0.403
## X X -2.18 -2 0.184
## X2 X2 2.90 3 0.101
## X3 X3 -3.93 -4 0.066
[분석]
결과를 보면 모든 변수가
통계적으로 유의하다는 것을 확인할 수 있으며, 특히 X2와 X3의 p-value가
2e-16으로 매우 작아 예측값에 매우 유의미한 영향을 준다고 볼 수
있다.
해당 모델에서는 전과 달리 feature 수를 줄여 모델의 복잡도는
줄었지만 오히려 예측 성능은 올라갔다. 이러한 결과를 통해 feature의 수를
잘 선택하는 것이 모델의 성능에 얼마나 영향을 미치는 지 느낄 수 있었다.
4.X, X2, X3…, X10의 10개 변수를 feature로, Y를
target으로 Lasso regression model을 만들어 본다. Cross validation을 통해
합리적인 모델을 찾아보자. 이 모델에는 어떤 변수가 포함되었는가?
regression coefficient 값을 실제 값과 비교해보자. 그리고 결과를 바탕으로
Lasso regression의 효과에 대해서 설명해보자.
# seed 설정
set.seed(30)
# Lasso 실행
lasso_model_4 <- train(Y~., data = data_2, method = "glmnet",
tuneGrid = data.frame(alpha = 1, lambda = seq(0, 2, length.out = 100)),
trControl = trainControl(method = "repeatedcv", number = 5, repeats = 5))
# 그래프 그리기
ggplot(lasso_model_4) + my_theme()
## Best Lambda: 0.3838384# 최적의 lambda 값에 대해 계수 추출
lasso_coef <- coef(lasso_model_4$finalModel, s = best_lambda)
lasso_coef <- as.data.frame(as.matrix(lasso_coef)) # 데이터를 쉽게 처리하기 위해 데이터프레임으로 변환
# 계수 확인
print(lasso_coef)## s=0.3838384
## (Intercept) 1.8249890549
## X -2.6193611458
## X2 2.4649672988
## X3 -3.5405204853
## X4 0.0000000000
## X5 0.0000000000
## X6 0.0000000000
## X7 0.0000000000
## X8 0.0000000000
## X9 -0.0006772809
## X10 0.0000000000# 추정된 계수
est_coef_4 <- lasso_coef[, 1] # 계수 추출 (1열에 계수가 저장됨)
# 계수 비교를 위한 데이터프레임 생성
coef_compare_4 <- data.frame(
Variable = rownames(lasso_coef),
Estimated = est_coef_4,
True = true_coef,
Difference = abs(est_coef_4 - true_coef))
# 비교 결과 출력
print(coef_compare_4)## Variable Estimated True Difference
## Intercept (Intercept) 1.8249890549 1 0.8249890549
## X1 X -2.6193611458 -2 0.6193611458
## X2 X2 2.4649672988 3 0.5350327012
## X3 X3 -3.5405204853 -4 0.4594795147
## X4 X4 0.0000000000 0 0.0000000000
## X5 X5 0.0000000000 0 0.0000000000
## X6 X6 0.0000000000 0 0.0000000000
## X7 X7 0.0000000000 0 0.0000000000
## X8 X8 0.0000000000 0 0.0000000000
## X9 X9 -0.0006772809 0 0.0006772809
## X10 X10 0.0000000000 0 0.0000000000
[분석]
그래프를 보면 최적의 lambda
값은 약 0.4 정도로 보이며, 정확한 값은 0.384이다. 해당 모델의 결과를
보면, X1, X2, X3 외에 X9도 선택되었음을 확인할 수 있다. 추정된 계수와
실제 값을 비교해봤을 때, 변수 3개만 사용한 모델보다는 성능이 조금
떨어지지만, 모든 변수를 사용한 모델보다는 훨씬 더 좋은 성능을
보였다.
따라서 변수의 개수를 모르는 상황에서 모델을 구성해야 할
경우, 단순히 선형 회귀를 적용하기보다는 Lasso 회귀를 활용하는 것이 더
좋아보인다. Lasso 회귀는 불필요한 변수를 자동으로 제거해주기 때문에,
변수 선택이 어려운 상황에서 매우 유용한 방법이라고 판단된다.