Ch 11. 지도시각화

11 - 1. 미국 주별 강력 범죄율 단계 구분도 만들기

구분 내용
Data USArrests(R base data, 1973년 미국 주별 강력 범죄율 데이터), state
Packages ggiraphExtra,tibble, dplyr, ggplot2

패키지 장착

library(ggiraphExtra)
library(tibble)

데이터 확인

head(USArrests)
##            Murder Assault UrbanPop Rape
## Alabama      13.2     236       58 21.2
## Alaska       10.0     263       48 44.5
## Arizona       8.1     294       80 31.0
## Arkansas      8.8     190       50 19.5
## California    9.0     276       91 40.6
## Colorado      7.9     204       78 38.7
str(USArrests)
## 'data.frame':    50 obs. of  4 variables:
##  $ Murder  : num  13.2 10 8.1 8.8 9 7.9 3.3 5.9 15.4 17.4 ...
##  $ Assault : int  236 263 294 190 276 204 110 238 335 211 ...
##  $ UrbanPop: int  58 48 80 50 91 78 77 72 80 60 ...
##  $ Rape    : num  21.2 44.5 31 19.5 40.6 38.7 11.1 15.8 31.9 25.8 ...
DT::datatable(USArrests)

데이터 전처리

head(), datatable() 결과를 통해서 USArrests 데이터의 첫 열에 변수명이 없음을 알 수 있다.

따라서 첫 열에 변수명 state를 부여하고 관측치의 이름을 소문자로 수정한다.

  • 책에 있는 코드
  • dplyr::%>%, mutate 이용한 코드
# Book Code

crime <- rownames_to_column(USArrests, var = "state")
crime$state <- tolower(crime$state)
# Using "dplyr" - USArrests class가 data.frame 이라 바로 가능함.

library(dplyr)

crime <- USArrests %>%
  rownames_to_column(var = "state") %>%
  mutate(state = tolower(state))

head(crime)
##        state Murder Assault UrbanPop Rape
## 1    alabama   13.2     236       58 21.2
## 2     alaska   10.0     263       48 44.5
## 3    arizona    8.1     294       80 31.0
## 4   arkansas    8.8     190       50 19.5
## 5 california    9.0     276       91 40.6
## 6   colorado    7.9     204       78 38.7

미국 주 지도 데이터 준비하기

library(ggplot2)
states_map <- ggplot2::map_data("state")

str(states_map)
## 'data.frame':    15537 obs. of  6 variables:
##  $ long     : num  -87.5 -87.5 -87.5 -87.5 -87.6 ...
##  $ lat      : num  30.4 30.4 30.4 30.3 30.3 ...
##  $ group    : num  1 1 1 1 1 1 1 1 1 1 ...
##  $ order    : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ region   : chr  "alabama" "alabama" "alabama" "alabama" ...
##  $ subregion: chr  NA NA NA NA ...

단계 구분도 만들기

  • 살인 범죄 건수에 따라 단계를 구분하고 색을 다르게 합니다.
library(ggiraphExtra)

ggChoropleth(data = crime,
             aes(fill = Murder,      # 색의 변화 기준이 되는 변수 
                 map_id = state),    # 지역 기준 변수 == 지역 명 변수
             map = states_map)       # 배경 지도 데이터

인터랙티브 단계 구분도 만들기

library(ggiraphExtra)

ggChoropleth(data = crime,
             aes(fill = Murder,      # 색의 변화 기준이 되는 변수 
                 map_id = state),    # 지역 기준 변수 == 지역 명 변수
             map = states_map,       # 배경 지도 데이터 
             interactive = T)        # 인터랙티브 유무

11 - 2. 대한민국 시도별 인구, 결핵 환자 수 단계 구분도 만들기

구분 내용
Data kormaps2014::korpop1 (2015년 시도별 센서스 데이터), kormap1(2014년 시도별 한국 행정 지도)
Packages stringi, devtools, kormaps2014, dplyr, ggplot2, ggiraphExtra
  • “kormaps2014”는 devtools 패키지를 이용해서 다운로드

