Chapter 6 Tables for public health 1

6.1 데이터표 준비

데이터 표를 만드는 실습은 6차 근로환경조사 자료를 통해 실습할 것입니다.. 자료는 안전보건공단, 근로환경조사 원시자료 사이트 (http://kosha.or.kr/kosha/data/primitiveData.do) 에서 신청할 수 있습니다.. 시각화는 표와 그래프 둘을 다룰 것이고, 표는 보건학에서 주로 사용하는 표1,2,3과 그래프는 막대 그래프부터 시계열 자료까지 나타내도록 하겠습니다. 실습에 필요한 라이브러리를 불러오겠습니다.

if(!require("tidyverse")) install.packages("tidyverse")
if(!require("htmlTable")) install.packages("htmlTable")
if(!require("haven")) install.packages("haven")

데이터를 불러오겠습니다. 안전보건공단 홈페이에서 자료를 다운 받는게 원칙입니다. 다만 실습을 빠르게 진행하기 위해서, dspubs.org 페이지에 있는 파일을 이용하겠습니다.
dspubs.org open_data

url <- "https://dspubs.org/webapps/forum/open_data/kwcs6th.sav"
download.file(url, "data/kwcs6th.sav")
kwcs <- read_sav("data/kwcs6th.sav")

unicode chart

표를 만들다 보면 크다, 작다, 같다 등은 표시가 쉬운 반면, 크거나 같다, 작거나 같다 등은 표시가 어렵습니다. 이 때 사용해야할 것이 unicode 입니다. 보건학에서 필요한 대표적 유니코드는 사실 2개입니다. 크거나 같다, 작거나 같다. 나머지는 키보드에 이미 있으니 이것을 사용하면 됩니다. 추가적인 unicode는 아래의 항목을 통해 살펴 볼수 있습니다.
Unicode Chart

tibble(
  "symbole    " = c("\u2264", "\u2265", "\u00b1"), 
  "unicode    " = c("u2264", "u2265", "u00b1")
) %>% 
  addHtmlTableStyle(css.cell = c("width: 100;","width: 100;")) %>%
  htmlTable(caption ="Simple unicode and symbols")
Simple unicode and symbols
symbole unicode
1 u2264
2 u2265
3 ± u00b1

다음과 같이 사용할 수 있습니다. 여기사 \는 “escape character”로 뒤에 오는 것이 문자가 아니라 약솓된 결과를 나타내 달라는 뜻입니다.

print("x \u2264 10")
## [1] "x ≤ 10"

그럼 어떻게 문자를 사용할까요? print 명령을 위와 같이 사용하는 것도 좋지만, 변수를 생성하기에는 좋은 방법이 따로 있습니다. 좋은 방법이란 처음에는 어렵지만, 알고나면 엄청난 것들을 말합니다. paste와 sprintf 입니다. 어떤 것이 좋은가요?. 우리의 목표가 데이터 표현을 자동화 하는 것인데, 자동화를 위해서는 구조를 만들고 거기에 값을 대입 시키는 것이 기본입니다. 그러니, sprintf 를 더 자주 사용하게 됩니다. paste는 모두 붙여 주는 것이니, 쉽게 이해될 것이고, sprintf는 %s 마다 무언가를 넣어 붙여 주는 방식입니다. paste와 sprintf는 정말 자주 사용되는 함수이니 꼭 익숙해 지시기 바랍니다.

xp1 = paste("x", "\u2264", "10")
xp2 = sprintf("x %s 10", "\u2264")
xp1
## [1] "x ≤ 10"
xp2
## [1] "x ≤ 10"

데이터 확인

데이터의 변수를 확인하는 방법은 colnames() , names() 또는 head()를 하용하는 것입니다. 약 400개의 변수가 있으므로, 처음부터 10번째까지 [1:10] 변수를 찾아 보겠습니다.

colnames(kwcs)[1:10]
##  [1] "id"           "wt"           "area"         "hh_num"       "hm_01_gender"
##  [6] "hm_01_year"   "hm_01_estat"  "hm_01_rel_t"  "hm_02_gender" "hm_02_year"

그럼 45번째 변수 부터 50번째 변수까지 찾아 보겠습니다. []안을 채워보세요.

colnames(kwcs)[     ]
## [1] "target"      "YEAR"        "ESTAT"       "AGE"         "country"    
## [6] "country_etc"

연령인 AGE가 있네요, 성별에 대한 AGE도 있습니다. 그럼 이것을 이용해서 실습을 해보겠습니다. 변수를 하나 선택하는 것은 데이터에 $ 표시를 하고 이후에 변수를 넣는 방식입니다. “데이터\(변수" 입니다. 10개만 보겠습니다. 숫자 처럼 보이네요, 확인하겠습니다. `class(kwcs\)AGE)`를 이용해봅니다. numberic 으로 숫자입니다. 숫자여야 평균 표준편차 등의 계산이 가능합니다.

kwcs$AGE[1:10]
##  [1] 74 54 64 78 78 65 61 57 38 47
class(kwcs$AGE)
## [1] "numeric"

이번에는 SEX 변수(variable)의 변수값(value)를 살펴 보겠습니다. class 가 무엇일까요? double 또는 interger 라고 나올 텐데요, R에서는 numberic 값에 interger와 double 속성을 사용한다고 생각하시면 됩니다. 다만 label을 붙여 놓아서 알기 쉽게 되어 있네요.

kwcs$TSEX[1:10]
## <labelled<double>[10]>: ▣ 분석 데모변수 : 응답자 성별
##  [1] 2 2 2 1 1 2 1 2 1 2
## 
## Labels:
##  value label
##      1  남성
##      2  여성
class(kwcs$TSEX)
## [1] "haven_labelled" "vctrs_vctr"     "double"

국적에 대해서도 알아보겠습니다. character 이네요. 한국에서 시행한 조사라, 한국 국적은 빈칸으로, 이외에는 국적을 적었습니다. kwcs$country_etc를 해보면, 대부분 빈칸입니다. 빈칸인 경 =="", 빈칸이 아닌 경우 !=""을 이용해서 어떤 국적이 있는지 살펴 보겠습니다.

class(kwcs$country_etc)
## [1] "character"
kwcs$country_etc[kwcs$country_etc != ""][1:10]
##  [1] "베트남"   "중국"     "캐나다"   "중국"     "중국"     "중국"    
##  [7] "중국"     "중국"     "파키스탄" "중국"

근로자의 지위를 알아 보겠습니다. 근로자 지위는 emp_stat 입니다. 1은 상용근로자, 2는 임시근로자, 3은 일용근로자 입니다.

kwcs$emp_stat 
kwcs$emp_stat %>% head()
## <labelled<double>[6]>: Q6.임금근로자_고용형태
## [1] NA  1  3 NA NA  1
## 
## Labels:
##  value      label
##      1 상용근로자
##      2 임시근로자
##      3 일용근로자

매우 중요한 개념이 나옵니다. 빈칸과 NA 입니다. 모두 값에 대한 정보가 없다는 것입니다. 값에 대한 정보가 없으면 어떻게 해야 할지는 매우 중요한 개념입니다. 우선 여기서는 값에 대한 정보가 없는 것을 제외하고 분석해 보겠습니다. 값에 대한 정보가 없는 데이터는 제거하여, 새로운 데이터를 만들겠습니다. 이때 filter 라는 것을 이용합니다. is.na라는 것은 NA라는 것을 의미하고, 앞에 !는 그 반대를 말합니다.

kwcs %>%
  filter(!is.na(emp_stat))
## # A tibble: 33,063 × 438
##       id    wt      area hh_num hm_01_gender hm_01_year  hm_01_estat hm_01_rel_t
##    <dbl> <dbl> <dbl+lbl>  <dbl>    <dbl+lbl>      <dbl>    <dbl+lbl>   <dbl+lbl>
##  1     2 0.137  3 [대구]      1     2 [여성]       1966 1 [임금 근… 0 [응답자 …
##  2     3 0.350  3 [대구]      1     2 [여성]       1956 1 [임금 근… 0 [응답자 …
##  3     6 0.156  7 [울산]      1     2 [여성]       1955 2 [일이 있… 0 [응답자 …
##  4     7 0.344  7 [울산]      1     1 [남성]       1959 1 [임금 근… 0 [응답자 …
##  5     9 0.168  7 [울산]      1     1 [남성]       1982 1 [임금 근… 0 [응답자 …
##  6    10 0.675  6 [대전]      4     2 [여성]       1973 1 [임금 근… 0 [응답자 …
##  7    18 2.71   6 [대전]      4     1 [남성]       1972 1 [임금 근… 0 [응답자 …
##  8    19 3.83   3 [대구]      5     1 [남성]       1969 1 [임금 근… 3 [부모]  …
##  9    20 2.90   3 [대구]      4     1 [남성]       1960 1 [임금 근… 3 [부모]  …
## 10    21 2.09   3 [대구]      3     2 [여성]       1940 1 [임금 근… 0 [응답자 …
## # … with 33,053 more rows, and 430 more variables: hm_02_gender <dbl+lbl>,
## #   hm_02_year <dbl>, hm_02_estat <dbl+lbl>, hm_02_rel_t <dbl+lbl>,
## #   hm_03_gender <dbl+lbl>, hm_03_year <dbl>, hm_03_estat <dbl+lbl>,
## #   hm_03_rel_t <dbl+lbl>, hm_04_gender <dbl+lbl>, hm_04_year <dbl>,
## #   hm_04_estat <dbl+lbl>, hm_04_rel_t <dbl+lbl>, hm_05_gender <dbl+lbl>,
## #   hm_05_year <dbl>, hm_05_estat <dbl+lbl>, hm_05_rel_t <dbl+lbl>,
## #   hm_06_gender <dbl+lbl>, hm_06_year <dbl>, hm_06_estat <dbl+lbl>, …

요통에 대한 변수는 heal_prob1 입니다. 이 변수를 확인해 보겠습니다. 1번은 있다, 2번은 업다, 8번은 무응답, 9번은 거절입니다. 이제 빈칸과 NA가 아니더라도 필요없는 정보가 있습니다. 모르거나/무응답했거나, 거절한 사람입니다. 이를 제거해 보겠습니다.

kwcs %>% 
  filter(!is.na(emp_stat)) %>%
  filter(!is.na(heal_prob1)) %>%
  filter(heal_prob1 !=8) %>%
  filter(heal_prob1 !=9)

이번에는 heal_prob1 에 1, 2 인 사람만 포함시켜 보겠습니다. 어떤 것이 더 편한가요?

kwcs %>% 
  filter(!is.na(emp_stat)) %>%
  filter(!is.na(heal_prob1)) %>%
  filter(heal_prob1 %in% c(1, 2))

이번에는 sleep1이라는 변수를 살펴 보겠습니다. sleep1이라는 변수는 잠들기 어려운 것이 매일(1), 한주에 여러번 (2), 한달에 여러번(3), 드물게(5), 전혀 없음(5) 의 5점 척도 입니다. 이데 모름과 거절, 그리고 NA 값을 제거해 보겠습니다.

그리고 dat라는 새로운 data를 만들고 kwcs를 제거 하겠습니다. dat를 이용해서 분석을 해보겠습니다. 지금까지 사용했던, TSEX, AGE, emp_stat, heal_prob1, sleep1 의 변수를 사용하겠습니다.

dat <- kwcs %>% 
  filter(!is.na(emp_stat)) %>%
  filter(!is.na(heal_prob1)) %>%
  filter(heal_prob1 %in% c(1, 2)) %>%
  filter(sleep1 %in% c(1:5)) %>%
  select(TSEX, AGE, emp_stat, heal_prob1, sleep1, sleep2, sleep3)
rm(kwcs)

dat를 살펴보겠습니다.

head(dat)
## # A tibble: 6 × 7
##        TSEX   AGE       emp_stat heal_prob1              sleep1   sleep2  sleep3
##   <dbl+lbl> <dbl>      <dbl+lbl>  <dbl+lbl>           <dbl+lbl> <dbl+lb> <dbl+l>
## 1  2 [여성]    54 1 [상용근로자]   2 [없다] 5 [전혀없음]        5 [전혀… 5 [전…
## 2  2 [여성]    64 3 [일용근로자]   2 [없다] 2 [한 주에 여러 번] 2 [한 … 2 [한 …
## 3  2 [여성]    65 1 [상용근로자]   2 [없다] 5 [전혀없음]        5 [전혀… 5 [전…
## 4  1 [남성]    61 3 [일용근로자]   1 [있다] 4 [드물게]          4 [드물… 4 [드…
## 5  1 [남성]    38 1 [상용근로자]   2 [없다] 5 [전혀없음]        5 [전혀… 5 [전…
## 6  2 [여성]    47 1 [상용근로자]   2 [없다] 5 [전혀없음]        4 [드물… 5 [전…

6.2 Central Tendency (mean, median, mode)

대표값중 가장 많이 사용하는 것은 mean 과 median 입니다. 이것을 나타내는 표를 만들어 봅시다. 이것의 평균과 표준 편차를 을 구해 보겠습니다.

mean(dat$AGE)
## [1] 46.564
sd(dat$AGE)
## [1] 14.45298

이번에는 “tidyverse”를 통해 pipe 코드를 짜 보겠습니다. 데이터를 변형하는 것은 다른 시간에 수행하겠지만, 여기서는 select와 filter, mutate, group_by를 사용하겠습니다. kwcs$TSEX 는 kwcs에서 TSEX를 select하라는 것으로 다름과 같이 사용할 수 있습니다.

dat %>% select(AGE) 
## # A tibble: 33,007 × 1
##      AGE
##    <dbl>
##  1    54
##  2    64
##  3    65
##  4    61
##  5    38
##  6    47
##  7    48
##  8    23
##  9    28
## 10    80
## # … with 32,997 more rows

6.2.1 하나씩 반복

이후 이것을 가지고 나와서 (pull), 이어 받고 (.) 평균을 구해보겠습니다. 이후 이것을 mean 과 sd 라는 변수에 assign 하겠습니다.

dat %>% select(AGE) %>% pull(.) %>% mean(.)
## [1] 46.564
dat %>% select(AGE) %>% pull(.) %>% sd(.)
## [1] 14.45298
dat %>% select(AGE) %>% pull(.) %>% mean(.) -> mean
dat %>% select(AGE) %>% pull(.) %>% sd(.)   -> sd

그럼 표현해 볼까요?

mean
## [1] 46.564
sd
## [1] 14.45298
paste("평균은 ", mean, "표준편차는 ", sd)
## [1] "평균은  46.5640015754234 표준편차는  14.4529846074525"

보기 불편하네요, 소숫점 2째 자리까지 표현하겠습니다. round를 이용합니다. 더 자세한것은 구글에게 물어 보세요.

mean <- dat %>% select(AGE) %>% pull(.) %>% mean(.) %>% round(., 2) 
sd   <- dat %>% select(AGE) %>% pull(.) %>% sd(.)   %>% round(., 2)
paste("평균은", mean, ", 표준편차는", sd)
## [1] "평균은 46.56 , 표준편차는 14.45"
sprintf("평균은 %s, 표준편차는 %s", mean, sd)
## [1] "평균은 46.56, 표준편차는 14.45"

unicode를 이용해서 약속된 표현을 사용해 보겠습니다.

paste(mean, "\u00b1", sd)
## [1] "46.56 ± 14.45"
sprintf("%s \u00b1 %s", mean, sd)
## [1] "46.56 ± 14.45"

sprintf 에서 %s 대신에 %.2f 를 사용할 텐데요, 어떤지 살펴 봅시다. 네 %2.f 라는 것은 소수 2째 자리까지 살려서 표현하는 방식입니다. %.2f 는 어떨까요? 공부하는 방법입니다. 무언가를 더 해보는 것!

sprintf("%.2f \u00b1 %.2f", mean, sd)
## [1] "46.56 ± 14.45"

그럼 median 을 구해볼까요? 해보세요. 무언가를 해보는 것!

dat %>% pull(AGE) %>% median(.)
## [1] 46
dat %>% pull(AGE) %>% quantile(., c(0.5))
## 50% 
##  46

그럼 quantile을 구해볼까요?

dat %>% pull(AGE) %>% quantile(., c(0.25, 0.5, 0.75)) 
## 25% 50% 75% 
##  35  46  57
p50 = dat %>% pull(AGE) %>% quantile(., c(0.50))
p25 = dat %>% pull(AGE) %>% quantile(., c(0.25)) 
p75 = dat %>% pull(AGE) %>% quantile(., c(0.75)) 
sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)
## [1] "46 (35-57)"

이것을 남녀를 나누어서 해보겠습니다. filter 명령문을 써보겠습니다. filter(TSEX==1)이라는 것은 TSEX==1 인 남자만을 고르라는 것입니다. 평균은 아래와 같이 구합니다. sd도 구해봅시다.

dat %>% filter(TSEX==1) %>% pull(AGE) %>% mean(.)
## [1] 45.30149
dat %>% filter(TSEX==2) %>% pull(AGE) %>% mean(.)
## [1] 47.67811

6.2.2 Group_by summary and Table 1

Group_by

남녀를 나누어서 평균을 구하는 다른 방법을 사용해 봅니다.

dat %>%
  group_by(TSEX) %>%
  summarise(avg = mean(AGE), 
            std = sd(AGE))
## # A tibble: 2 × 3
##        TSEX   avg   std
##   <dbl+lbl> <dbl> <dbl>
## 1  1 [남성]  45.3  13.9
## 2  2 [여성]  47.7  14.8

동일한 결과가 나오나요. 네 굉장합니다. 남녀를 나누어서 계산했네요. 지금은 2개의 집단을 나누지만, 만약 100개의 집단이라면 filter를 반복하면 어떻게 해야하나요ㅠㅠ. group_by 는 정말 대단한 명령어 입니다. group_by로 코드 파이프 안에서 다음과 같은 결과를 얻었습니다.

dat %>%
  group_by(TSEX) %>%
  summarise(avg = mean(AGE), 
            std = sd(AGE)) %>%
  mutate(smry = sprintf("%.2f \u00b1 %.2f", avg, std))
## # A tibble: 2 × 4
##        TSEX   avg   std smry         
##   <dbl+lbl> <dbl> <dbl> <chr>        
## 1  1 [남성]  45.3  13.9 45.30 ± 13.94
## 2  2 [여성]  47.7  14.8 47.68 ± 14.80

그럼 중간값과 p25-p75를 표현해 봅시다.

dat %>%
  group_by(TSEX) %>%
  summarise(avg = mean(AGE), 
            std = sd(AGE), 
            p25 = quantile(AGE, prob=c(0.25)), 
            p50 = quantile(AGE, prob=c(0.50)), 
            p75 = quantile(AGE, prob=c(0.75)), 
            ) %>%
  mutate(smry1 = sprintf("%.1f \u00b1 %.1f", avg, std)) %>%
  mutate(smry2 = sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)) 
