# ================================ << Pulse 데이터 분석 >> ================================
getwd() # 작업위치 확인
setwd("/media/sukjune/Windows/수업/R") # 작업위치 변경
setwd("E:\\R데이터분석\\분석")
load("./Pulse.Rdata") # 변수들을 읽기
save.image("./Pulse.Rdata") # 변수들을 저장
# install.packages("Stat2Data") # Pulse데이터가 있는 패키지를 다운
# install.packages("tidyverse")
# install.packages("reshape2")
# install.packages("caret")
# install.packages("rpart")
# install.packages("e1071")
# install.packages("ipred")
# install.packages("gbm")
# install.packages("nnet")
update.packages(ask = "FALSE") # 기존의 패키지 업데이트
library(nnet) # 다항분포의 로지스틱회귀분석을 실시하기 위한 패키지
library(gbm) # bootstrap을 위한 패키지
library(ipred) # 앙상블예측(bagging)을 하기위한 패키지
library(e1071) # SVM(서포트 벡터 머신)을 위한 패키지
library(rpart) # 트리모형의 패키지
library(caret) # 모형을 생성하는데 도움을 주는 패키지
library(tidyverse)
library(reshape2)
# 휴식기 심박수 기준
# https://ko.wikipedia.org/wiki/%EC%8B%AC%EB%B0%95%EC%88%98
data(Pulse, package = "Stat2Data") # Pulse변수에 데이터 저장
Pulse <- tbl_df(Pulse) # tibble로 자료형을 변환한다.
# ================================ 1. 데이터 탐색 ================================
# ---- (1) 데이터 자료형 확인 및 변형 ----
?Stat2Data::Pulse
str(Pulse) # Pulse 데이터 구조
# Action : 운동후 맥박수(분당 맥박수)
# Rest : 휴식시 맥박수(분당 맥박수)
# Smoke : 1 = 흡연자, 0 = 비흡연자
# Sex : 1 = 여자, 0 = 남자
# Excercise : 1주일 단위의 일반적인 운동시간
# Hgt : 키(단위 : inches)
# Wgt : 무게(단위 : pounds)
# Sex, Smoke는 범주형 변수여야 한다.
# Exercise도 1,2,3 시간밖에 없으므로 범주형으로 가정한다.
Pulse <- mutate(.data = Pulse,
Sex = factor(Sex, labels = c("M", "F")), # M = 여자, F = 남자
Smoke = factor(Smoke, levels = c(1,0), labels = c("Yes", "NO")),
Exercise = factor(Exercise, labels = c("1시간", "2시간", "3시간"))) # Sex와 Smoke를 factor로 변환
any(is.na(Pulse)) # 결측값이 있는지 확인
# 결측값이 존재하지 않는다.
summary(Pulse)# 기초통계량
head(Pulse, n = 10) # 처음 10개의 관측치
Pulse <- mutate(.data = Pulse,
Hgt = round(Hgt*2.54,2), # 1inch = 2.54cm
Wgt = round(Wgt/2.205,2) # 1kg = 2.205pound
)# Wgt와 Hgt를 kg과 cm로 변환한다.(소수 2자리까지 표현)
summary(Pulse)# 기초통계량
head(Pulse, n = 10) # 처음 10개의 관측치
# -- Rest의 범주화 --
# Rest를 기준에 따라 좋음, 보통, 나쁨으로 분류한다
# 데이터의 나이는 18 ~ 25세로 가정한다.
# 남성 : 운동선수 < 56 < =좋음 < 66 <= 보통 < 74 <= 나쁨
# 여성 : 운동선수 < 61 < =좋음 < 69 <= 보통 < 79 <= 나쁨
group_Rest <- function(sex, Rest){
man <- function(rest){ # 남성의 기준
cut(rest, breaks = c(0,56,62,66,70,74,82,150),
labels = c("운동선수", "뛰어남", "좋음","평균이상","평균","평균이하","나쁨"),
levels = c("운동선수", "뛰어남", "좋음","평균이상","평균","평균이하","나쁨"),
include.lowest = FALSE, right = TRUE)
}
woman <- function(rest){ # 여성의 기준
cut(rest, breaks = c(0,61,66,70,74,79,85,150),
labels = c("운동선수", "뛰어남", "좋음","평균이상","평균","평균이하","나쁨"),
levels = c("운동선수", "뛰어남", "좋음","평균이상","평균","평균이하","나쁨"),
include.lowest = FALSE, right = TRUE)
}
ifelse(sex == "M",man(Rest),woman(Rest) ) # 성별에따라 함수를 다르게 사용
}
Pulse <- Pulse %>% mutate(Rest_g = group_Rest(Sex, Rest) %>% factor(., labels = c("운동선수", "뛰어남", "좋음","평균이상","평균","평균이하","나쁨"))) # 새로 만든 변수를 factor로 바꾼다.
Pulse %>% group_by(Rest_g) %>% summarise(N = n()) # 빈도
# ---- (2) 데이터 분포 ----
# -- Active의 분포 --
ggplot(data = Pulse, aes(x = Active)) +
geom_histogram(binwidth = 10, colour = "black", fill = "yellow", aes(y = ..density..), alpha = 0.7) + #각 구간의 크기가 10인 히스토그램
stat_function(fun = dnorm, args = list(mean = mean(Pulse$Active), sd = sd(Pulse$Active)), aes(colour = "정규분포")) + # 정규분포의 곡선
geom_line(aes(colour = "커널밀도곡선"), size = 1.05, stat = "density") + # Active의 커널밀도곡선
ggtitle("Active의 분포") + # 제목
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size =20)) +
scale_color_discrete(name = "")
# Active는 정규분포가 아닌것 같다.
# 정규성 검정
qqnorm(y = Pulse$Active) # qqplot 직선이면 정규분포
qqline(y = Pulse$Active, col = "red") # 정규분포 직선
# 정규분포는 아니라 판단하고 감마분포와 비교한다.
# 감마분포로 qqplot을 그려 검정해본다.
# 감마분포의 2개의 모수는 평균과 분산으로 구할 수 있다.
theta = var(Pulse$Active)/mean(Pulse$Active) # 감마분포의 모수1
alpha = mean(Pulse$Active)/theta # 감마분포의 모수2
ggplot(data = Pulse, aes(x = Active)) + # 그래프
geom_histogram(binwidth = 10, colour = "black", fill = "yellow", aes(y = ..density..), alpha = 0.7) + #각 구간의 크기가 10인 히스토그램
stat_function(fun = dgamma, args = list(shape = alpha, rate = 1/theta), aes(colour = "감마분포")) + # 감마분포의 곡선
geom_line(aes(colour = "커널밀도곡선"), size = 1.05, stat = "density") + # Active의 커널밀도곡선
ggtitle("Active의 분포") + # 제목
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size =20)) +
scale_color_discrete(name = "")
# Active는 감마분포와 비슷하다.
# 감마분포 qqplot그래프 검정
qqplot(x = qgamma(ppoints(500), shape = alpha, rate = 1/theta), y = Pulse$Active) # 산점도
qqline(y = Pulse$Active, distribution = function(x){qgamma(x, shape = alpha, rate = 1/theta)}, col = "red") # 직선
# 감마분포를 따른다고 할 수 있다.
# -- Rest의 분포 --
ggplot(data = Pulse, aes(x = Rest)) +
geom_histogram(binwidth = 10, colour = "black", fill = "yellow", aes(y = ..density..), alpha = 0.7) + # 각 구간의 크기가 10인 히스토그램
geom_line(aes(colour = "blue"), size = 1.05, stat = "density", group = 2) + # 커널밀도곡선
stat_function(fun = dnorm, args = list(mean = mean(Pulse$Rest), sd = sd(Pulse$Rest)), aes(colour = "black"), size = 1.05) + # 정규분포곡선
ggtitle("Rest 와 정규분포") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))+
scale_color_discrete(name = "",labels = c("커널밀도곡선", "정규분포"))
# 정규분포에서 약간 왼쪽으로 치우쳐있다.
# 정규성 qqplot 검정
qqplot(x = qnorm(ppoints(500),mean = mean(Pulse$Rest), sd = sd(Pulse$Rest)), y = Pulse$Rest) # qqplot
qqline(y = Pulse$Rest, distribution = function(x){qnorm(x, mean = mean(Pulse$Rest), sd = sd(Pulse$Rest))}, col = "red") # 정규직선
# 끝의 값들이 직선보다 위에 있어 감마분포랑 비교해 본다.
# 평균과 분산으로 감마분포의 모수를 계산한다.
theta = var(Pulse$Rest)/mean(Pulse$Rest) # 감마분포의 모수1
alpha = mean(Pulse$Rest)/theta # 감마분포의 모수2
# 감마분포 qqplot 검정
qqplot(x = qgamma(ppoints(500), shape = alpha, rate = 1/theta), y = Pulse$Rest)
qqline(y = Pulse$Rest, distribution = function(x){qgamma(x, shape = alpha, rate = 1/theta)}, col = "red")
# 정규분포보다 감마분포에 더 근접해 있다. 따라서 감마분포를 따른다고 할 수 있다.
# -- Active와 Rest간의 분포차이 --
data <- melt(data = Pulse %>% select(Active, Rest), id.vars = NULL, variable.name = "group", value.name = "pulse")
str(data) # Active와 Rest를 스택형태의 데이터로 변환
ggplot(data = data, aes(x = pulse))+
geom_histogram(binwidth = 10, colour = "black", fill = "yellow") +
facet_grid(group ~ .) + # group별로 그래프를 위아래로 나눈다.
ggtitle("Rest 와 Active") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))
# -- Rest와 Active의 관계 --
ggplot(data = Pulse, aes(x = Rest, y = Active))+
geom_point() + # 산점도
ggtitle("Rest 와 Active") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))
# Rest와 Active는 상관관계가 있다.
# -- 키의 분포 --
ggplot(data = Pulse, aes(x = Hgt))+
geom_histogram(fill = "yellow", colour = "black", binwidth = 5) + # 구간의 크가가 5이다.
ggtitle("Hgt") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))
# 170~180의 사람이 많다.
# -- 무게의 분포 --
ggplot(data = Pulse, aes(x = Wgt))+
geom_histogram(fill = "yellow", colour = "black", binwidth = 7) +
ggtitle("Wgt") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))
# 60~80kg의 사람이 대부분이다.
# -- Rest와 Hgt간 관계 --
ggplot(data = Pulse, aes(x = Hgt, y = Rest)) +
geom_point(colour = "blue") +
ggtitle("Rest 와 Hgt") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))
# 상관관계가 있는지 알기 힘들다.
# 키를 150 부터 5cm씩 그룹을만들어 각 그룹의 평균 pulse를 비교한다.
data <- Pulse %>% mutate(Hgt_g = cut(Hgt, breaks = seq(150,200,5),
labels = (seq(150,200,5)[-11] + seq(150,200,5)[-1])/2 )) # Hgt로 그룹을 만든다
data <- data %>% group_by(Hgt_g) %>% summarise(Rest = mean(Rest)) # 각 그룹마다 평균Rest를 구한다.
data <- data %>% mutate(Hgt_g = as.numeric(as.character(Hgt_g))) # Hgt_g을 숫자로 바꿔준다.
head(data)
ggplot(data = data, aes(x = Hgt_g, y = Rest)) + # Rest와 Hgt_g
geom_point(colour = "blue") +
ggtitle("Rest x Hgt_g") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))
# Hgt_g과 Rest는 음의 상관관계가 있어보인다.
# -- Rest와 Wgt간의 관계 --
ggplot(data = Pulse, aes(x = Wgt, y = Rest)) +
geom_point(colour = "blue") +
ggtitle("Rest 와 Hgt") +
theme(plot.title = element_text(hjust = 0.5, size = 20, vjust = 1))
# 상관관계가 있는지 알기 힘들다.
# Wgt를 45에서 120까지 5kg단위로 그룹을 만들어 그룹당 Rest의 평균과 비교한다.
group_cut <- seq(45, 120, 5) # 각 구간의 끝점
length(group_cut)
data <- Pulse %>% mutate(Wgt_g = cut(Wgt, breaks = group_cut,
labels = (group_cut[-16] + group_cut[-1])/2)) # 그룹을 만든다.
data <- data %>% group_by(Wgt_g) %>% summarise(Rest = mean(Rest)) # 각 그룹당 Rest평균계산
ggplot(data = data, aes(x = Wgt_g, y = Rest))+ # Wgt_g와 Rest
geom_point(colour = "blue")+
ggtitle("Rest와 Wgt_g")+
theme(plot.title = element_text(size = 20, hjust = 0.5, vjust = 1))
# Wgt_g와 Rest는 음의 상관관계가 있어보인다.
# -- Rest와 Exercise -- #
ggplot(data = Pulse, aes(x = Exercise, y = Rest, fill = Exercise)) +
geom_boxplot(colour = "black", notch = TRUE) + # 중앙값의 신뢰구간도 표현
ggtitle("Rest와 Exercise") +
theme(plot.title = element_text(size = 20, vjust = 1, hjust = 0.5))
# 이상치가 존재한다.
# 운동시간이 늘어날수록 심박수가 줄어든다.
# -- Rest와 Sex -- #
ggplot(data = Pulse, aes(x = Sex, y = Rest, fill = Sex)) +
geom_boxplot(notch = TRUE) +
ggtitle("Rest와 Sex") +
theme(plot.title = element_text(size = 20, vjust = 1, hjust = 0.5))
# 이상치가 존재한다.
# 여성이 남성보다 심박수가 높다
# 크게 차이나지 않는다.
# -- Rest와 Smoke -- #
ggplot(data = Pulse, aes(x = Smoke, y = Rest, fill = Smoke)) +
geom_boxplot() +
ggtitle("Rest와 Smoke") +
theme(plot.title = element_text(size = 20, vjust = 1, hjust = 0.5))
# 이상치가 존재한다.
# 흡연자가 비흡연자보다 심박수가 높다.
# 크게 차이나지 않는다.
# -- 키와 몸무게의 관계 --
ggplot(data = Pulse, aes(x = Hgt, y = Wgt)) +
geom_point(colour = "blue")
# 키와 몸무게간의 상관관계가 양의 상관관계로 크다고 볼수있다.
# 키와 몸무게를 하나의 변수 체지방량으로 만든다.
# -- 체지방량과 Rest --
# BMI = Wgt/(Hgt/100)^2
# 출저 : https://ko.wikipedia.org/wiki/%EC%B2%B4%EC%A7%88%EB%9F%89_%EC%A7%80%EC%88%98
fat <- function(kg, cm){
kg/(cm/100)^2
}# BMI을 계산한다.
Pulse <- Pulse %>% mutate(fat = fat(Wgt, Hgt)) # Pulse데이터에 fat변수를 추가한다
summary(Pulse$fat)
ggplot(data = Pulse, aes(x = fat)) + # fat의 분포
geom_histogram(fill = "yellow", colour = "black", binwidth = 2) +
ggtitle("fat의 분포") +
theme(plot.title = element_text(size = 20, hjust = 0.5, vjust= 1))
# 왼쪽으로 기울어진 분포이다.
# -- Activet와 fat --
ggplot(data = Pulse, aes(x = fat, y = Active)) +
geom_point(colour = "blue") +
ggtitle("Active와 fat")+
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size =20))
# 상관관계가 있는지 알기 힘들다.
# 상관계수를 계산한다.
cor(Pulse$Active, Pulse$fat) # 0.041
# 0에 가까운 값으로 상관관계가 없다고 판단할 수 있다.
# -- Rest와 fat --
ggplot(data = Pulse, aes(x = fat, y = Rest)) +
geom_point(colour = "blue") +
ggtitle("Rest와 fat")+
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size =20))
# 상관관계를 알기어렵다
# fat를 16에서 36까지 2씩 그룹을 만들어 그룹별 Rest를 비교한다.
group_cut <- seq(14, 36, 2)
length(group_cut)
data <- Pulse %>% mutate(fat_g = cut(fat, breaks = group_cut, labels = (group_cut[-12] + group_cut[-1])/2)) # fat를 그룹화하여 fat_g를 만들었다.
data <- data %>% group_by(fat_g) %>% summarise(Rest = mean(Rest)) # 각 그룹별 평균구하기
data <- data %>% mutate(fat_g = as.numeric(as.character(fat_g))) # factor를 숫자로 변환
ggplot(data = data, aes(x = fat_g, y = Rest)) +
geom_point(colour = "blue") +
ggtitle("Rest와 fat_g") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# BMI는 체질량을 정확하게 계산하지 못하지만 실생활에서 많이 사용하므로 무게와 키를 대신하여 사용한다.
# Rset와 fat_g는 상관관계가 없어 보이지만 분석에 추가하여 영향을 주지않으면 제거한다.
# -- Sex별 Rest와 fat의 관계 --
ggplot(data = Pulse, aes(x = fat, y = Rest)) +
geom_point(size = 1.3) +
facet_grid( . ~ Sex) +
ggtitle("Sex별 Rest와 fat") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 여자들의 체지방량은 낮은곳에 모여있지만 심박수는 더 퍼져있다.
# -- Smoke별 Rest와 fat의 관계 --
ggplot(data = Pulse, aes(x = fat, y = Rest)) +
geom_point(size = 1.3) +
facet_grid(. ~ Smoke) +
ggtitle("Smoke별 Rest와 fat의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 유의한 차이가 없어 보인다.
# 관측치 개수의 차이가 크다.
# -- Exercise별 Rest와 fat의 관계 --
ggplot(data = Pulse, aes(x = fat, y = Rest)) +
geom_point(size = 1.3) +
facet_grid(. ~ Exercise) +
ggtitle("Exercise별 Rest와 fat의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 유의한 차이가 없어 보인다.
# -- Exercise와 Sex별 Rest와 fat의 관계 --
ggplot(data = Pulse, aes(x = fat, y = Rest)) +
geom_point(size = 1.3) +
facet_grid(. ~ Exercise ~ Sex) +
ggtitle("Exercise와 Sex별 Rest와 fat의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 유의한 차이가 없어 보인다.
# 남여모두 운동시간에따라 심박수가 달라진다.
# -- Exercise와 Smoke별 Rest와 fat의 관계 --
ggplot(data = Pulse, aes(x = fat, y = Rest)) +
geom_point(size = 1.3) +
facet_grid(. ~ Exercise ~ Smoke) +
ggtitle("Exercise와 Smoke별 Rest와 fat의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 유의한 차이가 없어 보인다.
# -- 각 그룹별 Rest와 fat의 관계 --
ggplot(data = Pulse, aes(x = fat, y = Rest)) +
geom_point(size = 1.3) +
facet_grid(. ~ Sex + Smoke ~ Exercise) +
ggtitle("각 그룹별 Rest와 fat의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 유의한 차이가 없어 보인다.
# -- Rest와 1/fat의 관계 --
ggplot(data = Pulse, aes(x = 1/fat, y = Rest)) +
geom_point(size = 1.3) +
ggtitle("Rest와 1/fat의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 상관관계가 없어본인다.
# -- 각 그룹별 Rest와 1/fat의 관계 --
ggplot(data = Pulse, aes(x = 1/fat, y = Rest)) +
geom_point(size = 1.3) +
facet_grid(. ~ Sex + Smoke ~ Exercise) +
ggtitle("각 그룹별 Rest와 1/fat의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 특정한 패턴이 보이지 않는다.
# -- Rest와 fat의 상관계수 --
cor(Pulse$Rest, Pulse$fat)
cor(Pulse$Rest, 1/Pulse$fat)
cor(Pulse$Rest, Pulse$Hgt)
cor(Pulse$Rest, Pulse$Wgt)
cor(Pulse$Rest, 1/Pulse$Hgt)
cor(Pulse$Rest, 1/Pulse$Wgt)
# 비만인 사람은 심박수가 높아 Hgt와 Wgt대신 BMI를 사용한다.
# fat는 BMI(체질량 지수)로 비만도를 나타내지만 보디빌더나 운동선수같은 사람들을 비만으로 판단하게되는 경우도 있어 Rest에 영향을 미치지 못하면 제거한다.
# -- fat의 그룹화 --
# 기준 : 저체중 < 18.5 <= 정상 < 25 <= 과체중 < 30 <= 비만
summary(Pulse$fat)
group_fat <- function(fat){
ifelse(fat < 18.5, "저체중",
ifelse(fat < 25, "정상",
ifelse(fat < 30, "과체중", "비만")))
}
Pulse <- Pulse %>% mutate(fat_g = factor(group_fat(fat), levels = c("저체중", "정상", "과체중", "비만") ))
table(Pulse$fat_g) # 비만과 저체중이 개수가 적어 통합을 한다.
# 정상과 과체중으로만 구분한다.
Pulse <- Pulse %>% mutate(fat_g = factor(ifelse(fat_g %in% c("저체중","정상"), "정상", "과체중"),
levels = c("정상", "과체중"))) # 범주를 2개씩 합친다.
table(Pulse$fat_g)
# -- Rest와 fat_g --
ggplot(data = Pulse, aes(x = fat_g, y = Rest, fill = fat_g)) +
geom_boxplot(notch = TRUE) +
ggtitle("Rest와 fat_g") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 큰 차이가 없어보인다.
# Exercie별로 분포를 확인한다.
# -- Exercise별 Rest와 fat_g --
ggplot(data = Pulse, aes(x = fat_g, y = Rest, fill = fat_g)) +
geom_boxplot(notch = TRUE) +
facet_grid(. ~ Exercise) +
ggtitle("Exercise별 Rest와 fat_g") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 큰 차이가 없어보인다.
# fat대신 fat_g를 사용한다.
# -- Active와 Rest_g의 관계 --
ggplot(data = Pulse, aes(x = Rest_g, y = Active, fill = Rest_g)) +
geom_boxplot() +
ggtitle("Active와 Rest_g의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# Rest_g별로 Active가 차이가 난다.
# 평균과 좋음이 크게 차이가 없다.
# -- Sex별 Avtive와 Rest_g의 분포 --
ggplot(data = Pulse, aes(x = Rest_g, y = Active, fill = Rest_g)) +
geom_boxplot() +
facet_grid(. ~ Sex) +
ggtitle("Active와 Rest_g의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 평균과 좋음을 구분할 수 없다.
# -- Exercise별 Avtive와 Rest_g의 분포 --
ggplot(data = Pulse, aes(x = Rest_g, y = Active, fill = Rest_g)) +
geom_boxplot() +
facet_grid(. ~ Exercise) +
ggtitle("Active와 Rest_g의 관계") +
theme(plot.title = element_text(hjust = 0.5, vjust = 1, size = 20))
# 평균과 좋음을 구분할 수 없다.
# Rest를 예측하는데 다른 변수가 필요할 거라 예상한다.
# -- Rest_g의 범주별 상관셩 --
model <- glm(Active ~ Rest_g:Sex, data = Pulse, family = Gamma) # Active의 분포는 Gamma분포였다.
summary(model) # "운동선수"와 다른 변수들 간의 차이가 있다고 할 수 있다.
# "뛰어남"과 비교
data <- Pulse %>% filter(Rest_g != "운동선수") # "운동선수" 제거
model <- glm(Active ~ Rest_g, data = data, family = Gamma) # "운동선수"를 제외하여 "뛰어남"과 비교한다.
summary(model) # "뛰어남"과 "평균", "평균이상"은 차이가 없는것을 알 수 있다.
# "좋음"과 비교
data <- data %>% filter(Rest_g != "뛰어남") # "뛰어남" 제거
model <- glm(Active ~ Rest_g, data = data, family = Gamma) # "뛰어남"를 제외하여 "좋음"과 비교한다.
summary(model) # "좋음"은 "평균", "평균이상"과 차이가 없다.
# "평균이상"과 비교
data <- data %>% filter(Rest_g != "좋음") # "좋음" 제거
model <- glm(Active ~ Rest_g, data = data, family = Gamma) # "좋음"를 제외하여 "평균이상"과 비교한다.
summary(model) # "평균이상"은 "평균"과 차이가 없다.
# "평균"과 비교
data <- data %>% filter(Rest_g != "평균이상") # "평균이상" 제거
model <- glm(Active ~ Rest_g, data = data, family = Gamma) # "평균이상"를 제외하여 "평균"과 비교한다.
summary(model) # "평균"은 남은 변수들과 차이가 있다.
# "평균이하"와 비교
data <- data %>% filter(Rest_g != "평균") # "평균" 제거
model <- glm(Active ~ Rest_g, data = data, family = Gamma) # "평균"을 제외하여 "평균이하"와 비교한다.
summary(model) # "평균이하"은 "나쁨"과 차이가 있다.
# "뛰어남", "좋음", "평균이상"을 하나의 범주 "좋음"으로 합친다.
Pulse <- Pulse %>% mutate(
Rest_g = ifelse(Rest_g %in% c("뛰어남", "좋음", "평균이상"),"좋음", as.character(Rest_g)) %>% # 범주 통합
factor(levels = c("운동선수", "좋음", "평균", "평균이하", "나쁨"))
)
# -- 데이터 분할 --
# 데이터의 80%는 모형을 만드는데 사용하고 20%는 모형을 평가하는데 사용한다.
# Hgt, Wgt, fat는 제외한다.
data <- Pulse %>% select(-Hgt, - Wgt, -fat)
set.seed(20130691) # 매번 같은 값을 얻기위한 코드
train.Pulse <- sample_frac(data, 0.8) # 80%만 추출
test.Pulse <- anti_join(data, train.Pulse) # Pulse로 부터 train.Pulse에 없는 데이터만 추출
str(Pulse)
str(train.Pulse)
str(test.Pulse)
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse")))
# ================================ 2. 회귀분석 ================================
# 반응변수 : Rest
# 설명변수 : Active, Smoke, Sex, Exercise, fat_g
# y의 분포를 정규성으로 가정한 모형과 gamma분포로 가정한 모형을 비교하여 더 잘 예측하는 모형을 사용한다.
train <- train.Pulse %>% select(-Rest_g) # Rest_g제외
train
# -- 변수선택 --
# 널모형에서 1차항과 모든 교호작용이 있는 모형까지 단계별선택법으로 변수를 선택한다.
# AIC를 기준으로 변수를 선택하며 step()를 이용한다.
model <- lm(Rest ~ 1, data = train) # 널 모형(정규성)
step(model, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 전진 선택법
# 첫번재로 Active가 선택되었다.
# Rest ~ Active모형에서 Exercise가 선택되었다.
# Rest ~ Active + Exercise모형에서 Active:Exercise모형이 선택되었다.
# 최종모형은 Rest ~ Active + Exercise + Active:Exercise모형이다.
model <- glm(Rest ~ 1, data = train, family = Gamma) # 널 모형(감마분포)
step(model, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 전진 선택법
# Active + Exercise + Active:Exercise 변수들이 선택되었다.
# -- 모형생성 및 유의성 검정 --
# 정규성 모형
model_reg <- lm(Rest ~ Active + Exercise + Active:Exercise, data = train) # 최선의 모형
anova(model_reg, test = "Chisq") # 각 변수들과 Rest와 관계가 있다.
summary(model_reg) # 각 변수들의 회귀계수의 검정
# Exercise2시간은 회귀계수가 0이라 판단된다.(p-value가 유의수준 0.05보다 크다.)
# 결정계수가 0.5이고 F-통계량의 p-value가 0.05보다 작다.
# Exercise2시간이 유의하지 않지만 교호작용에서 유의하므로 유지시킨다.
layout(matrix(1:4,ncol=2))
plot(model_reg) # 잔차분석
layout(matrix(1,ncol = 1))
# 잔차도 특별한 형태를 가지지않고 qqplot도 정규성을 띈다고 할 수있다.
# 영향력이 높은 값은 제거하지 않는다.
# Gamma모형
model_glm <- glm(Rest ~ Active + Exercise + Active:Exercise, data = train, family = Gamma) # 최선의 모형
anova(model_glm, test = "Chisq") # 각 변수의 유의성 검정
# 모든 변수가 유의한 것을 확인할 수 있다.
summary(model_glm) # Exercise변수가 유의수준 0.05에서 유의하지 않다.
# 따라서 Exercise변수를 제거한다.
model_glm <- update(model_glm, ~ . - Exercise) # Exercise제거
summary(model_glm) # 모든 변수가 유의하다.
# -- 오차 --
coef(model_reg) # 모형의 회귀계수 추정값
pred <- predict(model_reg, newdata = test.Pulse, type = "response") # 정규분포 예측값
mean((test.Pulse$Rest - pred)^2) # 오차
# 오차가 68.93으로 크다.
pred_glm <- predict(model_glm, newdata = test.Pulse, type = "response") # Gamma분포 예측값
mean((test.Pulse$Rest - pred_glm)^2) # 오차
# 77로 정규분포보다 오차가 크다.
error <- c("model_reg" = mean((test.Pulse$Rest - pred)^2)) # 오차저장
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error")))
# ================================ 3. 회귀나무 ================================
# rpart를 사용하여 회귀나무모형을 만든다.
# 회귀나무모형은 최대로 성장시킨후 가지치기를 한다.
train <- train.Pulse %>% select(-Rest_g)
train
# -- 회귀나무모형 생성 --
set.seed(20130691)
model_tree_reg <- rpart(Rest~ ., data = train, cp = -1, method = "anova") # 모형을 최대로 성장시킨다.
plot(model_tree_reg, main = "가지치기 전") # 나무모형그림
text(model_tree_reg)
# -- 가지치기 --
printcp(model_tree_reg) # cp값별 오차테이블
index <- which.min(model_tree_reg$cptable[,"xerror"]) # 오차를 최소로하는 위치를 찾는다.
model_tree_reg$cptable[index,] # 오차가 최소인 행
# 분할은 10번한다. cp 는 0.003으로 계산된다.
cp <- model_tree_reg$cptable[index,"CP"] # 오차를 최소로하는 CP
model_tree_reg <- prune(model_tree_reg, cp = cp) # 가지치기
par(mar = c(2, 0.5, 2, 0.5), oma = c(0, 0, 0, 0))
plot(model_tree_reg, main = "가지치기 후", margin = 0.2) # 나무모형그림
text(model_tree_reg, col = "red", use.n = TRUE)
summary(model_tree_reg) # 모형의 요약
# 변수의 중요도는 Active, Exercise순으로 Active가 영향을 많이준다.
# 가장 처음 분류하는 기준은 Active이다.
model_tree_reg # 회귀나무모형
# 대부분 Active와 Exercise로 구분한다.
# -- 오차 --
pred <- predict(model_tree_reg, newdata = test.Pulse, type = "vector")
mean((test.Pulse$Rest - mean(test.Pulse$Rest))^2) # 구분전의 오차
mean((test.Pulse$Rest - pred)^2) # 오차제곱평균
# 63.27로 회귀분석보다 오차가 작다.
error <- c(error, "model_tree_reg" = mean((test.Pulse$Rest - pred)^2)) # 오차추가
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error")))
# ================================ 4. 회귀분석(부트스트랩) ================================
train <- train.Pulse %>% select(-Rest_g)
train
# train데이터에서 150개의 데이터를 반복랜덤추출을 하여 모형을 만든다.
# 모형은 회귀분석에서 사용했던 step()을 이용해서 변수를 단계별 선택한다.
# 선택된 변수로 모형을 만들고 예측값을 계산해 각 모형의 예측값들의 평균을 구한다.
# 반복은 100회를 한다.
bootstrap <- function(){
random_index <- sample(1:nrow(train), replace = TRUE, size = 150) # 행번호 추출
temp <- train[random_index,] # 데이터 추출
model <- lm(Rest ~ 1,data = temp) # 널모형
# 변수선택
var_step <- step(model, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 단계별 선택
fm <- formula(var_step) # 선택된 모형의 함수식을 반환
model <- lm(fm, data = temp) # 선택된 모형을 생성
return(model) # 모형 반환
}# function end
# train데이터에서 표본을 추출하여 모형을 만든다.
set.seed(20130691) # 시드설정
pred <- sapply(1:100, function(x){ # 모형을 100개 만든든다.
model <- bootstrap() # 모형생성
return(predict(model, newdata = test.Pulse, type = "response")) # 예측값 계산
}) # warning은 무시한다
pred <-apply(pred,1,mean) # 각 행의 평균을 계산한다.
mean((test.Pulse$Rest - pred)^2) # 오차
# 65.61로 회귀분석보다 오차가 작다
error <- c(error, "reg_bootstrap" = mean((test.Pulse$Rest - pred)^2)) # 오차저장
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error")))
# ================================ 5. 회귀나무(bootstrap) ================================
train <- train.Pulse %>% select(-Rest_g)
train
# 회귀분석과 마찬가지로 150개의 데이터를 복원추출로 추출한 후 모형을 생성한다.
# 회귀나무모형은 가지치기를 실시한다.
# 결과값을 평균하여 에측값으로 사용한다.
bootstrap <- function(){ # 모형을 생성하는 함수
random <- sample(1:nrow(train), size = 150, replace = TRUE) # 복원추출
temp <- train[random,] # 데이터 추출
model <- rpart(Rest ~ ., data = temp, cp = -1) # 최초의 모형
index <- which.min(model$cptable[,"xerror"]) # 최소오차의 위치
cp <- model$cptable[index, "CP"] # cp값 추출
model <- prune(model, cp = cp) # 가지치기
return(model)
}# function end
set.seed(20130691)
pred <- sapply(1:100, function(y){
predict(bootstrap(), newdata = test.Pulse)
}) # 예측값 계산
pred <- apply(pred, 1, mean) # 예측값을 평균한다.
mean((test.Pulse$Rest - pred)^2) # 오차
# 53.98로 회귀분석보다 오차가 적다.
error <- c(error, "tree_reg_bootstrap"= mean((test.Pulse$Rest - pred)^2)) # 오차추가
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error")))
# ================================ 6. 신경망모형 ================================
# -----------------------------------------------------------------------#
# svm(Support Vector Machine) 방법
#
# 1990년대 후반 개발된 분석방법으로 회귀분석과 분류분석에 주로사용한다.
# 분류기법중 최상의 모델로 사용된다.
# 단점은 의사결정나무처럼 해석이 불가능하다.
# 따라서 결과를 해석하기 위해서는 의사결정나무를 사용하고 정확한 값을 위해서는 smv를 사용한다.
#
# -----------------------------------------------------------------------#
# -----------------------------------------------------------------------#
# 최대 마진 분류기
#
# 각 관찰값들을 선형경계로 구분하는 방법
# 장점 : 직관적이고 이해하기 쉽다.
# 단점 : 비선형 데이터는 적용할 수 없다.
#
# 초평면(Hyperplane) : 최대 마진 분류기가 선형경계로 사용하는 선, 데이터가 n차원이라면 초평면은 n-1차원을 가진다.
# 분리 초평면(Separating Hyperplane) : 관측치들을 클래스 라벨에 따라 완벽하게 분리하는 초평면들
# 최대 마진 분류기(Maximal Margin Classifier)
# 일반적으로 초평면을 이용하여 분류하게되면 분리 초평면은 무한히 존재한다.
# 따라서 무한한 분리 초평면중 어느 초평면을 사용할지 결정해야한다.
# 선택하는 초평면은 훈련관측치로 부터 가장 멀리 떨어진 최대 마진 초평면을 선택한다.
# 각 훈련데이터들에 초평면까지의 거리를 게산하여 이 값에 따라 최대 마진 초평면을 결정한다.
# 그리고 이때 관측치들에서 초평면까지의 가장 짧은 거리를 마진이라고 한다.
# 따라서 최적 마진 초평면은 분리 초평면으로부터 마진이 가장 큰 초평면이다.
# 검정 관측치가 이 최대 마진 초평면의 어느 쪽에 놓이는가로 분류하는데 이것이 최대 마진 분류기이다.
# 최대 마진 분류기는 분리 초평면이 존재한다면 분류를 수행하기 매우 좋은 방법이다.
# 하지만 어떤 데이터의 경우 초평면이 존재하지 않을 수 있다.
#
# -----------------------------------------------------------------------#
# -----------------------------------------------------------------------#
# 서포트 벡터 분류기(Support Vector Classifier)
#
# 서포트 벡터 분류기는 최대 마진 분류기를 확장한 개념으로 몇몇 훈련 관측치들을 잘못 분류하더라도 나머지 관측치들을 더 잘 분류할 수 있는 방법이다.
# 모든 관측치가 초평면뿐만 아니라 마진의 올바른 쪽(정확히 분류되는 경우)에 있도록 가능한 가장 큰 마진을 찾는 대신에, 일부 관측치들은 마진보다 가까이 있거나 마진의 옳지 않은 쪽(오류가 발생)에 있도록 허용된다.(이때 마진이 소프트하다고 하는 이유는 일부 훈련 관측치들이 마진에 위반되서 분류되기 때문이다.)
# 약간의 오류를 허용함으로 과적합을 방지할 수 있다.
# 일정 수준의 오류를 허용해 주었기 때문에, SVM에서 중요한것은 Cost(오류를 허용해주는 정도)를 정하는 일이다.
# Cost가 큰값을 가지면 margin의 폭은 작아지고, Cost의 값이 작으면 margin의 폭이 커진다.
#
# -----------------------------------------------------------------------#
# -----------------------------------------------------------------------#
# 서포트 벡터 머신(Support Vector Machine)
#
# 서포트 벡터 분류기를 확장하여 비선형의 클래스 경계를 수용하는 분류방법중 하나이다.
# 즉, 선형분류기를 비선형 구조로 변경해서 관측값들을 분류하는 것이다.
# 커널(Kernels)은 클래스들 사이의 비선형 경계를 수용하기 위행 변수공간을 확장하고자 할 때 사용하는 계산기법이다.
# 커널의 차원을 높임으로써 좀 더 다양한 결정경계가 만들어지고 이를 다항식 커널이라 한다.
# 즉, 변수공간에서가 아니라 차수가 d인 다항식들이 관련되는 더 높은 차원의 공간에서 서포트 벡터 분류기를 적합하는것이다.
# 서포트 벡터 분류기가 비선형 커널과 결합될 때 얻어지는 분류기를 서포트 벡터 머신이라고 한다.
# 비선형 커널 이외에도 방사커널(Radial kernel)이 널리 사용된다.
#
# -----------------------------------------------------------------------#
# -----------------------------------------------------------------------#
# 주요 파라미터 튜닝
#
# Cost : 오차 허요에 대한 파라미터
# Gamma : 커널과 관련있는 수치, svm함수에서 가우시안 함수의 표준편차를 조정하는 인자로써, 큰 값을 가질수록 작은 값의 표준편차를 가진다. 때문에 R에서는 svm함수의 kernel인자에 "radial"이 설정될 때 gamma값을 조정한다.(일반 초평면은 kernel인자에 "linear"가 입려되고, 이때는 cost만 조정하면 된다.)
# karnel인자에 "linear"를 입력하면 서포트 벡터 분류기가 실행되고 따라서 cost값으로 소프트 마진을 조정한다.
# kernel인자에 "polynomial"을 입력하면 서포트 벡터 머신이 되고 degree로 차수를 지정한다.(degree가 1을 초과하게 되면 서포트 벡터 머신이 된다.)
# kernel인자에 "radial"을 입력하면 방사커널이 수행되고 gamma와 cost값을 조절해줘야 한다.
#
# -----------------------------------------------------------------------#
# -----------------------------------------------------------------------#
# R에서 SVM수행시 고려사항
# 함수 자체에서 변수 scale(평균 0, 표준편차 1)을 수행해주기 때문에 별도로 scale을 할 필요가 없다.(단, scale이 필요없을 때 scale인자에 FALSE를 입력)
# kernel의 기본값은 radial이므로 cost와 gamma값을 동시에 조정해줘야 한다.
# 선형 초평면을 원한다면 kernel에 "linear"을 입력하고 cost값만 조정한다.
# 다항식 커널을 원하면 kernel에 "polynomial"을 입력하고 degree인자로 차수를 지정한다.
# 데이터 분포를 plot으로 나타냄으로써 선형 초평면, SVM에서 다항식 커널, SVM에서 방사커널 중 무엇을 이용해야 하는지 결정한다.
#
# -----------------------------------------------------------------------#
# SVM(서포트 벡터 머신)을 사용한다.
# kernel은 "radial", "linear", "polynomia"모두 생성하고 오차계산에서 가장작은 오차를 가지는 모형을 선택한다.
train <- train.Pulse %>% select(-Rest_g) # 훈련데이터
train
# -- kernel에 따른 조정이자 튜닝 --
set.seed(20130691) # 랜덤성분을 제어
model1 <- tune.svm(data = train, Rest ~ ., gamma = 2^(-5:0), cost = 2^(0:5), kernel = "radial") # 방사커널
model2 <- tune.svm(data = train, Rest ~ ., cost = 1:40, kernel = "linear") # 서포트 벡터 분류기
model3 <- tune.svm(data = train, Rest ~ ., degree = 1:5, cost = 2^(0:5), kernel = "polynomia") # 다항식 커널
# -- 적정 파라미터 확인 --
model1$best.parameters # gamma = 0.0125, cost = 32
model2$best.parameters # cost = 1
model3$best.parameters # degree = 2, cost = 16
# -- 모형생성 --
model1 <- svm(Rest ~ ., data = train, gamma = 0.0125, cost = 32, kernel = "radial")
model2 <- svm(Rest ~ ., data = train,cost = 1, kernel = "linear")
model3 <- svm(Rest ~ ., data = train, degree = 2, cost = 16, kernel = "polynomia")
# -- 모형 확인 --
summary(model1)
summary(model2)
summary(model3)
# 2개의 모형("linear", "polynomia")의 gamma값이 0.1428571 이다.
# -- 오차 --
pred_1 <- predict(model1, newdata = test.Pulse)
pred_2 <- predict(model2, newdata = test.Pulse)
pred_3 <- predict(model3, newdata = test.Pulse)
mean((test.Pulse$Rest - pred_1)^2) # 방사커널
mean((test.Pulse$Rest - pred_1)^2) # 서포트 벡터 분류기
mean((test.Pulse$Rest - pred_1)^2) # 다항식 커널
# 3개의 모형모두 오차가 60.20이다.
# 따라서 가장 간단한 서포트 벡터 분류기를 모형으로 선택한다.
# -- 변수정리 --
error <- c(error, "model_neural_reg" = mean((test.Pulse$Rest - pred_1)^2))
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error")))
# ================================ 7. 앙상블 예측(bagging) ================================
#-------------------------------------------------------------------#
# bagging
#
# boostrap aggregation의 약자, variance를 줄이는것이 목적(variance가 높으면 데이터셋에 따라 모형이 심하게 변동된다.)
# bagging은 표본에서 데이터를 복원추출하여 각 데이터의 특성을 평규하여 특성을 추정하는 부트스트랩이다.
# 이 방법은 데이터가 소량일 때 유용하다.
# bagging은 decision tree와 달리 prune을 하지않는다.
# OOB(Out-of-Bag) : Bagging에서 모형을 만들때 사용되지 않은 데이터
#
#-------------------------------------------------------------------#
# ipred패키지의 bagging함수를 사용한다.
# nbagg를 100으로 설정하고 정확도를 확인한 후 조정한다.
train <- train.Pulse %>% select(-Rest_g)
train
# -- 모형생성 --
set.seed(20130691) # 시드값 설정
model <- bagging(Rest ~ ., data = train, nbagg = 100) # 100회의 반복으로 bagging을 실행
# nbagg = 부트스트랩의 반복횟수
# rpart의 설정을 사용할 수있다
model$mtrees[1] # bagging의 첫번째 트리모형, 많은 터미널 노드를 가진다.
# 가지치기를 실시하지 않는다.
# -- 오차 --
pred <- predict(model, newdata = test.Pulse) # 예측값
mean((test.Pulse$Rest - pred)^2) # 오차
# 약 54.78으로 오차가 회귀분석보다 작다.
error <- c(error, "tree_reg_bagging" = mean((test.Pulse$Rest - pred)^2)) # 오차추가
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error")))
# ================================ 8. 앙상블 예측(boosting) ================================
#-------------------------------------------------------------------#
# boosting
#
# 다수의 분류기를 생성하는 기법중 하나
# bagging과는 다르게 순차적으로 분류기를 생성
# 생성된 분류기에서 잘못 분류된 관측값이나 훈련에 사용되지 않은 값들은 다음 분류기 생성 데이터 추출때 가중치를 부여
#
#-------------------------------------------------------------------#
#-------------------------------------------------------------------#
# R에서의 인자
#
# n.trees : 생성할 나무의 개수
# interaction.depth : 각 노드에서 뻗어나갈 가지의 개수를 지정(나무의 깊이지정), 이 값이 N이면 2*N+1개의 terminal node가 생성
# shrinkage : Learning Rate로, gradient boosting을 학습하는데 걸리는 시간을 조정, 만약 이값이 크다면 빨리 학습을 수행하지만, 수행도중 중요한점을 놓칠 수 있다. 때문에 이 값이 크면 n.trees인자로 나무를 많이 생성함으로 보완한다. shrinkage가 작고 n.trees가 크면 과적합이 발생 할 수 있다.
# n.minobsinnode : 트리의 terminal node의 최소 관측수, 정수로 입력
# 종속변수가 범주형이면 distributio = "bernoulli"를 입력
#
#-------------------------------------------------------------------#
# gbm패키지의 함수를 사용한다.
# n.tree와 interaction.depth두 변수의 파라미터를 조정해야한다.
train <- train.Pulse %>% select(-Rest_g)
train
# -- 튜닝 --
set.seed(20130691)
control <- caret::trainControl(method = "boot", search = "grid") # 함수 수행을 위한 셋팅
tunegrid <- expand.grid(n.trees = seq(20,100,10), interaction.depth = 1:10, shrinkage = 0.1, n.minobsinnode = 1:10) # gbm인자들을 비교할 목록
gbm_tuning <- caret::train(Rest ~ ., data = train, method = "gbm", tuneGrid = tunegrid, trControl = control) # train함수로 모형생성
print(gbm_tuning) # 결과
# n.trees = 30, interaction.depth = 2, shrinkage = 0.1, n.minobsinnode = 10이 가장 최선의 모형 파라미터이다.
# -- 모형생성 --
set.seed(20130691)
reg_boosting_model <- gbm(Rest ~ ., data = train, n.tree = 30, interaction.depth = 2, shrinkage = 0.1, n.minobsinnode = 10) # 튜닝의 결과를 사용
summary(reg_boosting_model) # 변수별 중요도를 출력한다.
# Active가 가장중요하고 Exercise가 두번째로 중요하다.
# -- 오차(모델의 정확도) --
pred <- predict(reg_boosting_model, newdata = test.Pulse, type = "response", n.trees = 30) # 예측값 계산
mean((test.Pulse$Rest - pred)^2) # 오차
# 57.12 회귀분석보다 오차가 작다.
error <- c(error , "tree_reg_boosting" = mean((test.Pulse$Rest - pred)^2)) # 오차저장
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error")))
# ================================ 9. Rest의 최종결론 ================================
error # 각 모형의 오차
sort(error) # 오름차순으로 정렬
# 앙상블 모형이 가장좋고 그 다음이 svm, 나무모형, 회귀분석순으로 오차가 작다.
# bagging의 모형이 가장 좋은 모형이다.
# bootstrap모형은 bagging과 유사하게 만들었다.
# 10.번부터는 Rest_g을 예측한다.
# ================================ 10. 로지스틱 회귀분석 ================================
#---------------------------------------------------------------------------#
# 일반화 선형모형
#
# 일반화 선형모형(GLM)은 분산이 일정하지 않고 오차가 정규분포를 따르지 않는 경우 사용할 수 있다.
# 비율로 표현되는 가운트 데이터 - 로지스틱 회귀분석
# 비율이 아닌 카운트 데이터 - 카운트에 대한 로그선형 모형
# 이진 반응 변수(예 : 사망 또는 생존)
# 평균과 함께 분산이 선형증가율보다 빠르게 증가하는 생존시간(예 : 감마오차를 갖는 시간데이터)
# 1. 오차구조
# - 포아송 오차 : 카운트 데이터
# - 이항 오차 : 비율 데이터
# - 감마 오차 : 변화의 일정한 계수를 보여주는 데이터
# - 지수 오차 : 생존시간 데이터(생존분석)
#
# 2. 선형예측 모형
# - 예측값은 선형 예측치에서 발생한 값의 변환으로 계산 할 수있다.
#
# 3. 연결함수
# - 선형 예측치와 y의 평균과의 관계인 함수
# - normal -> identiry
# - poisson -> log
# - binomial -> logit
# - Gamma -> reciprocal : R에서는 감마오차만 대문자이다.
#
#---------------------------------------------------------------------------#
# 로지스틱 회귀석은 nnet패키지의 multinom()과 glm()을 모두 사용한다.
# glm()의 경우 "나쁨"," 평균이하", "평균", "좋음", "운동선수"모형을 많들어 예측력이 좋은 4개의 모형으로 5개의 범주를 예측한다.
# multinom()의 인자는 glm()의 인자와 비슷하다.
# glm()의 랜덤성분은 binomial이고 연결은 logit으로 설정한다.
# 변수의 선택은 step()를 이용하여 널 모형에서 최대모형까지 단계별 선택법을 사용합니다.
train <- train.Pulse %>% mutate(Rest_g_1 = factor(ifelse(Rest_g == "나쁨", "Yes", "No"), levels = c("Yes", "No")), # 나쁨 변수
Rest_g_2 = factor(ifelse(Rest_g == "평균이하", "Yes", "No"), levels = c("Yes", "No")), # 평균이하 변수
Rest_g_3 = factor(ifelse(Rest_g == "평균", "Yes", "No"), levels = c("Yes", "No")), # 평균 변수
Rest_g_4 = factor(ifelse(Rest_g == "좋음", "Yes", "No"), levels = c("Yes", "No")), # 좋음 변수
Rest_g_5 = factor(ifelse(Rest_g == "운동선수", "Yes", "No"), levels = c("Yes", "No")) # 운동선수 변수
) %>% select(-Rest) # 각 항목별 변수를 생성
train
#++++++++++++++++++++++++++++++++ 1. "나쁨" +++++++++++++++++++++++++++++++++++++++++++++++++++++
# -- 변수선택 --
model_0 <- glm(Rest_g_1 ~ 1, data = train, family = binomial("logit")) # 널모형
var_step <- step(model_0, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 모든 교호작용이 있고 1차항 까지 있는 모형이 최대모형이다.
# 단게별 선택을 한다.
# -- 모형 생성 --
model_1 <- glm(formula(var_step), data = train, family = binomial("logit")) # 최선의 모형
# -- 모형의 유의성 --
anova(model_0, model_1, test = "Chisq") # 모든 회귀계수가 0인지 검정
# 유의수준 0.05보다 작으므로 귀무가설(모든 회귀계수가 0)을 기각한다.
anova(model_1, test="Chisq") # 각 변수의 유의성 검정.
summary(model_1) # 회귀계수가 0인지에 대한 검정
# Exercie3시간의 회귀계수가 0에 가깝다. 하지만 Exercie2시간은 유의하므로 그대로 둔다.
# -- 결정계수 --
# 결정계수로 모형을 선택한다.
# 결정계수 = (1- exp((D - null_D)/null_df))/(1 - exp(-null_D/null_df))
D <- model_1$deviance
null_D <- model_1$null.deviance
df <- model_1$df.null
R_1 <- (1 - exp((D-null_D)/df))/(1 -exp(-null_D/df)) ; R_1
#++++++++++++++++++++++++++++++++ 2. "평균이하" +++++++++++++++++++++++++++++++++++++++++++++++++++++
# -- 변수선택 --
model_0 <- glm(Rest_g_2 ~ 1, data = train, family = binomial("logit")) # 널모형
var_step <- step(model_0, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 모든 교호작용이 있고 1차항 까지 있는 모형이 최대모형이다.
# 단게별 선택을 한다.
# -- 모형 생성 --
model_2 <- glm(formula(var_step), data = train, family = binomial("logit")) # 최선의 모형
# -- 모형의 유의성 --
anova(model_0, model_2, test = "Chisq") # 모든 회귀계수가 0인지 검정
# 유의수준 0.05보다 작으므로 귀무가설(모든 회귀계수가 0)을 기각한다.
anova(model_2, test="Chisq") # 각 변수의 유의성 검정
# Active:Exercise 변수가 유의하지 않다.따라서 제거한다.
model_2 <- update(model_2, ~ . -Exercise) # 변수제거
anova(model_2, test = "Chisq") # 각 변수의 유의성 검정
# 모든 변수가 유의하다
summary(model_2) # 회귀계수가 0인지에 대한 검정
# Exercie2시간의 회귀계수가 0에 가깝다. 하지만 Exercie3시간은 유의하므로 그대로 둔다.
# -- 결정계수 --
# 결정계수로 모형을 선택한다.
# 결정계수 = (1- exp((D - null_D)/null_df))/(1 - exp(-null_D/null_df))
D <- model_2$deviance
null_D <- model_2$null.deviance
df <- model_2$df.null
R_2 <- (1 - exp((D-null_D)/df))/(1 -exp(-null_D/df)) ; R_2
#++++++++++++++++++++++++++++++++ 3. "평균" +++++++++++++++++++++++++++++++++++++++++++++++++++++
# -- 변수선택 --
model_0 <- glm(Rest_g_3 ~ 1, data = train, family = binomial("logit")) # 널모형
var_step <- step(model_0, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 모든 교호작용이 있고 1차항 까지 있는 모형이 최대모형이다.
# 단게별 선택을 한다.
# -- 모형 생성 --
model_3 <- glm(formula(var_step), data = train, family = binomial("logit")) # 최선의 모형
# -- 모형의 유의성 --
anova(model_0, model_3, test = "Chisq") # 모든 회귀계수가 0인지 검정
# 유의수준 0.05보다 작으므로 귀무가설(모든 회귀계수가 0)을 기각한다.
anova(model_3, test="Chisq") # 각 변수의 유의성 검정
# Sex변수가 유의하지 않다.
model_3 <- update(model_3, ~ . -Sex) # Sex변수 제거
anova(model_3, test ="Chisq") # 각 변수의 유의성 검정
summary(model_3) # 회귀계수가 0인지에 대한 검정
# Exercie2시간의 회귀계수가 0에 가깝다. 하지만 Exercie3시간은 유의하므로 그대로 둔다.
# -- 결정계수 --
# 결정계수로 모형을 선택한다.
# 결정계수 = (1- exp((D - null_D)/null_df))/(1 - exp(-null_D/null_df))
D <- model_3$deviance
null_D <- model_3$null.deviance
df <- model_3$df.null
R_3 <- (1 - exp((D-null_D)/df))/(1 -exp(-null_D/df)) ; R_3
#++++++++++++++++++++++++++++++++ 4. "좋음" +++++++++++++++++++++++++++++++++++++++++++++++++++++
# -- 변수선택 --
model_0 <- glm(Rest_g_4 ~ 1, data = train, family = binomial("logit")) # 널모형
var_step <- step(model_0, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "both") # 모든 교호작용이 있고 1차항 까지 있는 모형이 최대모형이다.
# 단게별 선택을 한다.
# -- 모형 생성 --
model_4 <- glm(formula(var_step), data = train, family = binomial("logit")) # 최선의 모형
# -- 모형의 유의성 --
anova(model_0, model_4, test = "Chisq") # 모든 회귀계수가 0인지 검정
# 유의수준 0.05보다 작으므로 귀무가설(모든 회귀계수가 0)을 기각한다.
summary(model_4) # 회귀계수가 0인지에 대한 검정
# Exercie2시간의 회귀계수가 0에 가깝다. 하지만 Exercie3시간은 유의하므로 그대로 둔다.
# -- 결정계수 --
# 결정계수로 모형을 선택한다.
# 결정계수 = (1- exp((D - null_D)/null_df))/(1 - exp(-null_D/null_df))
D <- model_4$deviance
null_D <- model_4$null.deviance
df <- model_4$df.null
R_4 <- (1 - exp((D-null_D)/df))/(1 -exp(-null_D/df)) ; R_4
#++++++++++++++++++++++++++++++++ 5. "운동선수" +++++++++++++++++++++++++++++++++++++++++++++++++++++
# -- 변수선택 --
model_0 <- glm(Rest_g_5 ~ 1, data = train, family = binomial("logit")) # 널모형
var_step <- step(model_0, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 모든 교호작용이 있고 1차항 까지 있는 모형이 최대모형이다.
# 단게별 선택을 한다.
# -- 모형 생성 --
model_5 <- glm(formula(var_step), data = train, family = binomial("logit")) # 최선의 모형
# -- 모형의 유의성 --
anova(model_0, model_5, test = "Chisq") # 모든 회귀계수가 0인지 검정
# 유의수준 0.05보다 작으므로 귀무가설(모든 회귀계수가 0)을 기각한다.
anova(model_5, test="Chisq") # 각 변수의 유의성 검정
# 변수가 유으하다.
summary(model_5) # 회귀계수가 0인지에 대한 검정
# Exercise의 두 더비변수의 회귀계수가 모두 0이므로 Exercise를 제거한다.
model_5 <- update(model_5, ~ . -Exercise) #Exercise제거
summary(model_5) # 회귀계수가 0인지에 대한 검정
# 모든 변수의 회귀계수가 0이 아니다.
# -- 결정계수 --
# 결정계수로 모형을 선택한다.
# 결정계수 = (1- exp((D - null_D)/null_df))/(1 - exp(-null_D/null_df))
D <- model_5$deviance
null_D <- model_5$null.deviance
df <- model_5$df.null
R_5 <- (1 - exp((D-null_D)/df))/(1 -exp(-null_D/df)) ; R_5
#++++++++++++++++++++++++++++++++ 5. 다항분포 +++++++++++++++++++++++++++++++++++++++++++++++++++++
# -- 변수 선택 --
model_0 <- multinom(Rest_g ~ 1, data = train) # 널모형
var_step <- step(model_0, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "both")
model_multi <- multinom(formula(var_step), data = train) # 최선의 모형
# -- 모형의 유의성 --
anova(model_0, model_multi, test = "Chisq") # 모든 회귀계수에 대한 검정
# p-value가 유의수준 0.05보다 작으므로 모든 회귀계수는 0이 아니라 할 수 있다.
summary(model_multi) # 각 회귀계수에 대한 검정
# 검정 통계량이 계산되지 않으므로 회귀계수/오차로 검정톨계량을 계산하여 z-test를 실시한다.
a <- summary(model_multi)$coefficients/ summary(model_multi)$standard.errors # 검정통계량
a <-pnorm(abs(a), mean = 0, sd = 1, lower.tail = FALSE)*2; a # p-value
a > 0.05 # 0.05보다 큰지 확인
# '좋음'에서 Exercise3시간의 회귀계수가 유의하고 3시간은 유의하지 않다.
# Sex는 "운동선수"의 모형에서만 유의하다.
#++++++++++++++++++++++++++++++++ 6. 모형평가 +++++++++++++++++++++++++++++++++++++++++++++++++++++
# -- 로지스틱모형 예측 --
R_1; R_2; R_3; R_4; R_5 # 각 모형의 결정계수
# 1번, 2번, 3번, 5번 모형을 사용한다.
pred_1 <- predict(model_1, newdata = test.Pulse, type = "response") # 모형1의 예측값
pred_2 <- predict(model_2, newdata = test.Pulse, type = "response") # 모형1의 예측값
pred_3 <- predict(model_3, newdata = test.Pulse, type = "response") # 모형3의 예측값
pred_5 <- predict(model_5, newdata = test.Pulse, type = "response") # 모형4의 예측값
levels(train$Rest_g_1) # 예측값의 숫자는 "No"일 확률의 예상값이다.
# 0.5보다 낮으면 Yes로 예측한다.
pred_logit <- sapply(1:nrow(test.Pulse),function(x){ # 예측값의 확률로 범주형으로 예측한다.
a <- min(pred_1[x],pred_2[x],pred_3[x],pred_5[x]) # 가장 낮은 확률
if(a >= 0.5){ # 모든 예측값이 No이면 "보통"으로 반환
return("좋음")
}else if(pred_1[x] == a){ # 가장 낮은 예측값이 모형1이면 "나쁨"
return("나쁨")
}else if(pred_2[x] == a){ # 가장 낮은 예측값이 모형2이면 "평균이하"
return("평균이하.")
}else if(pred_3[x] == a){ # 가장 낮은 예측값이 모형3이면 "평균"
return("평균.")
}else { # 가장 낮은 예측값이 모형5이면 "운동선수"
return("운동선수")
}
})
pred_logit <- factor(pred_logit, levels = c("운동선수", "좋음", "평균", "평균이하", "나쁨")) ; pred_logit
confusionMatrix(pred_logit, test.Pulse$Rest_g) # 59%의 예측률
# -- 다항분포 모형의 예측 --
pred_multi <- predict(model_multi, newdata = test.Pulse, type = "class") # 다항분포의 로지스틱 모형
confusionMatrix(pred_multi, test.Pulse$Rest_g) # 57%의 예측률
error_rate <- c("binomial_logistic" = 0.587)
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", 'error', "error_rate")))
# ================================ 11. 로지스틱 회귀분석(bootstrap) ================================
# 로지스틱 회귀분석은 multinom()을 사용한다.
# 100개의 표본으로 모형을 만든다. 각 표본의 크기는 150으로 한다.
# 각 모형은 step()함수를 이용하여 단계별 선택법을 사용한다.
# 각 회귀계수에 대한 검정은 하지않는다.
train <- train.Pulse %>% select(-Rest)
train
# -- 모형생성 --
bootstrap <- function(){ # 표본을 추출하여 모형을 생성하여 예측값을 반환
index = sample(1:nrow(train), replace = TRUE, size = 150) # 데이터의 행을 추출
temp <<- train[index,] # 데이터 추출
model <- multinom(Rest_g ~ 1, data = temp) # 널모형
var_step <- step(model, scope = list(lower = ~ 1, upper = ~ Active*Smoke*Sex*Exercise*fat_g), direction = "forward") # 단계별 선택법으로 변수를 선택
model <- multinom(formula(var_step), data= temp) # 최선의 모형생성
return(model)
}
set.seed(20130691)
pred <- sapply(1:100, function(x){ # 100개 모형의 예측값을 계산한다.
model <- bootstrap()
predict(model, newdata = test.Pulse, type = "class")
})
# -- 예측값 --
pred <- apply(pred, 1, function(x){ # 각 행마다 계산
# 각 행의 예측값중 가장 많이 나온 범주를 사용한다.
# 같은 개수이면 levels이 낮은 값을 사용한다.
a <- table(x) # 개수를 개산
return(names(which.max(a))) # 빈도가 가장높은 첫번째 값의 이름을 출력
}) # 예측값들을 통합하여 가장 많이 나오는 예측값을 추출
pred <- factor(pred, levels = c("운동선수", "좋음", "평균", "평균이하", "나쁨")) # factor로 변환
# -- 오류율 --
confusionMatrix(pred, test.Pulse$Rest_g) # 63%의 정확도
error_rate <- c(error_rate, "multinomial_logistic_bootstrap" = 0.6304) # 정분류율저장
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error", "error_rate")))
# ================================ 12. 의사결정나무 ================================
# 의사결정나무는 rpart패키지를 사용한다.
# 튜닝은 cp값만 조정한다.
train <- train.Pulse %>% select(-Rest)
train
# -- 튜닝 --
# 의사결정나무를 최대로 성장시킨후 오차가 최소가 되는 cp값을 찾는다.
mycontrol <- rpart.control(minsplit = 12, maxdepth = 30) # 기본 나무모형의 설정
# 최소 분할하는 갯수가 12개, 최대깊이가 30
set.seed(20130691)
model <- rpart(Rest_g ~ ., data = train, cp = -1, control = mycontrol) # 최대 성장 나무모형
# 나무모형 그림
plot(model, main = "가지치기 전")
text(model)
printcp(model) # cp값과 오차 테이블
index <- which.min(model$cptable[,"xerror"]) # 최소오차의 위치
cp <- model$cptable[index, "CP"] # 최소오차인 cp
model$cptable[index,] # 총 6번의 분할을 실시
# -- 가지치기 --
model <- prune(model, cp = cp) # 가지치기
model # 나무모형
par(mar = c(2,0.5,2,0.5))
plot(model, margin = 0.5, main = "가지치기 후") # 나무모형 그림
text(model, col = "red")
model$variable.importance
# -- 오분류율 --
pred <- predict(model, newdata = test.Pulse, type = "class") # 예측값
confusionMatrix(pred, test.Pulse$Rest_g) # 59%의 정확도
error_rate <- c(error_rate, "decision_Tree" = 0.587) # 오차추가
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error", "error_rate")))
# ================================ 13. 의사결정나무(bootstrap) ================================
# 100개의 표본을 추출하여 모형을 생성한다.
# 각 모형마다 예측값을 구하고 가장 많이 나온 값을 예측값으로 사용한다.
train <- train.Pulse %>% select(-Rest)
train
# -- 모형 만들기 --
bootstrap <- function(){
# - 데이터 추출 - #
index <- sample(1:length(train.Pulse), replace= TRUE, size = 150) # 각 표본의 크기는 150이다.
temp <- train.Pulse[index,] # 데이터 추출
# - 튜닝 - #
mycontrol <- rpart.control(minsplit = 10) # 나무모형의 기본설정
model <- rpart(Rest_g ~ ., data = temp, control = mycontrol, cp = -1) # 최대 성장 나무모형
table_index <- which.min(model$cptable[,"xerror"]) # 오차가 최소인 위치
cp <- model$cptable[table_index, "CP"] # 오차가 최소인 cp값
# - 가지치기 - #
model <- prune(model, cp = cp) # cp값으로 가지치기
return(model) # 모형반환
}
# -- 예측값 계산 --
set.seed(20130691) # 랜덤성분 제어
pred <- sapply(1:100, function(x){ # 예측값 생성
model <- bootstrap() # 모형생성
predict(model, newdata = test.Pulse, type = "class")
})
pred <- apply(pred, 1, function(x){ # 각 행마다 계산
a <- table(x) # 빈도계산
return(names(which.max(a))) # 최빈 항목을 반환한다.
})
pred <- factor(pred, levels = c("운동선수", "좋음", "평균", "평균이하", "나쁨"))
# -- 오분류율 --
confusionMatrix(pred, test.Pulse$Rest_g) # 52%의 정확도
error_rate <- c(error_rate, "decision_tree_bootstrap" = 0.5217) # 오차저장
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error", "error_rate")))
# ================================ 14. 신경망 모형(SVM) ================================
# 신경망모형은 e1071패키지의 svm()를 사용한다.
# 커널은 다항식 커널을 사용한다.
train <- train.Pulse %>% select(-Rest)
train
# -- 튜닝 --
set.seed(20130691)
result <- tune.svm(Rest_g ~ ., data = train, degree = 1:10, cost = seq(0.01, 5, 0.01), kernel = "polynomia") # 여러 종류의 모형을 생성하여 파라미터를 찾는다.
result$best.parameters # degree = 6, cost = 4.21
# -- 모형생성 --
set.seed(20130691)
model <- svm(Rest_g ~ ., data = train, degree = 6, cost = 4.21, kernel = "polynomia") # 서포트 벡터 분류기 생성
summary(model)
# -- 예측 --
pred <- predict(model, newdata = test.Pulse, type = "class")
pred
# -- 오차 --
confusionMatrix(pred, test.Pulse$Rest_g) # 57%의 정확도
error_rate <- c(error_rate, "SVM" = 0.5652) # 오차저장
rm(list = setdiff(ls(), c("Pulse","train.Pulse", "test.Pulse", "error", "error_rate")))
# ================================ 15. 앙상블 예측(bagging) ================================
# package "ipred"를 사용한다.
# 반복횟수는 100회
# bagging(fomula, data, mfinal = 100, control)
# nbagg : Tree반복생성 횟수, control : 통제(rpart.control을 사용)
# rpart.control(minsplit = 5, minbucket = round(minsplit/3), cp = 0.01, maxcompete = 4, maxsurrogate = 5, usesurrogate = 2, xval = 10, surrogatestyle = 0, maxdepth = 30, ..)
# minsplit = Tree의 최소 Node수, minbucket = Terminal Node의 최소 자료수
# cp = 복잡성 계수, max.compete = 분기후보, maxsurrogate = 최대 분기후보, usesurrogate = 허용 분기후보
# xval = 교차인정값, maxdepth = 최대깊이
train <- train.Pulse %>% select(-Rest) # 모형을 생성하는 데이터
train
# -- 모형생성 --
# bagging은 동일한 가중치로 복원추출하여 여러 데이터를 만들어 각 데이터마다 분류기를 만드는 방법이다.
# bagging방법은 각각의 분류기의 분산을 감소시킴으로 오차를 줄이는 방법이다.
# 따라서 adabag에서 사용되는 bagging의 나무모형은 가지치기를 하지않는다.
set.seed(20130691)
model <- bagging(Rest_g ~ ., data = train, nbagg = 100, coob = TRUE) # 모형생성
model # 모형
# Out-of-Bag들의 오차가 0.5323이 된다.
# -- 예측값 --
pred <- predict(model, newdata = test.Pulse, type = "class") # 예측값
pred <- factor(pred, levels = c("운동선수", "좋음", "평균", "평균이하", "나쁨"))
table(pred) # 예측값들의 빈도
# -- 오차 --
confusionMatrix(pred, test.Pulse$Rest_g) # 오분류율
# 52%의 정확도를 가진다.
error_rate <- c(error_rate, "classificaton_bagging" = 0.5217) # 오차저장
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error", "error_rate")))
# ================================ 16. 앙상블 예측(boosting) ================================
# boosting은 gradien boosting을 사용한다.(gbm을 사용)
# n.trees : 생성할 나무의 개수
# interaction.depth : 각 노드에서 뻗어나갈 가지의 개수를 지정(나무의 깊이지정), 이 값이 N이면 2*N+1개의 terminal node가 생성
# shrinkage : Learning Rate로, gradient boosting을 학습하는데 걸리는 시간을 조정, 만약 이값이 크다면 빨리 학습을 수행하지만, 수행도중 중요한점을 놓칠 수 있다. 때문에 이 값이 크면 n.trees인자로 나무를 많이 생성함으로 보완한다. shrinkage가 작고 n.trees가 크면 과적합이 발생 할 수 있다.
# n.minobsinnode : 트리의 terminal node의 최소 관측수, 정수로 입력
# 종속변수가 범주형이면 distributio = "bernoulli"를 입력
train <- train.Pulse %>% select(-Rest)
train
# -- 튜닝 --
control <- caret::trainControl(method = "boot", search = "grid") # 함수 수행을 위한 셋팅
# method = "boot" : 샘플링방버을 설정 - 반복학습/ 테스트 분할할 경우 "boot" 사용
# search = "grid" 옵션은 그리드 검색 루틴을 사용한다. - "random"일 경우 무작위 검색절차
tunegrid <- expand.grid(n.trees = seq(20,100,10), interaction.depth = 1:10, shrinkage = 0.1, n.minobsinnode = 1:10) # gbm인자들을 비교할 목록
set.seed(20130691)
gbm_tuning <- caret::train(Rest_g ~ ., data = train, method = "gbm",tuneGrid = tunegrid, trControl = control, metric = "Accuracy") # train함수로 모형생성
# tuneGrid : 가능한 튜닝값을 가진 데이터 프레임을 받아들인다.
# metric : 최적모형을 구분하는 기준, 회귀분석의 경우 "RMSE"나 "Rsquared"를 상요하고 분류의 경우 "Accuracy"나 "Kappa"를 사용한다.
print(gbm_tuning)
# n.trees = 20, interaction.depth = 1, shrinkage = 0.1, n.minobsinnode = 9
# -- 모형 생성 --
model <- gbm(Rest_g ~ . ,data = train, n.trees = 20, interaction.depth = 1, shrinkage = 0.1, n.minobsinnode = 9) # boosting 모형생성
summary(model) # Active변수를 중요하게 사용했다.
# -- 오차 --
set.seed(20130691)
pred <- predict(model, newdata = test.Pulse, type = "response", n.trees = 20)
pred # 각 항목이 나타난 비율을 나타낸다.
pred <- apply(pred, 1, function(x){
index <- which.max(x) # 최대 비율을 갖는 위치를 찾는다.
name <- switch(index, "운동선수", "좋음", "평균", "평균이하", "나쁨") # 최대비율을 갖는 항목을 출력한다.
return(name)
}) # 최대비율을 갖는 항목을 예측값으로 변환
pred <- factor(pred, levels = c("운동선수", "좋음", "평균", "평균이하", "나쁨")) # factor로 변환
confusionMatrix(pred, test.Pulse$Rest_g) # 59%의 정확도를 갖는다.
error_rate <- c(error_rate, "classification_boosting" = 0.587) # 오차저장
rm(list = setdiff(ls(), c("Pulse", "train.Pulse", "test.Pulse", "error", "error_rate")))
# ================================ 17. 결과 ================================
error_rate
sort(error_rate, decreasing = TRUE) # 의사결정 나무가
# SVM > logistic > Tree > boosting > bagging 순으로 정확하다.
'데이터 마이닝' 카테고리의 다른 글
| 모형의 평가방법( 예측모형) (0) | 2019.08.29 |
|---|---|
| Kaggle 대회 데이터 분석 (0) | 2019.08.28 |
| 머신러닝 시작 (0) | 2019.08.09 |
| WineQuality 분석 (0) | 2019.07.26 |
| 통계상담과제 (0) | 2019.07.22 |