패키지 설치

install.packages(c("stringi", "devtools"))
devtools::install_github("cardiomoon/kormaps2014")

데이터 확인

library(DT)
library(kormaps2014)

# 시도별 인구 데이터 
str(changeCode(korpop1))
## 'data.frame':    17 obs. of  25 variables:
##  $ C행정구역별_읍면동     : chr  "'11" "'21" "'22" "'23" ...
##  $ 행정구역별_읍면동      : chr  "서울특별시" "부산광역시" "대구광역시" "인천광역시" ...
##  $ 시점                   : chr  "2015" "2015" "2015" "2015" ...
##  $ 총인구_명              : chr  "9904312" "3448737" "2466052" "2890451" ...
##  $ 남자_명                : chr  "4859535" "1701347" "1228511" "1455017" ...
##  $ 여자_명                : chr  "5044777" "1747390" "1237541" "1435434" ...
##  $ 내국인_계_명           : chr  "9567196" "3404667" "2436770" "2822601" ...
##  $ 내국인_남자_명         : chr  "4694317" "1675339" "1211219" "1414793" ...
##  $ 내국인_여자_명         : chr  "4872879" "1729328" "1225551" "1407808" ...
##  $ 외국인_계_명           : chr  "337116" "44070" "29282" "67850" ...
##  $ 외국인_남자_명         : chr  "165218" "26008" "17292" "40224" ...
##  $ 외국인_여자_명         : chr  "171898" "18062" "11990" "27626" ...
##  $ 가구_계_가구           : chr  "3914820" "1348315" "937573" "1066297" ...
##  $ 일반가구_가구          : chr  "3784490" "1335900" "928528" "1045417" ...
##  $ 집단가구_가구          : chr  "2261" "686" "574" "713" ...
##  $ 외국인가구_가구        : chr  "128069" "11729" "8471" "20167" ...
##  $ 주택_계_호             : chr  "2793244" "1164352" "738100" "942244" ...
##  $ 단독주택_호            : chr  "355039" "225697" "155801" "102914" ...
##  $ 아파트_호              : chr  "1636896" "738068" "509068" "577346" ...
##  $ 연립주택_호            : chr  "117235" "32120" "9381" "21589" ...
##  $ 다세대주택_호          : chr  "654372" "154253" "53098" "232346" ...
##  $ 비거주용_건물내_주택_호: chr  "29702" "14214" "10752" "8049" ...
##  $ 주택이외의_거처_호     : chr  "150951" "50810" "15304" "39964" ...
##  $ C행정구역별            : chr  "11" "21" "22" "23" ...
##  $ code                   : chr  "11" "21" "22" "23" ...
# 지도 데이터 
str(changeCode(kormap1))
## 'data.frame':    8831 obs. of  15 variables:
##  $ id       : chr  "0" "0" "0" "0" ...
##  $ long     : chr  "137.774352627938" "137.779270931415" "137.780545929866" "137.814504843261" ...
##  $ lat      : chr  "50.6883045072662" "50.6899249663447" "50.6900586920365" "50.6937941360883" ...
##  $ order    : chr  "1" "2" "3" "4" ...
##  $ hole     : chr  "FALSE" "FALSE" "FALSE" "FALSE" ...
##  $ piece    : chr  "1" "1" "1" "1" ...
##  $ group    : chr  "0.1" "0.1" "0.1" "0.1" ...
##  $ SP_ID    : chr  "0" "0" "0" "0" ...
##  $ SIDO_CD  : chr  "11" "11" "11" "11" ...
##  $ SIDO_NM  : chr  NA NA NA NA ...
##  $ BASE_YEAR: chr  "2014" "2014" "2014" "2014" ...
##  $ name     : chr  "서울특별시" "서울특별시" "서울특별시" "서울특별시" ...
##  $ name1    : chr  NA NA NA NA ...
##  $ region   : chr  "11" "11" "11" "11" ...
##  $ code     : chr  "11" "11" "11" "11" ...