## # A tibble: 2 × 8
##        TSEX   avg   std   p25   p50   p75 smry1       smry2     
##   <dbl+lbl> <dbl> <dbl> <dbl> <dbl> <dbl> <chr>       <chr>     
## 1  1 [남성]  45.3  13.9    34    44    55 45.3 ± 13.9 44 (34-55)
## 2  2 [여성]  47.7  14.8    36    48    57 47.7 ± 14.8 48 (36-57)

이제 필요한 것만 남겨 보겠습니다.

dat %>%
  group_by(TSEX) %>%
  summarise(avg = mean(AGE), 
            std = sd(AGE), 
            p25 = quantile(AGE, prob=c(0.25)), 
            p50 = quantile(AGE, prob=c(0.50)), 
            p75 = quantile(AGE, prob=c(0.75)), 
            ) %>%
  mutate(smry1 = sprintf("%.1f \u00b1 %.1f", avg, std)) %>%
  mutate(smry2 = sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)) %>%
  select(TSEX, smry1, smry2)
## # A tibble: 2 × 3
##        TSEX smry1       smry2     
##   <dbl+lbl> <chr>       <chr>     
## 1  1 [남성] 45.3 ± 13.9 44 (34-55)
## 2  2 [여성] 47.7 ± 14.8 48 (36-57)

