데이터마이닝의 이론과 실제 기말고사 연습문제4

연습문제4

Toyota 중고차 가격을 결정하는 모델을 만들고자 한다.

주행거리에 따라서 중고차 가격이 어떻게 결정되는지 가격모델을 만들어 보세요

price: 가격 (유로)

km: 주행킬로미터 (kilometers driven)

단순 선형회귀분석(Simple linear regression)

1.기본 package 설정

1.1 library 로드

library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.5     v purrr   0.3.4
## v tibble  3.1.6     v dplyr   1.0.8
## v tidyr   1.2.0     v stringr 1.4.0
## v readr   2.1.2     v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(tidymodels)
## -- Attaching packages -------------------------------------- tidymodels 0.2.0 --
## v broom        0.8.0     v rsample      0.1.1
## v dials        0.1.1     v tune         0.2.0
## v infer        1.0.0     v workflows    0.2.6
## v modeldata    0.1.1     v workflowsets 0.2.1
## v parsnip      0.2.1     v yardstick    0.0.9
## v recipes      0.2.0
## -- Conflicts ----------------------------------------- tidymodels_conflicts() --
## x scales::discard() masks purrr::discard()
## x dplyr::filter()   masks stats::filter()
## x recipes::fixed()  masks stringr::fixed()
## x dplyr::lag()      masks stats::lag()
## x yardstick::spec() masks readr::spec()
## x recipes::step()   masks stats::step()
## * Dig deeper into tidy modeling with R at https://www.tmwr.org
library(rstatix)
## 
## 다음의 패키지를 부착합니다: 'rstatix'
## The following objects are masked from 'package:infer':
## 
##     chisq_test, prop_test, t_test
## The following object is masked from 'package:dials':
## 
##     get_n
## The following object is masked from 'package:stats':
## 
##     filter
library(skimr)

2.데이터 불러오기

toyota_tb <- read_csv('data4.csv', 
                      col_names = TRUE,
                      locale=locale('ko', encoding='euc-kr'), # 한글
                      na=".") %>%
  mutate_if(is.character, as.factor)
## Rows: 112 Columns: 3
## -- Column specification --------------------------------------------------------
## Delimiter: ","
## dbl (3): id, price, km
## 
## i Use `spec()` to retrieve the full column specification for this data.
## i Specify the column types or set `show_col_types = FALSE` to quiet this message.
str(toyota_tb)
## tibble [112 x 3] (S3: tbl_df/tbl/data.frame)
##  $ id   : num [1:112] 1 2 3 4 5 6 7 8 9 10 ...
##  $ price: num [1:112] 12995 13900 13950 14250 13450 ...
##  $ km   : num [1:112] 29198 22000 20019 20000 17003 ...
skim(toyota_tb)
Data summary
Name toyota_tb
Number of rows 112
Number of columns 3
_______________________
Column type frequency:
numeric 3
________________________
Group variables None

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
id 0 1 56.50 32.48 1 28.75 56.5 84.25 112 ▇▇▇▇▇
price 0 1 9908.37 2139.27 6250 8400.00 9722.5 11485.00 14900 ▆▇▆▅▂
km 0 1 47735.07 17098.59 12500 34448.75 49955.5 60834.25 74963 ▃▅▆▇▇

3.기본통계치 확인

toyota_tb %>%
  get_summary_stats(price, km)
## # A tibble: 2 x 13
##   variable     n   min   max median     q1     q3    iqr    mad   mean     sd
##   <chr>    <dbl> <dbl> <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>
## 1 km         112 12500 74963 49956. 34449. 60834. 26386. 18729. 47735. 17099.
## 2 price      112  6250 14900  9722.  8400  11485   3085   2294.  9908.  2139.
## # ... with 2 more variables: se <dbl>, ci <dbl>

4.그래프 그리기(산점도)

toyota_tb %>%
  ggplot(mapping = aes(x = km,
                       y = price)) +
  geom_jitter() +
  geom_smooth(method = "lm", se = TRUE)
## `geom_smooth()` using formula 'y ~ x'

5.단순회귀분석

표준화 회귀계수

install.packages(“lm.beta”)

library(lm.beta)

모델 생성

toyota_fit <- lm(price ~ km, data = toyota_tb) %>%
  lm.beta()      # 표준화 회귀계수

ANOVA 분석

anova(toyota_fit)
## Analysis of Variance Table
## 
## Response: price
##            Df    Sum Sq   Mean Sq F value    Pr(>F)    
## km          1 428569798 428569798  593.58 < 2.2e-16 ***
## Residuals 110  79420788    722007                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

회귀계수

tidy(toyota_fit, conf.int = TRUE) %>%
  mutate_if(is.numeric, round, 5)