데이터 전처리

  • 분석에 사용할 변수 : 총인구_명, 행정구역별_읍면동 -> pop, name 으로 변경
  • 책의 코드와 11 - 1 과 유사하게 dplyr 패키지를 이용한 코드를 첨부했습니다.
# Book Code

korpop1 <- rename(korpop1,
                  pop  = 총인구_명,
                  name = 행정구역별_읍면동)

str(korpop1)
# Using "dplyr"

library(dplyr)

korpop1 <- korpop1 %>%
  rename(pop  = 총인구_명,
         name = 행정구역별_읍면동)

인터랙티브 그래프 - 대한민국 시도별 인구 밀도

  • korpop1, kormap1 사용
library(kormaps2014)
library(ggiraphExtra)
library(ggplot2)

# options(encoding = "UTF-8") 한글이 깨지므로 인코딩을 UTF-8 로 변경 

ggChoropleth(data = korpop1,       # 지도에 표현할 데이터
             aes(fill    = pop,    # 색깔로 표현할 변수
                 map_id  = code,   # 지역 기준 변수
                 tooltip = name),  # 지도 위에 표시할 지역명
             map = kormap1,        # 지도 데이터
           interactive = T)        # 인터랙티브

인터랙티브 그래프 - 대한민국 시도별 결핵 환자 수

  • tbc, kormap1 데이터 사용
library(kormaps2014)
library(ggiraphExtra)
library(ggplot2)

# options(encoding = "UTF-8") 한글이 깨지므로 인코딩을 UTF-8 로 변경 

ggChoropleth(data = tbc,            # 지도에 표현할 데이터
             aes(fill    = NewPts,  # 색깔로 표현할 변수
                 map_id  = code,    # 지역 기준 변수
                 tooltip = name),   # 지도 위에 표시할 지역명
             map = kormap1,         # 지도 데이터
             interactive = T)       # 인터랙티브

Ch 12. 인터랙티브 그래프 (Interactive Graph)

구분 내용
Data ggplot2::mpg, ggplot2::diamonds, ggplot2::economics
Packages ggplot2, plotly, dygraphs

인터랙티브 그래프 생성 과정

    1. ggplot2로 그래프 생성
    1. 위 그래프를 plotly::ggplotly에 적용

12 - 1. plotly 패키지로 인터랙티브 그래프 만들기

library(ggplot2)
library(plotly)

p1 <- mpg %>%
  ggplot(aes(displ, hwy,         # x축 배기량, y축 고속도로연비
             col = drv)) +       # drv(구동방식)에 따라 색을 다르게
  geom_point()                   # 산점도

plotly::ggplotly(p1)
p2 <- diamonds %>%
  ggplot(aes(cut,                # 컷팅 방식
             fill = clarity)) +  # 가격
  geom_bar(position = "dodge")   # 컷팅, 가격별로 막대그래프 배열 

plotly::ggplotly(p2)

12 - 2. dygraphs 패키지로 인터랙티브 시계열 데이터 만들기

library(ggplot2)
library(dygraphs)
library(xts)

eco <- xts(economics$unemploy,          # 실업자 수
           order.by = economics$date)   # 날짜순으로 순서 부여

eco %>%
  dygraph() %>%
  dyRangeSelector()

참고1) dygraphs 패키지를 이용한 비트코인 시계열 데이터 인터랙티브 그래프

xwMOOC 금융(Finance) - 암호화폐

Ch 13. 가설검정

13 - 1. 통계적 가설 검정이란? (Statistical Hypothesis Testing)

Do it R 교재 내용에 부가설명을 추가해서 정리했습니다.

통계 분석

  • 기술 통계 : 요약 통계량을 이용한 통계 분석기법(평균월급 비교 etc)
  • 추론 통계 : 집단의 특성, 분포(이항, 포아송 등)에 맞게 통계 분석 실시하여 가설 검정