이번에는 같은 내용을 수면에 대해서 해보겠습니다. 수면 점수가 높다믄 것은 잠들기 어렵다는 것이 전혀 없음(5점)에 가깝다는 것입니다. 그러니 역의 점수를 만들겠습니다. 1, 2, 3, 4, 5로 기록된 것을 5, 4, 3, 2, 1로 바꾸고 싶은 것입니다. 아래를 이용할 것입니다.

x= 1:5
y= 6-x
y
## [1] 5 4 3 2 1
dat %>%
  mutate(sleep1in = 6-sleep1)
## # A tibble: 33,007 × 8
##         TSEX   AGE       emp_stat heal_prob1     sleep1  sleep2  sleep3 sleep1in
##    <dbl+lbl> <dbl>      <dbl+lbl>  <dbl+lbl>  <dbl+lbl> <dbl+l> <dbl+l>    <dbl>
##  1  2 [여성]    54 1 [상용근로자]   2 [없다] 5 [전혀없… 5 [전… 5 [전…        1
##  2  2 [여성]    64 3 [일용근로자]   2 [없다] 2 [한 주… 2 [한 … 2 [한 …        4
##  3  2 [여성]    65 1 [상용근로자]   2 [없다] 5 [전혀없… 5 [전… 5 [전…        1
##  4  1 [남성]    61 3 [일용근로자]   1 [있다] 4 [드물게… 4 [드… 4 [드…        2
##  5  1 [남성]    38 1 [상용근로자]   2 [없다] 5 [전혀없… 5 [전… 5 [전…        1
##  6  2 [여성]    47 1 [상용근로자]   2 [없다] 5 [전혀없… 4 [드… 5 [전…        1
##  7  1 [남성]    48 1 [상용근로자]   1 [있다] 5 [전혀없… 5 [전… 5 [전…        1
##  8  1 [남성]    23 2 [임시근로자]   2 [없다] 4 [드물게… 4 [드… 4 [드…        2
##  9  2 [여성]    28 1 [상용근로자]   2 [없다] 4 [드물게… 4 [드… 4 [드…        2
## 10  2 [여성]    80 3 [일용근로자]   1 [있다] 4 [드물게… 4 [드… 4 [드…        2
## # … with 32,997 more rows

