2022-06-10-Data-mining-theory-and-reality-midtestExample03

3. 대학의 이교수는 영어를 가르치고 있다. 이교수는 영어 수업이 효과가 있는지를 검증하려고 한다. 영어 수업을 듣는 학생들을 대상으로 첫주에 영어시험을 보고, 한학기가 끝난 후에 다시 영어시험을 보았다. 과연 이교수의 영어 수업은 학생들의 영어점수를 향상시켜 주었는가?

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()
## * Use suppressPackageStartupMessages() to eliminate package startup messages
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)
brain_tb <- read_csv('data3.csv', 
                     col_names = TRUE,
                     locale=locale('ko', encoding='euc-kr'), # 한글
                     na=".") %>%
  mutate_if(is.character, as.factor)
## Rows: 40 Columns: 2
## -- Column specification --------------------------------------------------------
## Delimiter: ","
## dbl (2): 사전, 사후
## 
## 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(brain_tb)
## tibble [40 x 2] (S3: tbl_df/tbl/data.frame)
##  $ 사전: num [1:40] 71 56 76 76 82 83 64 70 72 73 ...
##  $ 사후: num [1:40] 68 54 74 74 80 81 63 69 71 72 ...
brain_tb
## # A tibble: 40 x 2
##     사전  사후
##    <dbl> <dbl>
##  1    71    68
##  2    56    54
##  3    76    74
##  4    76    74
##  5    82    80
##  6    83    81
##  7    64    63
##  8    70    69
##  9    72    71
## 10    73    72
## # ... with 30 more rows

long형으로 변형

brain_tb_long <- brain_tb %>%
  select(사전, 사후) %>% # 필요한 변수만 선택
  pivot_longer(c("사전","사후"), 
               names_to = "구분",
               values_to = "영어점수")
brain_tb_long
## # A tibble: 80 x 2
##    구분  영어점수
##    <chr>    <dbl>
##  1 사전        71
##  2 사후        68
##  3 사전        56
##  4 사후        54
##  5 사전        76
##  6 사후        74
##  7 사전        76
##  8 사후        74
##  9 사전        82
## 10 사후        80
## # ... with 70 more rows

차이 추가

brain_tb <- brain_tb %>%
  mutate(차이 = 사후-사전) # mutate -> 새로운 칼럼 생성

brain_tb
## # A tibble: 40 x 3
##     사전  사후  차이
##    <dbl> <dbl> <dbl>
##  1    71    68    -3
##  2    56    54    -2
##  3    76    74    -2
##  4    76    74    -2
##  5    82    80    -2
##  6    83    81    -2
##  7    64    63    -1
##  8    70    69    -1
##  9    72    71    -1
## 10    73    72    -1
## # ... with 30 more rows

기술통계분석

skim(brain_tb)
Data summary
Name brain_tb
Number of rows 40
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
사전 0 1 74.15 6.43 54 72 75.0 77 85 ▁▁▂▇▂
사후 0 1 74.00 6.54 54 72 75.5 77 85 ▁▁▃▇▂
차이 0 1 -0.15 1.31 -3 -1 0.0 1 3 ▃▇▇▅▂
brain_tb %>% 
  get_summary_stats(사전, 사후, 차이)
## # A tibble: 3 x 13
##   variable     n   min   max median    q1    q3   iqr   mad  mean    sd    se
##   <chr>    <dbl> <dbl> <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 사전        40    54    85   75      72    77     5  3.71 74.2   6.43 1.02 
## 2 사후        40    54    85   75.5    72    77     5  3.71 74     6.54 1.03 
## 3 차이        40    -3     3    0      -1     1     2  1.48 -0.15  1.31 0.207
## # ... with 1 more variable: ci <dbl>

그래프 그리기(박스그래프,히스토그램)

brain_tb %>% 
  ggplot(aes(y = 차이)) +
  geom_boxplot() +
  labs(y = "차이")

brain_tb %>% 
  ggplot(aes(x = 차이)) +
  geom_histogram(binwidth = 0.05, 
                 color = "white", 
                 fill = "steelblue")

이상치 제거

이상치 확인

brain_tb %>%
  identify_outliers(차이)
## [1] 사전       사후       차이       is.outlier is.extreme
## <0 행> <또는 row.names의 길이가 0입니다>

정규분포 검정

brain_tb %>%
  shapiro_test(차이)
## # A tibble: 1 x 3
##   variable statistic      p
##   <chr>        <dbl>  <dbl>
## 1 차이         0.952 0.0912

0.05 보다 크므로 정규분포 만족

paired t-검정

brain_tb_long %>% 
  t_test(formula = 영어점수 ~ 구분,
         ref.group = "사후",
         paired = TRUE,
         alternative = "two.sided",
         detailed = TRUE)