13 - 2. t-test - 두 집단의 평균 비교 (Two-Sample t-test)

  • t-검정(t-test)은 검정통계량이 스튜던트의 t-분포를 따르는 통계절차를 말한다.
  • “두 집단의 평균 비교”는 분산분석(ANOVA)의 특이한 예시중에 하나이다.

ggplot2::mpg 데이터를 사용한 t-test 과정

    1. mpg 데이터에서 class, cty 변수만 추출
    1. class==“compact”, “suv” 인 데이터 추출 : 두 집단을 compact와 suv로 선택
    1. t.test() 실행
library(ggplot2)
library(dplyr)
library(descr)

# 데이터 취득
mpg <- tbl_df(ggplot2::mpg)

mpg_diff <- mpg %>%
  select(class, cty) %>%
  filter(class %in% c("compact", "suv"))

# 데이터 확인
head(mpg_diff)
## # A tibble: 6 x 2
##     class   cty
##     <chr> <int>
## 1 compact    18
## 2 compact    21
## 3 compact    20
## 4 compact    21
## 5 compact    16
## 6 compact    18
descr::CrossTable(mpg_diff$class)
##    Cell Contents 
## |-------------------------|
## |                       N | 
## |           N / Row Total | 
## |-------------------------|
## 
## | compact |     suv |
## |---------|---------|
## |      47 |      62 |
## |   0.431 |   0.569 |
## |---------|---------|

t-test 실습1 : 교재에 생략되어있는 가설 설정부터 실행

H0(귀무가설) : class(차종)에 따른 cty(도시연비)차이가 통계적으로 유의하지 않다.

(class가 cty에 영향을 주지 않는다.)

H1(대립가설) : class(차종)에 따른 cty(도시연비)차이가 통계적으로 유의하다

(class가 cty에 영향을 준다.)

t-test 를 통한 가설 검정

t.test(data = mpg_diff,   # 검정에 사용할 데이터 지정
       cty ~ class,       # 종속변수 ~ 독립변수 / y = x 형태라고 이해하면 된다.
       var.equal = T)     # 집단 간 분산이 같다고 가정함.
## 
##  Two Sample t-test
## 
## data:  cty by class
## t = 11.917, df = 107, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  5.525180 7.730139
## sample estimates:
## mean in group compact     mean in group suv 
##              20.12766              13.50000

p-value < 2.2e-16 이므로 “집단 간 차이가 통계적으로 유의”하다

혹은 “귀무가설 반박, 대립가설 채택”이라고 결론을 내림.

즉 class가 cty에 영향을 준다

단, p-value 만을 보고 판단하는 것은 좋지 않다.

t-test 실습2

  • 일반 휘발유와 고급 휘발유 두 집단의 cty t-test

H0(귀무가설) : fl(휘발유)에 따른 cty(도시연비)차이가 통계적으로 유의하지 않다.

(fl가 cty에 영향을 주지 않는다.)

H1(대립가설) : fl(휘발유)에 따른 cty(도시연비)차이가 통계적으로 유의하다

(fl가 cty에 영향을 준다.)

mpg_diff2 <- mpg %>%
  select(fl, cty) %>%
  filter(fl %in% c("r", "p"))   # r : regular / p : premium

head(mpg_diff2)
## # A tibble: 6 x 2
##      fl   cty
##   <chr> <int>
## 1     p    18
## 2     p    21
## 3     p    20
## 4     p    21
## 5     p    16
## 6     p    18
descr::CrossTable(mpg_diff2$fl)
##    Cell Contents 
## |-------------------------|
## |                       N | 
## |           N / Row Total | 
## |-------------------------|
## 
## |       p |       r |
## |---------|---------|
## |      52 |     168 |
## |   0.236 |   0.764 |
## |---------|---------|
t.test(data = mpg_diff2,
       cty ~ fl,
       var.equal = T)
## 
##  Two Sample t-test
## 
## data:  cty by fl
## t = 1.0662, df = 218, p-value = 0.2875
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -0.5322946  1.7868733
## sample estimates:
## mean in group p mean in group r 
##        17.36538        16.73810