이번에는 mutate을 이용해 새로운 변수를 생성해 보겠습니다. 변형된 dat는 dat1에 할당합니다.

dat1<-dat %>%
  mutate(sleep1in = 6-sleep1, 
         sleep2in = 6-sleep2, 
         sleep3in = 6-sleep3 
         ) %>%
  mutate(sexgp = case_when(
    TSEX ==1 ~ "Men", 
    TRUE ~ "Women"
  ))

이번에는 sexgp를 group으로 하여 성별, 연령의 대표값을 나나타내 봅시다.

dat1 %>%
  group_by(sexgp) %>%
  summarise(
            avg = mean(AGE), 
            std = sd(AGE), 
            p50 = quantile(AGE, prob=c(0.50)), 
            p25 = quantile(AGE, prob=c(0.25)), 
            p75 = quantile(AGE, prob=c(0.75))
            ) %>%
  mutate(smry1= sprintf("%.1f \u00b1 %.1f", avg, std),
         smry2= sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)
         ) %>%
    select(sexgp,smry1, smry2)
## # A tibble: 2 × 3
##   sexgp smry1       smry2     
##   <chr> <chr>       <chr>     
## 1 Men   45.3 ± 13.9 44 (34-55)
## 2 Women 47.7 ± 14.8 48 (36-57)