## # A tibble: 1 x 13
##   estimate .y.      group1 group2    n1    n2 statistic     p    df conf.low
## *    <dbl> <chr>    <chr>  <chr>  <int> <int>     <dbl> <dbl> <dbl>    <dbl>
## 1    -0.15 영어점수 사후   사전      40    40    -0.723 0.474    39   -0.569
## # ... with 3 more variables: conf.high <dbl>, method <chr>, alternative <chr>
brain_tb_long %>% 
  cohens_d(formula = 영어점수 ~ 구분, 
           paired = TRUE)
## # A tibble: 1 x 7
##   .y.      group1 group2 effsize    n1    n2 magnitude 
## * <chr>    <chr>  <chr>    <dbl> <int> <int> <ord>     
## 1 영어점수 사전   사후     0.114    40    40 negligible

추론(infer)을 이용한 가설검정 및 그래프

표본평균(x)을 이용한 검정그래프

표본평균(x) 계산

x_bar <- brain_tb %>%
  specify(response = 차이) %>%    # hypothesize 없음
  calculate(stat = "mean") %>%      # stat = "mean"
  print()
## Response: 차이 (numeric)
## # A tibble: 1 x 1
##    stat
##   <dbl>
## 1 -0.15

Bootstrapping을 이용한 귀무가설 분포 생성

set.seed(123) 
null_dist_x <- brain_tb %>%
  specify(response = 차이) %>%
  hypothesize(null = "point", 
              mu = 0) %>% 
  generate(reps = 1000, 
           type = "bootstrap") %>%
  calculate(stat = "mean") %>%
  print()
## Response: 차이 (numeric)
## Null Hypothesis: point
## # A tibble: 1,000 x 2
##    replicate      stat
##        <int>     <dbl>
##  1         1 -2.50e-17
##  2         2  2.5 e- 1
##  3         3  7.5 e- 2
##  4         4  1.75e- 1
##  5         5 -1.25e- 1
##  6         6 -2.50e- 2
##  7         7  3.5 e- 1
##  8         8  2.25e- 1
##  9         9 -2   e- 1
## 10        10  2.50e- 2
## # ... with 990 more rows

신뢰구간 생성

null_dist_ci <- null_dist_x %>%
  get_ci(level = 0.95, 
         type = "percentile") %>%
  print()
## # A tibble: 1 x 2
##   lower_ci upper_ci
##      <dbl>    <dbl>
## 1   -0.376    0.425

그래프 그리기

null_dist_x %>%
  visualize() +                                       # method 없음
  shade_p_value(obs_stat = x_bar,
                direction = "two-sided") +            # x_bar
  shade_confidence_interval(endpoints = null_dist_ci) # CI

p_value

null_dist_x %>%
  get_p_value(obs_stat = x_bar,           
              direction = "two-sided") %>%
  print()
## # A tibble: 1 x 1
##   p_value
##     <dbl>
## 1   0.498

P값이 0.05 보다 크므로 영어 수업은 학생들의 영어점수를 향상시켜 주었다고 볼수 없다.

t값을 이용한 검정그래프

t_cal <- brain_tb %>%
  specify(response = 차이) %>%
  hypothesize(null = "point",         # hypothesize 필요
              mu = 0) %>%  
  calculate(stat = "t") %>%           # stat = "t"              
  print()
## Response: 차이 (numeric)
## Null Hypothesis: point
## # A tibble: 1 x 1
##     stat
##    <dbl>
## 1 -0.723

Bootstrapping을 이용한 귀무가설 분포 생성

set.seed(123) 
null_dist_t <- brain_tb %>%
  specify(response = 차이) %>%
  hypothesize(null = "point", 
              mu = 0) %>% 
  generate(reps = 1000, 
           type = "bootstrap") %>%
  calculate(stat = "t") %>%
  print()
## Response: 차이 (numeric)
## Null Hypothesis: point
## # A tibble: 1,000 x 2
##    replicate      stat
##        <int>     <dbl>
##  1         1 -1.47e-16
##  2         2  1.40e+ 0
##  3         3  4.24e- 1
##  4         4  7.33e- 1
##  5         5 -9.02e- 1
##  6         6 -1.22e- 1
##  7         7  1.91e+ 0
##  8         8  1.10e+ 0
##  9         9 -1.06e+ 0
## 10        10  1.27e- 1
## # ... with 990 more rows

신뢰구간 생성

null_dist_ci <- null_dist_t %>%
  get_ci(level = 0.95, 
         type = "percentile") %>%
  print()
## # A tibble: 1 x 2
##   lower_ci upper_ci
##      <dbl>    <dbl>
## 1    -1.97     2.06
null_dist_t %>%
  visualize(method = "both") +           #method = "both": 이론분포+boot분포
  shade_p_value(obs_stat = t_cal,
                direction = "two-sided") +
  shade_confidence_interval(endpoints = null_dist_ci)
## Warning: Check to make sure the conditions have been met for the theoretical
## method. {infer} currently does not check these for you.

Leave a comment