p-value < 2.2e-16 이므로 “집단 간 차이가 통계적으로 유의하다”

혹은 “귀무가설 반박, 대립가설 채택”이라고 결론을 내림.

13 - 3. 상관분석(Correlation Analysis) : 두 연속형 변수의 관계성 분석

unemploy : 실업자 수, pce : 개인 소비 지출

economics <- tbl_df(ggplot2::economics)
cor.test(economics$unemploy, economics$pce)
## 
##  Pearson's product-moment correlation
## 
## data:  economics$unemploy and economics$pce
## t = 18.605, df = 572, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.5603164 0.6625460
## sample estimates:
##       cor 
## 0.6139997

cor.test의 p-vlaue < 2.2e-16 이므로 두 연속형 변수는 서로 연관이 있음을 알 수 있고

cor(상관계수) = 0.61 이므로 정비례 관계임을 알 수 있다.

(실업자 수가 증가하면 개인 소비 지출도 증가한다.)

참고로 unemploy 와 pce 는 각각 num, int 데이터임.

따라서 산점도를 그려서 시각적으로 상관성을 확인 할 수 있음.

economics 데이터는 시계열 데이터이므로 일정 기간을 잡아서 시각화 해야함.

library(ggplot2)
library(dplyr)

# 원하는 기간 데이터 생성 
s <- as.Date("2012-01-01")
e <- as.Date("2012-12-01")

date_month <- seq(from = s,
                  to   = e,
                  by   = 1)

# 매칭되는 데이터만 필터링 후 산점도로 시각화 
economics %>%
  filter(date %in% date_month) %>%
  ggplot(aes(unemploy, pce)) +
  geom_point()

상관행렬 히트맵 만들기

  • R base data 인 mtcars 사용
head(mtcars)
##                    mpg cyl disp  hp drat    wt  qsec vs am gear carb
## Mazda RX4         21.0   6  160 110 3.90 2.620 16.46  0  1    4    4
## Mazda RX4 Wag     21.0   6  160 110 3.90 2.875 17.02  0  1    4    4
## Datsun 710        22.8   4  108  93 3.85 2.320 18.61  1  1    4    1
## Hornet 4 Drive    21.4   6  258 110 3.08 3.215 19.44  1  0    3    1
## Hornet Sportabout 18.7   8  360 175 3.15 3.440 17.02  0  0    3    2
## Valiant           18.1   6  225 105 2.76 3.460 20.22  1  0    3    1
# 교재는 head() 만 사용하여 데이터를 확인했지만
# 모든 변수가 연속형인지 판단하기 위해 str()을 이용하여 한 번 더 확인한다.
str(mtcars)
## 'data.frame':    32 obs. of  11 variables:
##  $ mpg : num  21 21 22.8 21.4 18.7 18.1 14.3 24.4 22.8 19.2 ...
##  $ cyl : num  6 6 4 6 8 6 8 4 4 6 ...
##  $ disp: num  160 160 108 258 360 ...
##  $ hp  : num  110 110 93 110 175 105 245 62 95 123 ...
##  $ drat: num  3.9 3.9 3.85 3.08 3.15 2.76 3.21 3.69 3.92 3.92 ...
##  $ wt  : num  2.62 2.88 2.32 3.21 3.44 ...
##  $ qsec: num  16.5 17 18.6 19.4 17 ...
##  $ vs  : num  0 0 1 1 0 1 0 1 1 1 ...
##  $ am  : num  1 1 1 0 0 0 0 0 0 0 ...
##  $ gear: num  4 4 4 3 3 3 3 4 4 4 ...
##  $ carb: num  4 4 1 1 2 1 4 2 2 4 ...

위의 결과를 통해 mtcars의 속성이 data.frame, 모든 변수가 numeric 임을 알 수 있다.

따라서 아래 두 가지 방법을 이용하여 상관행렬을 생성할 수 있다.

# Book Code