html Table로 만들면 복사해서 붙여 넣기 매우 편리합니다.

dat1 %>%
  group_by(TSEX) %>%
  summarise(
            avg = mean(AGE), 
            std = sd(AGE), 
            p50 = quantile(AGE, prob=c(0.50)), 
            p25 = quantile(AGE, prob=c(0.25)), 
            p75 = quantile(AGE, prob=c(0.75))
            ) %>%
  mutate(smry1= sprintf("%.1f \u00b1 %.1f", avg, std),
         smry2= sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)
         ) %>%
    select(TSEX,smry1, smry2) %>%
  htmlTable(caption ="descritive statistics of study population")
descritive statistics of study population
TSEX smry1 smry2
1 1 45.3 ± 13.9 44 (34-55)
2 2 47.7 ± 14.8 48 (36-57)

그렇다면, sleep1 에 대해서도 만들어 보십시오.

연령과 sleep1 에 대한 대표값을 함께 볼 수 있을까요? 지금까지 배운 방법으로 해보도록 하겠습니다.

tab1 = dat1 %>%
  group_by(TSEX) %>%
  summarise(
            avg = mean(AGE), 
            std = sd(AGE), 
            p50 = quantile(AGE, prob=c(0.50)), 
            p25 = quantile(AGE, prob=c(0.25)), 
            p75 = quantile(AGE, prob=c(0.75))
            ) %>%
  mutate(smry1= sprintf("%.1f \u00b1 %.1f", avg, std),
         smry2= sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)
         ) %>%
    select(TSEX,smry1, smry2) %>%
  mutate(variable = "AGE")