## # A tibble: 2 x 8
##   term      estimate std_estimate std.error statistic p.value conf.low conf.high
##   <chr>        <dbl>        <dbl>     <dbl>     <dbl>   <dbl>    <dbl>     <dbl>
## 1 (Interce~  1.54e+4       NA     239.           64.4       0   NA        NA    
## 2 km        -1.15e-1       -0.919   0.00472     -24.4       0   -0.928    -0.909

설명력R2 : r.squared 값이 1에 가까울수록 실제값과 예측값의 차이가 거의 없는것이다.

glance(toyota_fit)
## Warning: The `glance()` method for objects of class lm.beta is not maintained
## by the broom team, and is only supported through the lm tidier method. Please be
## cautious in interpreting and reporting broom output.
## # A tibble: 1 x 12
##   r.squared adj.r.squared sigma statistic  p.value    df logLik   AIC   BIC
##       <dbl>         <dbl> <dbl>     <dbl>    <dbl> <dbl>  <dbl> <dbl> <dbl>
## 1     0.844         0.842  850.      594. 3.90e-46     1  -913. 1833. 1841.
## # ... with 3 more variables: deviance <dbl>, df.residual <int>, nobs <int>

6. 회귀분석 가정 검정

등분산성: Scale-Location, ncvTest

정규성: Nomal Q-Q, shapiro.test

선형성: Residuals vs Fitted,

독립성: durbinWatsonTest

이상치검정 : Residuals vs Leverage(cook’s distance) 4/n-k-1

그림으로 가정 검정

opar <- par(no.readonly = TRUE)
  par(mfrow=c(2,2))
  par("mar")
## [1] 5.1 4.1 4.1 2.1
  par(mar=c(3,3,3,3))
  plot(toyota_fit)

  par(opar)

잔차의 정규분포 검정 :P값이 0.05보다 크므로 정규분포를 따른다

shapiro_test(toyota_fit$residuals)
## # A tibble: 1 x 3
##   variable             statistic p.value
##   <chr>                    <dbl>   <dbl>
## 1 toyota_fit$residuals     0.987   0.374

수치로 가정 검정

library(car)
## 필요한 패키지를 로딩중입니다: carData
## 
## 다음의 패키지를 부착합니다: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
## The following object is masked from 'package:purrr':
## 
##     some

잔차의 등분산성 검정

유의값(p-value)이 .05 이상인 경우에는 귀무가설을 지지하여 ’분산에 차이가 없다.’고 말할 수 있다.

car::ncvTest(toyota_fit)
## Non-constant Variance Score Test 
## Variance formula: ~ fitted.values 
## Chisquare = 1.729766, Df = 1, p = 0.18844

이상치 확인

car::outlierTest(toyota_fit)
## No Studentized residuals with Bonferroni p < 0.05
## Largest |rstudent|:
##     rstudent unadjusted p-value Bonferroni p
## 45 -2.385319            0.01879           NA

#이상치 검정, sd, hat, d 통합검정

car::influencePlot(toyota_fit, id.method="identify")
## Warning in plot.window(...): "id.method"는 그래픽 매개변수가 아닙니다
## Warning in plot.xy(xy, type, ...): "id.method"는 그래픽 매개변수가 아닙니다
## Warning in axis(side = side, at = at, labels = labels, ...): "id.method"는 그래
## 픽 매개변수가 아닙니다

## Warning in axis(side = side, at = at, labels = labels, ...): "id.method"는 그래
## 픽 매개변수가 아닙니다
## Warning in box(...): "id.method"는 그래픽 매개변수가 아닙니다
## Warning in title(...): "id.method"는 그래픽 매개변수가 아닙니다
## Warning in plot.xy(xy.coords(x, y), type = type, ...): "id.method"는 그래픽 매개
## 변수가 아닙니다

##         StudRes        Hat        CookD
## 13  -0.01681014 0.04577181 6.839490e-06
## 14   1.13783401 0.04718519 3.197152e-02
## 45  -2.38531877 0.01133391 3.127959e-02
## 55  -1.86253610 0.01930106 3.338755e-02
## 86  -2.18997344 0.01696279 3.999824e-02
## 112 -2.21119165 0.01072818 2.560606e-02

이상치 제거

toyota_tb <- toyota_tb[-c(13:14),]
#toyota_tb <- toyota_tb[-c(47,66,42,4,8),]

5.회귀분석 다시 실시

07.모델을 이용한 예측

콜레스테롤이 130, 150일 경우 예측값

toyota_tb_new <- data.frame(km=c(75000))
predict(toyota_fit, newdata = toyota_tb_new)
##        1 
## 6775.129

Leave a comment