car_cor <- cor(mtcars)
round(car_cor, 2)
# Using "dplyr"

library(dplyr)

mtcars %>%
  cor() %>%
  round(2)
##        mpg   cyl  disp    hp  drat    wt  qsec    vs    am  gear  carb
## mpg   1.00 -0.85 -0.85 -0.78  0.68 -0.87  0.42  0.66  0.60  0.48 -0.55
## cyl  -0.85  1.00  0.90  0.83 -0.70  0.78 -0.59 -0.81 -0.52 -0.49  0.53
## disp -0.85  0.90  1.00  0.79 -0.71  0.89 -0.43 -0.71 -0.59 -0.56  0.39
## hp   -0.78  0.83  0.79  1.00 -0.45  0.66 -0.71 -0.72 -0.24 -0.13  0.75
## drat  0.68 -0.70 -0.71 -0.45  1.00 -0.71  0.09  0.44  0.71  0.70 -0.09
## wt   -0.87  0.78  0.89  0.66 -0.71  1.00 -0.17 -0.55 -0.69 -0.58  0.43
## qsec  0.42 -0.59 -0.43 -0.71  0.09 -0.17  1.00  0.74 -0.23 -0.21 -0.66
## vs    0.66 -0.81 -0.71 -0.72  0.44 -0.55  0.74  1.00  0.17  0.21 -0.57
## am    0.60 -0.52 -0.59 -0.24  0.71 -0.69 -0.23  0.17  1.00  0.79  0.06
## gear  0.48 -0.49 -0.56 -0.13  0.70 -0.58 -0.21  0.21  0.79  1.00  0.27
## carb -0.55  0.53  0.39  0.75 -0.09  0.43 -0.66 -0.57  0.06  0.27  1.00

상관행렬 히트맵을 이용한 시각화

  • 교재와 다르게 상관행렬을 저장하지않고 dplyr의 파이프 연산자를 이용하여 바로 시각화
library(corrplot)

mtcars %>%
  cor() %>%
  corrplot()

원 대신 상관계수 출력

mtcars %>%
  cor() %>%
  corrplot(method = "number")

색상, 상관계수 출력

col <- colorRampPalette(c("#BB4444", "#EE9988",
                          "#FFFFFF", "#77AADD",
                          "#4477AA"))

mtcars %>%
  cor() %>%
  corrplot(method      = "color",    # 색깔로 표현
           col         = col(200),   # 색상 200개 선정
           type        = "lower",    # 왼쪽 아래 행렬만 표시
           order       = "hclust",   # 유사한 상관계수끼리 군집화
           addCoef.col = "black",    # 상관계수 색깔
           tl.col      = "black",    # 변수명 색깔
           tl.srt      = 45,         # 변수명 45도 기울임
           diag        = F)          # 대각 행렬 제외

참고 1) 데이터 시각화 - pairwise 산점도 + 상관계수

library(ggplot2)
library(mycor)
library(dplyr)

pairs(mtcars,
      lower.panel = function(x,y){ points(x,y); abline(0, 1, col = "red")},
      upper.panel = panel.cor)

# mpg 데이터에서 연속형 변수들만 select,
# 하삼각행렬에 산점도와 상관계수 1인 직선 그래프
# 상삼각행렬에 각 변수들의 상관계수 출력
pairs(mpg %>%
        select(displ, year, cyl, cty, hwy),
      lower.panel = function(x,y){ points(x,y); abline(0, 1, col = "red")},
      upper.panel = panel.cor)

참고 2) Statistical Hypothesis Testing Algorithm - 통계적 가설 검정 절차

모집단의 특성에 따른 검정 절차 및 방법은 아래 사진과 같다.

출처 : POSTECH MOOC - 데이터 분석을 위한 R 프로그래밍

1. Do it 쉽게 배우는 R 데이터 분석

2. 따라하며 배우는 데이터 과학

3. 초보자를 위한 RStudio 마스터

4. POSTECH MOOC - 데이터 분석을 위한 R 프로그래밍