tab2 = dat1 %>%
  group_by(TSEX) %>%
  summarise(
            avg = mean(sleep1), 
            std = sd(sleep1), 
            p50 = quantile(sleep1, prob=c(0.50)), 
            p25 = quantile(sleep1, prob=c(0.25)), 
            p75 = quantile(sleep1, prob=c(0.75))
            ) %>%
  mutate(smry1= sprintf("%.1f \u00b1 %.1f", avg, std),
         smry2= sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)
         ) %>%
    select(TSEX,smry1, smry2) %>%
  mutate(variable = "sleep1")

이 둘을 합쳐 보겠습니다.

rbind(tab1, tab2)
## # A tibble: 4 × 4
##        TSEX smry1       smry2      variable
##   <dbl+lbl> <chr>       <chr>      <chr>   
## 1  1 [남성] 45.3 ± 13.9 44 (34-55) AGE     
## 2  2 [여성] 47.7 ± 14.8 48 (36-57) AGE     
## 3  1 [남성] 4.5 ± 0.8   5 (4-5)    sleep1  
## 4  2 [여성] 4.3 ± 0.8   5 (4-5)    sleep1

sleep1까지 했는데요, sleep2, sleep3 까지 하려면 어떻게 해야 하나요? 이것을 2번더 반복해야 합니다. 만약에 변수가 100개라면 어떻게 해야할까요, 네 1000번 하면됩니다. 좀 익숙한 표를 만들기 위해서 반복하는 것이 가장 기본이 됩니다. 그래도 이제는 이 방법을 응용해서, 컴퓨터가 쉽게 작동하는 방식으로 생각해 보겠습니다.

6.2.3 Long File and Table 1

long file

이제 group 대신에 group이 될 파일을 선택해서 long file로 만들어 보겠습니다. 원하는 파일만 골라 봅니다. 성에 따른 연령, sleep1, 2, 3를 정리해 보겠습니다.

dat1 %>% select(sexgp, AGE, sleep1in, sleep2in, sleep3in)
## # A tibble: 33,007 × 5
##    sexgp   AGE sleep1in sleep2in sleep3in
##    <chr> <dbl>    <dbl>    <dbl>    <dbl>
##  1 Women    54        1        1        1
##  2 Women    64        4        4        4
##  3 Women    65        1        1        1
##  4 Men      61        2        2        2
##  5 Men      38        1        1        1
##  6 Women    47        1        2        1
##  7 Men      48        1        1        1
##  8 Men      23        2        2        2
##  9 Women    28        2        2        2
## 10 Women    80        2        2        2
## # … with 32,997 more rows

이러한 파일을 wide라고 부릅니다. 이제 우리는 TSEX별, AGE와 Sleep1의 값을 얻고자 합니다. 그러면 TSEX와 AGE, TSEX와 Sleep이 하나의 group이 됩니다. group을 반복한것과 같습니다. 즉 group의 기본은 -로 제외하고, 이와 상관되는 변수는 모두 포함하여 variables라는 이름으로 만들고, 값은 values라는 변수로 변환시키는 것입니다.

dat1 %>% select(sexgp, AGE, sleep1in, sleep2in, sleep3in) %>%
  pivot_longer(-c(sexgp), names_to ="variables", values_to = "values")
## # A tibble: 132,028 × 3
##    sexgp variables values
##    <chr> <chr>      <dbl>
##  1 Women AGE           54
##  2 Women sleep1in       1
##  3 Women sleep2in       1
##  4 Women sleep3in       1
##  5 Women AGE           64
##  6 Women sleep1in       4
##  7 Women sleep2in       4
##  8 Women sleep3in       4
##  9 Women AGE           65
## 10 Women sleep1in       1
## # … with 132,018 more rows

평균에 대해서만 먼저 해보겠습니다.

dat1%>% select(sexgp, AGE, sleep1, sleep2, sleep3) %>%
  pivot_longer(-c(sexgp), names_to ="variables", values_to = "values") %>%
  group_by(sexgp, variables) %>%
  summarise(avg = mean(values, na.rm =TRUE), 
            std = sd(values, na.rm =TRUE), 
            p50 = quantile(values, prob = 0.50, na.rm =TRUE),
            p25 = quantile(values, prob = 0.25, na.rm =TRUE), 
            p75 = quantile(values, prob = 0.75, na.rm =TRUE) 
            ) %>%
  mutate(mean_std  = sprintf("%.1f\u00b1%.1f", avg, std), 
         median_IQR= sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)) %>%
  select(sexgp, variables, mean_std)
## # A tibble: 8 × 3
## # Groups:   sexgp [2]
##   sexgp variables mean_std 
##   <chr> <chr>     <chr>    
## 1 Men   AGE       45.3±13.9
## 2 Men   sleep1    4.5±0.8  
## 3 Men   sleep2    4.4±0.8  
## 4 Men   sleep3    4.2±0.9  
## 5 Women AGE       47.7±14.8
## 6 Women sleep1    4.3±0.8  
## 7 Women sleep2    4.3±0.9  
## 8 Women sleep3    4.2±1.0

무언가 복잡하지만 쉬워졌죠? 이말에 고개를 끄덕이셨다면 여러분은 이제 R coding에 빠져 들고 있는 것입니다.

이제 이것을 사람이 보기 편한 상태인 wide로 바꾸겠습니다. pivot_wider를 사용하고 우리가 원하는 가로로 필요한 정보를 names로 하고 원하는 값을 values로 하겠습니다. 이것을 tab1으로 정의하겠습니다.

dat1 %>% select(sexgp, AGE, sleep1, sleep2, sleep3) %>%
  pivot_longer(-c(sexgp), names_to ="variables", values_to = "values") %>%
  group_by(sexgp, variables) %>%
  summarise(avg = mean(values, na.rm =TRUE), 
            std = sd(values, na.rm =TRUE)
            ) %>%
  mutate(mean_std  = sprintf("%.1f\u00b1%.1f", avg, std), 
         median_IQR= sprintf("%.0f (%.0f-%.0f)", p50, p25, p75)) %>%
  select(sexgp, variables, mean_std) %>%
  pivot_wider(names_from = sexgp, values_from = c(mean_std)) -> tab1
tab1 %>% htmlTable(caption = "Table 1. Basic Characteristics of study population")
Table 1. Basic Characteristics of study population
variables Men Women
1 AGE 45.3±13.9 47.7±14.8
2 sleep1 4.5±0.8 4.3±0.8
3 sleep2 4.4±0.8 4.3±0.9
4 sleep3 4.2±0.9 4.2±1.0

6.3 Count and Distribution

heal_prob1 이 요통 변수라고 했었는데요, 1이 요통이 있다. 2가 요통이 없다 입니다. 이것을 표로 만들어 보겠습니다. ### 하나씩 반복 하나씩 반복하는게 가장 기본입니다. count를 사용하겠습니다. 요통이 몇몇 있나요?

dat1 %>% count(heal_prob1)
## # A tibble: 2 × 2
##   heal_prob1     n
##    <dbl+lbl> <int>
## 1   1 [있다]  8930
## 2   2 [없다] 24077

분율, 비율, 유병률을 나타내면 어떻게 될까요. 전체 합을 sum으로 만들고 이때 비율을 만들어야 합니다.

dat1 %>% count(heal_prob1) %>%
  mutate(total = sum(n)) %>%
  mutate(prob = n/total)
## # A tibble: 2 × 4
##   heal_prob1     n total  prob
##    <dbl+lbl> <int> <int> <dbl>
## 1   1 [있다]  8930 33007 0.271
## 2   2 [없다] 24077 33007 0.729

남자만 해보겠습니다.

dat1 %>% 
  filter(TSEX==1) %>%
  count(heal_prob1) %>%
  mutate(prob= n/sum(n)) 
## # A tibble: 2 × 3
##   heal_prob1     n  prob
##    <dbl+lbl> <int> <dbl>
## 1   1 [있다]  3613 0.234
## 2   2 [없다] 11860 0.766

여자만 해보세요.

6.3.1 Group_by summary and Table 1

group_by를 이용해 보겠습니다. 이면 위에서 연령을 이용한 방법을 실습했는데, 그것과 동일합니다

dat1 %>%
  group_by(TSEX) %>%
  count(heal_prob1) %>%
  mutate(prob = n/sum(n))
## # A tibble: 4 × 4
## # Groups:   TSEX [2]
##        TSEX heal_prob1     n  prob
##   <dbl+lbl>  <dbl+lbl> <int> <dbl>
## 1  1 [남성]   1 [있다]  3613 0.234
## 2  1 [남성]   2 [없다] 11860 0.766
## 3  2 [여성]   1 [있다]  5317 0.303
## 4  2 [여성]   2 [없다] 12217 0.697

익숙한 형태로 표시해 보겠습니다. 그리고 필요한 것만 남겨 보겠습니다 .

dat1 %>%
  group_by(TSEX) %>%
  count(heal_prob1) %>%
  mutate(prob = n/sum(n)) %>%
  mutate(smry1 = sprintf("%s (%.1f%%)", n, prob*100)) %>%
  select(TSEX, heal_prob1, smry1)
## # A tibble: 4 × 3
## # Groups:   TSEX [2]
##        TSEX heal_prob1 smry1        
##   <dbl+lbl>  <dbl+lbl> <chr>        
## 1  1 [남성]   1 [있다] 3613 (23.4%) 
## 2  1 [남성]   2 [없다] 11860 (76.6%)
## 3  2 [여성]   1 [있다] 5317 (30.3%) 
## 4  2 [여성]   2 [없다] 12217 (69.7%)

이번에는 연령을 5세단위로 바꾸로 각 연령의 분포를 확인해 보겠습니다.

dat1 <-dat1 %>% 
  mutate(agegp = case_when(
    AGE <25 ~ "<25",
    AGE <30 ~ "<30", 
    AGE <35 ~ "<35", 
    AGE <40 ~ "<40", 
    AGE <45 ~ "<45", 
    AGE <50 ~ "<50",
    AGE <55 ~ "<55", 
    AGE <60 ~ "<60",
    TRUE ~ "\u226560" # 나머지는 모두 >65 (\u2265는 크거나 같다는 symbol)
  )) 

dat1 %>%
  group_by(TSEX) %>%
  count(agegp) %>%
  mutate(prob = n/sum(n)) %>%
  mutate(smry1 = sprintf("%s (%.1f%%)", n, prob*100)) %>%
  select(TSEX, agegp, smry1) %>%
  arrange(TSEX, agegp)
## # A tibble: 18 × 3
## # Groups:   TSEX [2]
##         TSEX agegp smry1       
##    <dbl+lbl> <chr> <chr>       
##  1  1 [남성] <25   718 (4.6%)  
##  2  1 [남성] <30   1489 (9.6%) 
##  3  1 [남성] <35   1677 (10.8%)
##  4  1 [남성] <40   1998 (12.9%)
##  5  1 [남성] <45   1988 (12.8%)
##  6  1 [남성] <50   1872 (12.1%)
##  7  1 [남성] <55   1634 (10.6%)
##  8  1 [남성] <60   1479 (9.6%) 
##  9  1 [남성] ≥60   2618 (16.9%)
## 10  2 [여성] <25   833 (4.8%)  
## 11  2 [여성] <30   1390 (7.9%) 
## 12  2 [여성] <35   1613 (9.2%) 
## 13  2 [여성] <40   1707 (9.7%) 
## 14  2 [여성] <45   1861 (10.6%)
## 15  2 [여성] <50   2163 (12.3%)
## 16  2 [여성] <55   2315 (13.2%)
## 17  2 [여성] <60   2105 (12.0%)
## 18  2 [여성] ≥60   3547 (20.2%)

6.3.2 Long File and Table 1

그럼 이번에는 요통과 연령집단을 동시에 바꾸어 보겠습니다.

dat1 %>%
  mutate(backpain= case_when(heal_prob1==1 ~ "pain", 
                             TRUE ~ "no-pain")) %>%
  select(sexgp, agegp, backpain) %>%
  pivot_longer(-c(sexgp), names_to ="variables", values_to = "values")  %>%
  group_by(sexgp, variables) %>%
  count(values) %>%
  mutate(prob = n/sum(n)) %>%
  mutate(smry1 = sprintf("%s (%.1f%%)", n, prob*100)) %>%
  select(-n, -prob) %>%
  pivot_wider(names_from = sexgp, values_from = smry1) -> tab2

html 테이블로 살펴 보겠습니다.

tab2 %>% htmlTable()
variables values Men Women
1 agegp <25 718 (4.6%) 833 (4.8%)
2 agegp <30 1489 (9.6%) 1390 (7.9%)
3 agegp <35 1677 (10.8%) 1613 (9.2%)
4 agegp <40 1998 (12.9%) 1707 (9.7%)
5 agegp <45 1988 (12.8%) 1861 (10.6%)
6 agegp <50 1872 (12.1%) 2163 (12.3%)
7 agegp <55 1634 (10.6%) 2315 (13.2%)
8 agegp <60 1479 (9.6%) 2105 (12.0%)
9 agegp ≥60 2618 (16.9%) 3547 (20.2%)
10 backpain no-pain 11860 (76.6%) 12217 (69.7%)
11 backpain pain 3613 (23.4%) 5317 (30.3%)

무슨 생각이 드시죠? tab1과 tab2를 합치면 좋겠다는 생각이 드시죠, tab1에는 values 라는 변수가 없습니다 .그래서 합치기 어렵습니다. values 라는 변수를 생성하고 합쳐 보겠습니다.

tab1 = tab1 %>% mutate(values = "") %>% select(variables, values, Men, Women) 

rbind(tab1, tab2) %>%
  htmlTable()
variables values Men Women
1 AGE 45.3±13.9 47.7±14.8
2 sleep1 4.5±0.8 4.3±0.8
3 sleep2 4.4±0.8 4.3±0.9
4 sleep3 4.2±0.9 4.2±1.0
5 agegp <25 718 (4.6%) 833 (4.8%)
6 agegp <30 1489 (9.6%) 1390 (7.9%)
7 agegp <35 1677 (10.8%) 1613 (9.2%)
8 agegp <40 1998 (12.9%) 1707 (9.7%)
9 agegp <45 1988 (12.8%) 1861 (10.6%)
10 agegp <50 1872 (12.1%) 2163 (12.3%)
11 agegp <55 1634 (10.6%) 2315 (13.2%)
12 agegp <60 1479 (9.6%) 2105 (12.0%)
13 agegp ≥60 2618 (16.9%) 3547 (20.2%)
14 backpain no-pain 11860 (76.6%) 12217 (69.7%)
15 backpain pain 3613 (23.4%) 5317 (30.3%)

6.4 정리 1

  • 원하는 데이터를 불러옴
  • 원하는 변수를 선정
    • 관심 변수 (종속, 독립)
  • 변수 값을 살펴 보기
    • 숫자 인지 아닌지
    • missing value 가 얼마 인지
    • 제외할 변수가 얼마인지
  • 대표값 생성
    • count 함수 사용, prob 변수 생성
  • 표 생성
    • 반복
    • group_by
    • long file

6.4.1 과제

원하는 집단 변수 (남녀, 고용관계 등)과 연속형변수 (연령 등) 5개, 비연속변수 4개 를 선정하여 Table 1을 만드시오.