📜 제목으로 보기

✏마지막 댓글로

  • R 고급 실습(cohort study with dplyr-pipeline, broom-tidy) 포함됨

학습개요

학습목표

  • Cohort study의 정의를 설명할 수 있다.

  • Cohort study에서 RR과 OR을 추정하고, 위험요소와 질병의 독립성을 검정할 수 있다.

  • 로지스틱 회귀분석을 이용하여 confounding과 interaction에 대한 분석을 수행할 수 있다.

  • Cross-sectional study의 정의를 설명하고, cohort study와의 차이점을 설명할 수 있다.

주요 용어

  • Cohort study : 연구자가 위험요소에 노출된 그룹과 노출되지 않은 그룹을 일정기간 동안 추적관찰하여, 두 그룹 간 질병 발생률(incidence)이나 질병에 의한 사망률 등을 비교하는 연구

  • Mantel-Haenszel 방법 : Confounder값 별로 나눈 stratum에서 각각 구한 OR(또는 RR)을, 그 OR의 분산의 역수를 weight로 준 weighted average(가중평균)으로 종합하는 방법

  • Logistic regression : 결과변수(반응변수, 종속변수)가 이분형인 회귀분석

  • Cross-sectional study : 위험요소 노출여부와 질병유무를 한 시점에 동시에 조사하는 연구

Cohort study

연구디자인과 cohort study

image-20220417132324635

  • 연구종류에는
    • 대상자에게 뭔가를 가하는 Experimental(익스페리멘탈) study과
    • 대상자에게 뭔가를 하지않고 관찰만하는 Observational(옵져베이셔널) study가 있다.
  • 대상자에게 뭔가를 가하지 않는 옵져베이셔널 스터디(관찰연구)에는
    • 전향적인 연구인 Prospective study와 후향적인 연구인 Retrospective study가 있다.
      • 전향적인 연구: 연구시점부터 무슨일이 일어나는지 관찰한다.
      • 후향적인 연구: 다 끝난 뒤에서 무슨일이 일어났었는지 관찰한다.
    • 전향적인 관찰연구(→) 중 적어도 2시점 이상을 관찰하는 연구를 Cohort Study(=Longitudinal Study)라 한다.
    • 전향이든, 후향이든 관찰연구(→, ←) 중 딱 1시점에서만 관찰하는 연구를 Cross-sectional Study(=? Case-control study?)라 한다.

Cohort study

image-20220417133356197

  1. 연구자가 위험요소 1개에 대해 노출O(group)/ 노출X(group)의 그룹을 나누고

    • 위험요소 뿐만 아니라 치료법 O/X or 환자 특정 특성 O/X(ex나이)의 2그룹으로 나눌 수 있다.
      • ex> 남성/여성
    • 실제로는 2그룹 뿐만 아니라 3그룹 이상도 나눌 수 있다.
      • ex> 20-40세 / 40-60세 / 60세이상의 3그룹

  2. 일정기간 동안 추적 관찰하여, 그룹간 질병 발생률(incidence)나 질병 사망률(또는 다른 이벤트) 등을 비교하는 연구

cohort 수행 2가지 방법

방법1: 연구시작시점에서 (위험요소노출 O/X 등) 그룹을 먼저 나눈 뒤, 질병발생률만 관찰
  • 연구 시작 시점에 위험요소(치료법, 특성) 노출O 그룹 vs 노출X 그룹을 먼저 나누고 -> 각 그룹별 질병 발생/사망률 추적 관찰
    • 예> 방사능 노출이라는 위험요소에 대해, 노출시 아웃컴이 어떻게 변하는가?
      1. 원전사고가 있는 지역 = exposed group과 원전사고가 없지만, 비슷한 지역 = not exposed group 그룹을 먼저 나누고
      2. 질병 발생률이나 사망률을 follow up

image-20220417143658209

방법2: 연구시작시점에서 위험요소노출 X 집단을 선택하여, 노출여부 -> 빌병발생률 2가지를 순서대로 관찰
  • 연구 시작 시점에 아직 위험요소(치료법, 특성)에 노출X 집단만 선택해놓고 -> 위험요소 노출여부O/X 관찰 -> 질병발생률 다시 관찰(2가지를 관찰)
    • 이 때, 노출여부O/X는 랜덤화 하지않고, 자연스럽게 그룹이 나뉘게 한다.

image-20220417144449374

Cohort study의 장점

  • 연구 시작시점(방법1) or 연구기간(방법2) 중 위험요소 노출여부가 먼저 관측 -> 그 이후에 새로 발생한 질병 관측하기 때문에
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 질병발생 시점이 명확 -> incidence(질병발생 속도와 관련) 추정이 가능
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 인과관계 추론에 결정적 단서
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 위험요소 노출O/X이후 나올 수 있는 여러 질병들을 동시에 연구 가능하다.(타겟질병말고도 다른 질병이 나타나면 연구가능)

image-20220417145235716

Cohort study의 단점

  • 질병 종류에 따라서 너무 오래걸리거나 비용이 많이 들 수 있다.
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 질병 발생할때까지 계속 기다려야한다.
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 희귀병은 기다려도 몇 명 발생안하므로, 많은 수의 연구대상자가 확보되어야하고, 계속 follow up해야한다
      • 대안: retrospective cohort study
        • 전통적인 코호트 스터디는 전향적 연구( 연구시작 -> 노출O/X그룹 -> 질병발생률 관찰)인데
        • 요즘은 레트로스펙티브 코호트 스터디도 소개되고 있다.
        • 연구시작 시점 이후 발병하는 질병을 follow 하는 것은 동일하나,
        • 연구시작 시점 이전 위험요소 노출O/X를 <기록(데이터),기억(인터뷰)>를 통해 결정하는 것
          • 질병발생률은 prospective / 위험요소 노출여부는 retrospective로서 2가지가 섞여있는 디자인이라고도 볼 수 있다.
        • 위험요소 노출O/X를 기다렸다가 나누었다가 할 필요가 없이 바로 질병발생률만 관찰하므로 전체적인 연구기간이 짧아진다.
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 질병 발생할때까지 계속 기다리다가 위험요소 노출여부가 중간에 바뀌어버릴 수 있다.
      • 예> 흡연으로 위험요소 노출O였다가, 연구 중간에 금연해서 not/un exposed group으로 가버릴 수 있다
      • 대안: 위험요소 노출여부를 여러번 측정 -> longitudinal study(여러번 측정시)

image-20220417145655473

코호트 대표 예제

  • 전형적이지 않게, 1948년 연구 시작으로 지금까지 연구중
  • 위험요소가 무엇인지 모를 때 -> 심혈관계 질환 발생에 어떤 것들이??
    • 위험요소라는 결정적인 증거가 없었을 때, 위험요소들을 모아보기 위해 연구 시작
    • 식생활, 비만 등
  • 맨 첨 시작할 때는 위험요소노출O/X 이후 -> 질병발생까지 기다리는 시간이 오래 걸리니 30세 이상부터 연구 -> 5000명 정도 시작
    • 자식세대인 offspring cohort는 5000명 정도 follow 중
    • 3세대는 80세 미만으로 400명 정도가 follow 중 image-20220417152222107

코호트의 목적: 위험요소 - 질병발생 상관관계 탐구

  • 이 상관관계를 탐구하기 위해 사용하는 메져가 무엇이 있었는가?

    1. RR 랠러티브 리스크 for 상관관계
    2. OR 오즈 레이시오 for 상관관계
    3. attributable risk for 상관관계

  • RR과 OR을 코호트 스터디에선 어떻게 측정하는지 알아보자

RR & OR in cohort study

Review : RR & OR

image-20220417152833247

  • RR을 구하려면
    1. exposed risk(위험요소 노출 O group에서의 질병발생 확률)을 p1로 구한다.
    2. unexposed risk(위험요소 노출 X group에서의 질병발생 확률)을 p0로 구한다.
    3. RR는 p1/p0로 위험요소 노출X 그룹의 risk 분에 노출O 그룹의 risk로 위험요소노출 그룹별로 구한 risk들의 비율이다.
  • OR을 구하려면
    1. exposed risk p1을 구한다
    2. unexposed risk p0을 구한다
    3. exposed Odds를 p1/1-p1으로 구한다.
    4. unexposed Odds를 p0/1-p0으로 구한다.
    5. OR는 위험요소 노출X 그룹의 odds 분에 노출O 그룹의 odds로 위험요소노출 그룹별로 구한 odds의 비율이다.
RR in cohort study

image-20220417154258531

  1. 표에서 row/col을 보고 질병발생여부보다 위험요소 노출 여부가 있는 row별로 확인한다.

    • 위험요소 노출O -> 질병발생 a명 / 질병발생안한 사람 b명
    • 위험요소 노출X -> 질병발생 c명 / 질병발생안한 사람 d명

  2. RR은 각 위험요소 노출O/X group별 risk의 비율이므로 각 group별 risk를 구한다.

    1. exposed risk(p1) : a/a+b
    2. unexposed risk(p0) : c/c+d
    3. rr : p1/p0 = a/a+b / c/c+d image-20220417155327297
OR in cohort study

image-20220417155530944

  1. 표에서 row/col을 보고 질병발생여부보다 위험요소 노출 여부가 있는 row별로 확인한다. 2.OR는 각 위험요소 노출O/X group별 odds의 비율이므로 각 group별 odds를 구해야한다.

    1. 각 위험요소 노출여부 group별 odds는 그룹별 risk(p)부터 구한뒤 -> 그룹별 odds(p/1-p)를 구해야한다.
      1. exposed odds
        1. exposed risk(p1) : a/a+b
        2. exposed odds(p1/1-p1) : a/a+b / (1 - a/a+b) = a/b
      2. unexposed odds
        1. unexposed risk(p0) : c/c+d
        2. unexposed odds(p0/1-p0) : c/c+d / (1 - c/c+d) = c/d
    2. OR = exposed odds / unexposed odds = a/b / c/d = ad/bc image-20220417155613768

  2. 구하고 보니 위험요소여부-질병발생여부 표에서 cross-product하면, odds ratio를 바로 구할 수 도 있다.

    • 이 때, a는 위험요소노출o->질병o, d는 위험요소노출x->질병x 이다. image-20220417160310053 image-20220417160344988

예제로 보는 cohort study 속 위험요소-질병발생 상관관계 by RR, OR

상관관계의 크기를 보고 싶어 RR, OR의 값 추정(점 추정치)

image-20220417160458677

  1. 표를 보고 위험요소 노출O/X(row for group) -> 질병 발생여부(col) 순서로 확인한다.
    • 아~ 위험요소 2개의 interaction문제는 아니구나
    • 아~ 위험요소 1개 O/X별 -> 질병발생 O/X 의 관계이므로 -> 상관관계가 목적이고 -> RR, OR로 계산하겠구나. 그럴려면 위험요소 O/X를 group으로 나눠서 risk/odds를 계산하겠구나

  1. RR를 각 위험요소 노출여부 그룹별(→row별) risk 비율로 구한다.

    • RR = exposed risk(p1) / unexposed risk(p0) = (9/9+291) / (10/10+490) = 9/300 / 10/500 = 1.5

  2. OR는 위험요소노출여부-질병여부 표에서 ad/bc의 cross-product로 빠르게 구한다

    • 원래는 각 위험요소 노출여부 그룹별(→row별) odds(risk/1-risk) 비율로 구해야한다.
    • OR = a(oo)d(xx) / bc = 9 x 490 / 291 x 10 = 1.52 image-20220417161612247

상관관계 크기를 점 추정한 RR, OR에 신뢰구간 더해주기

  • RR, OR의 정확한 값 추정(점 추정치)는 95%신뢰구간으로 불확실성을 표시해줘야한다.
  • 신뢰구간의 여러계산법이 있지만, 거의다 비슷하게 나온다.
  • 다른 소프트웨어를 쓴다고 다른 값으로 나온다고 해도 의심안해도 된다.
  • 가장 많이 쓰는 방법은 Wald method(왈드 방법)인데, 좋은 방법이 아니라 계산이 쉬워서 흔하게 사용되는 것 image-20220417161744879

위험요소 노출-질병과의 상관관계가 통계적으로 있는 = 2개가 독립이 아닌지 통계적 검정까지 by RR,OR

  • RR, OR을 계산하는 이유는 위험요소 - 질병의 상관관계의 크기를 보기 위함이고,
  • RR, OR의 점추정에 대한 불확실성 때문에 상관관계의 크기 신뢰구간까지 더해준다.
  • 위험요소 - 질병의 상관관계가 통계적으로 유의한지를 2개가 독립이 아닌지로 통계적 검정까지 해줘야한다
    • 위험요소 - 질병 상관관계는 RR, OR로 구한다며... 근데.. 신뢰구간과 독립아닌지 검정까지 해줘야해?

  • 상관관계의 크기(RR, OR) 및 신뢰구간을 구했다면, 그래서 이것이 유의한 상관관계이냐를 독립이 아니냐의 통계적 검정(하이퍼 시시스 테스트)을 통해 증명한다.

    • 일반: 위험요소 - 질병간 상관관계가 있느냐?
    • 통계적 표현: 위험요소노출여부라는 변수와 질병발생여부라는 변수가 서로 독립이냐?(독립아니다의 대립가설을 채택해야함)

  • 위험요소 노출여부는 2그룹 이외에 3~4그룹으로도 나눌수 있지만, 질병발생여부는 무조건 2분형이어야한다.

  • 위험요소여부, 질병여부 -> 2개의 범주형 변수에 대한 독립성 검정은 카이스퀘어 검정을 이용한다

    • H0(귀무가설): 독립이다. RR=1 or OR=1
    • H1(대립가설): 독립이 아니다(상관관계 있다) RR != 1 or OR != 1
    • 카시스퀘어 검정은 샘플작을 때 작동을 잘 안하므로 Fisher's exact Test를 이용한다 image-20220417161854119
독립성검정(자세히 설명X)
  • 서로 독립이다를 수식으로 나타내면 Pij = Pi x Pj로 나타낸다.
    • 첫째 범주 i, 2번째 범주 j image-20220417163526215
  • 테스트하는 검정통계량
    • 귀무가설 하에서의 기대도수를 이용하여 계산된 검정통계량을 사용한다(공식) image-20220417163452119

R 실습(상관관계를 위한 RR,OR,신뢰구간,독립성검정까지 반영한 상관관계)

  • 여러 방법이 있지만, epiR 패키지로 계산한다.

    • 업데이트에 민감하니, epiR패키지를 업데이를 하고 R을 최신버전으로 업그레이드까지.. 안되면 epiR 제거후 다시 설치

  • 2by2 table을 입력할 때, 이미 형식이 정해져있다.

    • 위험요소가 row / 질병여부가 col에 있다고 가정하고
      • 위험o : 질병o(a) 질병x(b)
      • 위험x : 질병o(c) 질병x(d)
    • 의 순서대로 c(a, b, c, d)순서대로 입력해줘야한다.
      • my) 앞이 위험요소/뒤가 질병여부라 치면 row1: oo(a) ox(b) //row2: xo(c) xx(d)순으로 입력 image-20220417164351082 
# 아나콘다에서 epiR 검색후 r-epir 설치

library(epiR)

Loading required package: survival
Package epiR 1.0-4 is loaded
Type help(epi.about) for summary information

image-20220417164901121

#     9, 291, 
#     10, 490
# )

# tab <- c(9, 291, 10, 490)

# 바로 벡터만 넣어주면 안되더라... matrix nrow=2 byrow 써서 matrix를 만들어줘야함
# ?epi.2by2

tab <- matrix(c(9, 291, 10, 490), nrow = 2, byrow = TRUE)
# rownames(dat) <- c("Expose+", "Expose-"); 
# colnames(dat) <- c("Disease+", "Disease-");

# epi.2by2(dat = as.table(dat), method = "cross.sectional", 
#    conf.level = 0.95, units = 100, outcome = "as.columns")

tab
9 291
10 490 
epi.2by2(dat=tab, method="cohort.count") # method 코호트.카운트가 default긴 하다.

# case-control인 경우, method옵션을 바꾸면 된다.

             Outcome +    Outcome -      Total        Inc risk *        Odds
Exposed +            9          291        300              3.00      0.0309
Exposed -           10          490        500              2.00      0.0204
Total               19          781        800              2.38      0.0243

Point estimates and 95% CIs:
-------------------------------------------------------------------
Inc risk ratio                               1.50 (0.62, 3.65)
Odds ratio                                   1.52 (0.61, 3.77)
Attrib risk *                                1.00 (-1.29, 3.29)
Attrib risk in population *                  0.38 (-1.24, 1.99)
Attrib fraction in exposed (%)               33.33 (-62.19, 72.60)
Attrib fraction in population (%)            15.79 (-28.34, 44.75)
-------------------------------------------------------------------
 Test that odds ratio = 1: chi2(1) = 0.809 Pr>chi2 = 0.369
 Wald confidence limits
 CI: confidence interval
 * Outcomes per 100 population units 
# -> RR ( 95% 신뢰구간이 자동 제공된다.)

# Odds ratio                                   1.52 (0.61, 3.77)
# -> OR ( 95% 신뢰구간)

# 신뢰구간은 wald방법으로 뽑아진 것


# 아래 attributable risk로 계산 한 것
# Attrib risk *                                1.00 (-1.29, 3.29)
# Attrib risk in population *                  0.38 (-1.24, 1.99)
# Attrib fraction in exposed (%)               33.33 (-62.19, 72.60)
# Attrib fraction in population (%)            15.79 (-28.34, 44.75)

# Test that odds ratio = 1: chi2(1) = 0.809 Pr>chi2 = 0.369
# 0.05보다 높아서 -> H0 독립성 보장 -> **상관관계 유의하지 않다**

my) 카이제곱검정 -> 독립성검정 -> 상관관계로 해석해도 됨 h0:(독립성)상관관계 유의하지 않다 h1:(독립성X)상관관계 유의하다

confounding까지 고려한 코호트 속 위험요소-질병간 상관관계

  • 데이터분석 단계에서 confouning해결책 3가지 중에 2가지만 제시(후속강의에서 3번째 방법 진행)
  1. Stratification : 각 confounder 값 별로 strata를 나누어, 각 stratum별로 상관관계를 따로 추정한다
    • 각 confouder값 별로 따로 추정한 상관관계를 종합하는 방법이 Mantel-Haenszel 방법
      • 각 stratum별 RR이나 OR을 종합

  2. Adjustment : 기본적으로 Multivariable logistic regression를 하는 것으로 confounder를 공변량으로 넣은 회귀분석을 통해, confounder 값이 고정되어 있을 때의 X와 Y의 상관관계를 추정한다

    • Multivariable logistic regression: adjusted OR

      image-20220417170848737

Mantel-Haenszel 방법

image-20220417171320215

  1. confounder 값별로 따로 2by2 table을 만들고 -> OR(RR)을 계산한다
  2. 각 confounder값별로 나온 OR1,OR2 .. 를 평균을 내는데, 자신이 있는 지에 따라 분산의 역수를 weight로 준 가중 평균을 구한다.
    • 좀 더 계산에 자신 있는 OR에 비중을 많이 둬서 평균을 구하는데
    • (추정값에 대한) 자신 있다 없다는 분산(불확실성)으로 표현된다. 분산이 크면 -> 불확실성 크다 -> 자신 없다 -> 분산의 역수 = 자신 있다 -> 분산의 역수를 weight로 줘서 자신있는 곳에 비중을 둔다.
  • 예
    1. confoudner가 성별이라고 치면, 여성에 대한 2by2 table , 남성에 대한 2by2 table을 따로 만들고, OR도 따로 구한다.
    2. OR의 평균을 분산의 역수로 가중평균을 구한다. -> 추정값에 대한 변동이 적은=추정값에 대해 불확실성이 적은 = 추정값에 대해 내가 더 자신있는 = 대체로 샘플사이즈가 더 큰 confounder의 OR을 더 많이 반영할 수 있다. image-20220417172039782
  • 계산과정을 생략하고, 전체적으로 보면
    • confounder 값의 갯수만큼(I개)로 2by2 table -> OR이 나오며, 계산방법을 따라서 OR를 종합할 수 있다. image-20220417172352912

R 실습(confounding or interaction까지 고려한 상관관계) 

# install.packages("Epi")

also installing the dependencies 'RcppArmadillo', 'cmprsk', 'etm'


  There are binary versions available but the source versions are later:
                  binary     source needs_compilation
RcppArmadillo 0.10.4.0.0 0.11.0.0.0              TRUE
cmprsk            2.2-10     2.2-11              TRUE
Epi                 2.44       2.46              TRUE

  Binaries will be installed
package 'RcppArmadillo' successfully unpacked and MD5 sums checked
package 'cmprsk' successfully unpacked and MD5 sums checked
package 'etm' successfully unpacked and MD5 sums checked
package 'Epi' successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	C:\Users\cho_desktop\AppData\Local\Temp\RtmpgfxhqH\downloaded_packages

library(epiR)
library(Epi) # 데이터용
library(dplyr) # 데이터정리용

Warning message:
"package 'dplyr' was built under R version 3.6.3"
Attaching package: 'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union


?nickel # lung cancer (ICDs 162 and 163), nasal cancer (ICD 160)

data(nickel) # 예제데이터

dat1 <- nickel
head(dat1)

# exposure: 위험요소 노출된 정도 -> 위험요소 그룹은 2개를 넘어선다
# age1st : 위험요소에 노출된 첫번째 나이
# icd : 어떤이유로 사망했는지 알 수 있는데, 죽지 않았으면 0 
# - lung cancer는 162, 163으로 코딩되어있고 / nasal cancer 는 160으로 코딩
id icd exposure dob age1st agein ageout
3 0 5 1889.019 17.4808 45.2273 92.9808
4 162 5 1885.978 23.1864 48.2684 63.2712
6 163 10 1881.255 25.2452 52.9917 54.1644
8 527 9 1886.340 24.7206 47.9067 69.6794
9 150 0 1879.500 29.9575 54.7465 76.8442
10 163 2 1889.915 21.2877 44.3314 62.5413

dplyr

매핑(ifelse)한 칼럼추가(mutate) 
# ifelse()로 삼항연산자를 이용한 매핑을 한다.

dat1 <- nickel %>% 
    mutate(lung.cancer = ifelse(icd==162 | icd == 163, 1, 0))

head(dat1)
id icd exposure dob age1st agein ageout lung.cancer
3 0 5 1889.019 17.4808 45.2273 92.9808 0
4 162 5 1885.978 23.1864 48.2684 63.2712 1
6 163 10 1881.255 25.2452 52.9917 54.1644 1
8 527 9 1886.340 24.7206 47.9067 69.6794 0
9 150 0 1879.500 29.9575 54.7465 76.8442 0
10 163 2 1889.915 21.2877 44.3314 62.5413 1 
dat1 <- nickel %>% mutate(
    lung.cancer = ifelse(icd==162 | icd == 163, 1, 0),
    nasal.cancer = ifelse(icd==160, 1, 0),
)

head(dat1)
id icd exposure dob age1st agein ageout lung.cancer nasal.cancer
3 0 5 1889.019 17.4808 45.2273 92.9808 0 0
4 162 5 1885.978 23.1864 48.2684 63.2712 1 0
6 163 10 1881.255 25.2452 52.9917 54.1644 1 0
8 527 9 1886.340 24.7206 47.9067 69.6794 0 0
9 150 0 1879.500 29.9575 54.7465 76.8442 0 0
10 163 2 1889.915 21.2877 44.3314 62.5413 1 0
위험요소 변수를 연속형 -> 특정값을 기준으로 0or1의 이분형으로 매핑(ifelse) 칼럼 추가(mutate) 
dat1 <- nickel %>% mutate(
    lung.cancer = ifelse(icd==162 | icd == 163, 1, 0),
    nasal.cancer = ifelse(icd==160, 1, 0),
    # 연속형 변수로 표현되어있던 위험요소 노출 정도를 -> 특정값 기준으로 위험요소 노출여부(1or0)로 매핑
    exposure.grp = ifelse(exposure >=  10, 1, 0),
    # 위험요소노출 나이 또한 2분형 범주칼럼으로 매핑
    age1st.grp = ifelse(age1st >=  25, 1, 0))

head(dat1)
id icd exposure dob age1st agein ageout lung.cancer nasal.cancer exposure.grp age1st.grp
3 0 5 1889.019 17.4808 45.2273 92.9808 0 0 0 0
4 162 5 1885.978 23.1864 48.2684 63.2712 1 0 0 0
6 163 10 1881.255 25.2452 52.9917 54.1644 1 0 1 1
8 527 9 1886.340 24.7206 47.9067 69.6794 0 0 0 0
9 150 0 1879.500 29.9575 54.7465 76.8442 0 0 0 1
10 163 2 1889.915 21.2877 44.3314 62.5413 1 0 0 0
boolean형 칼럼들 -> 일괄 factor형으로 type변환 by mutate_at()
  • 여러개의 작업에 대해 똑같은 작업을 할 때 mutate_at( vars( 변수들 나열), list(~ 개별로 적용할 작업(.))
    • 참고) %>% 직전메서드 바로 뒤에 달아주고, 다음 메서드는 다음줄에 가능하다
dat1 <- nickel %>% mutate(
    lung.cancer = ifelse(icd==162 | icd == 163, 1, 0),
    nasal.cancer = ifelse(icd==160, 1, 0),
    exposure.grp = ifelse(exposure >=  10, 1, 0),
    age1st.grp = ifelse(age1st >=  25, 1, 0)) %>%
                    mutate_at(vars(lung.cancer, nasal.cancer, 
                                   exposure.grp, age1st.grp),
                              list(~as.factor(.)))

head(dat1)
id icd exposure dob age1st agein ageout lung.cancer nasal.cancer exposure.grp age1st.grp
3 0 5 1889.019 17.4808 45.2273 92.9808 0 0 0 0
4 162 5 1885.978 23.1864 48.2684 63.2712 1 0 0 0
6 163 10 1881.255 25.2452 52.9917 54.1644 1 0 1 1
8 527 9 1886.340 24.7206 47.9067 69.6794 0 0 0 0
9 150 0 1879.500 29.9575 54.7465 76.8442 0 0 0 1
10 163 2 1889.915 21.2877 44.3314 62.5413 1 0 0 0
이제 만들어진 범주형변수들 중 필요한 위험요소1 - 질병1을 뽑아 2by2 table을 만든다. 
# 질병 여부 -> nasal.cancer

# 정해진 범주형 칼럼을 table()로 만든다. -> 알아서 범주종류별로 n by m을 만들어주나보다.
# -> 만약, 값만 있다면? as.matrix( a,b,c,d , row=2, byrow=true)

# 먼저 적어준 범주형칼럼이 row로 오게 된다?
table(dat1$exposure.grp, dat1$nasal.cancer)

   
      0   1
  0 595  41
  1  28  15
tab2 <- table(dat1$exposure.grp, dat1$nasal.cancer,
              dnn = c("exposure.grp", "nasal.cancer"))

tab2 

# 출력해보니, 위험요소/질병여부 oo -> ox 순이 아니라  xx -> xo 순으로 나온다.
# -> 코호트나 분석을 위해선, oo, ox xo xx 순으로 나타나도록 
# --> 데이터 dat1 속 factor의 level순서를 재설정해줘야한다.

            nasal.cancer
exposure.grp   0   1
           0 595  41
           1  28  15
범주형 0,1 factor의 레벨 순서를 바꿔주기 by relevel() + mutate() 칼럼추가 외 칼럼 덮어쓰기를 이용해서
  • relevel(, ref="1"): reference 그룹 = 첫번째 기준그룹을 직접 지정해준다. default로 0, 1중 작은 수인 "0"으로 지정되어있다.
dat2 <- dat1 %>% mutate(lung.cancer=relevel(lung.cancer, ref="1"), 
                        nasal.cancer=relevel(nasal.cancer, ref="1"), 
                        exposure.grp=relevel(exposure.grp, ref="1"), 
                        )
#                         age1st.grp=relevel(age1st.grp, ref="1")) 
                        # 이놈은 나중에 confounder로 사용될 변수라서.. 딱히 1의 그룹을 첫번째그룹으로 안주어도 된다.

head(dat2)
id icd exposure dob age1st agein ageout lung.cancer nasal.cancer exposure.grp age1st.grp
3 0 5 1889.019 17.4808 45.2273 92.9808 0 0 0 0
4 162 5 1885.978 23.1864 48.2684 63.2712 1 0 0 0
6 163 10 1881.255 25.2452 52.9917 54.1644 1 0 1 1
8 527 9 1886.340 24.7206 47.9067 69.6794 0 0 0 0
9 150 0 1879.500 29.9575 54.7465 76.8442 0 0 0 1
10 163 2 1889.915 21.2877 44.3314 62.5413 1 0 0 0 
tab2 <- table(dat2$exposure.grp, dat2$nasal.cancer,
              dnn = c("exposure.grp", "nasal.cancer"))

tab2

            nasal.cancer
exposure.grp   1   0
           1  15  28
           0  41 595
epi.2by2(tab2)

             Outcome +    Outcome -      Total        Inc risk *        Odds
Exposed +           15           28         43             34.88      0.5357
Exposed -           41          595        636              6.45      0.0689
Total               56          623        679              8.25      0.0899

Point estimates and 95% CIs:
-------------------------------------------------------------------
Inc risk ratio                               5.41 (3.27, 8.96)
Odds ratio                                   7.77 (3.85, 15.69)
Attrib risk *                                28.44 (14.06, 42.81)
Attrib risk in population *                  1.80 (-1.01, 4.62)
Attrib fraction in exposed (%)               81.52 (69.40, 88.84)
Attrib fraction in population (%)            21.84 (9.26, 32.67)
-------------------------------------------------------------------
 Test that odds ratio = 1: chi2(1) = 43.042 Pr>chi2 = < 0.001
 Wald confidence limits
 CI: confidence interval
 * Outcomes per 100 population units

confounder를 고려 및 보정(핸들)해주며 코호트 속 상관관계 실습

핸들1: Stratification(confounder를 선택후 값별(strrata)로 나눠서 2by2 테이블을 보는 층화)
  • age1st.grp도 0or1로 매핑했었다 -> 이걸 confounder라고 가정하고, confounder값별로 나눠서 2by2 테이블을 그리는 것을 Stratification라 한다
위험요소1 + 질병여부1 + confounder변수까지 추출해서 table을 만든다. 
tab3 <- table(dat2$exposure.grp, dat2$nasal.cancer, dat2$age1st.grp,
              dnn = c("exposure.grp", "nasal.cancer", "age1st.grp"))

tab3 
# table()에 3번째 인자를 주면, 그 값별로 1,2번째 칼럼으로 만든 table을 나눠서 보여준다.
# 3번째칼럼의 범주종류별로 3차원을 table이 완성된다고 생각하자.

, , age1st.grp = 0

            nasal.cancer
exposure.grp   1   0
           1   1   9
           0  12 269

, , age1st.grp = 1

            nasal.cancer
exposure.grp   1   0
           1  14  19
           0  29 326

epi.2by2( tab3 )

# 데이터 자체는 위에서 만든 table과 유사하나, 아래 요약글을 보면
# RR, OR이 여러줄로 나타난다.( cf) rsk ratio = relative risk


# Inc risk ratio (crude)   -> crude: confounder(age1st.grp)을 무시한 rr + 신뢰구간

# Inc risk ratio (M-H) -> [맨틀-핸젤 방법으로 confounder별 (다른 table->) 다른 RR을 가중평균한] RR + 신뢰구간

# Odds ratio (crude)                           7.77 (3.85, 15.69)
# Odds ratio (M-H)                             7.08 (3.45, 14.54)
# -> confounder 무시한 OR   +   cofounder별 OR구한 뒤 맨틀-핸젤방법으로 가중평균한 OR

             Outcome +    Outcome -      Total        Inc risk *        Odds
Exposed +           15           28         43             34.88      0.5357
Exposed -           41          595        636              6.45      0.0689
Total               56          623        679              8.25      0.0899


Point estimates and 95% CIs:
-------------------------------------------------------------------
Inc risk ratio (crude)                       5.41 (3.27, 8.96)
Inc risk ratio (M-H)                         4.78 (2.87, 7.97)
Inc risk ratio (crude:M-H)                   1.13
Odds ratio (crude)                           7.77 (3.85, 15.69)
Odds ratio (M-H)                             7.08 (3.45, 14.54)
Odds ratio (crude:M-H)                       1.10
Attrib risk (crude) *                        28.44 (14.06, 42.81)
Attrib risk (M-H) *                          27.34 (-32.35, 87.03)
Attrib risk (crude:M-H)                      1.04
-------------------------------------------------------------------
 M-H test of homogeneity of incidence risk ratios: chi2(1) = 0.603 Pr>chi2 = 0.438
 M-H test of homogeneity of odds ratios: chi2(1) = 0.769 Pr>chi2 = 0.38
 Test that M-H adjusted odds ratio = 1: chi2(1) = 34.954 Pr>chi2 = <0.001
 Wald confidence limits
 M-H: Mantel-Haenszel; CI: confidence interval
 * Outcomes per 100 population units 
# [m]easure of [assoc] iation의 detail이라는 뜻.
epi.2by2( tab3 )$massoc.detail$OR.strata.wald
epi.2by2( tab3 )$massoc.detail$RR.strata.wald

NULL
NULL

  • 나는 안되네.. image-20220417214555170

    image-20220417214605495

핸들2: adjustment
  • confounder를 어저스트먼트(adjustment)로 보정하기 위해선 Multivaraible logistic regression이 필요하므로 logisitic regression에 대해 먼저 알아보자
실습전 코호트-confounder 핸들을 위한 Logistic regression
  • 결과 변수가 이분형(binary)인 회귀분석
  • 우리는 지금 코호트 속 설명변수 위험요소 노출여부에 대한 질병여부로서 결과변수를 이분형을 사용중이니 이분형 결과변수를 사용하는 logisitic regression을 적절한 회귀분석으로 사용할 수 있다
    • 위험요소 노출여부를 x(설명변수) , 질병여부를 y결과변수(종속변수)로 로지스틱 회귀분석 모델을 잡을 수 있는데
    • 회귀식의 좌항은 ln(질병발생O 그룹의 odds) = 1차 회귀식(독립변수의 linear combination)으로 성립한다.
      • 즉, 로지스틱 회귀 모델은 질병발생O의 odds에 자연로그 취한 값을 = 독립변수들의 linear combination으로 설명하겠다는 뜻이 된다.

image-20220417220454071

image-20220417215046958

Logistic regression에서 가장 중요한 OR(odds rations)와 B(회귀계수)의 관계
  • ln(질병O의 odds) = b0 + b1*X의 로지스틱 회귀 식에서
    1. X=1을 대입 후, 양변에 exponentiate(e, 지수함수)을 취해주자.
      • X=1이면, 위험요소 노출O의 상황이므로 risk를 p1이라고 표현한다(앞에서는 위험여부O일 때의 risk를 p1로 취급함) image-20220417221216636
    2. X=0을 대입후, 위험요소 노출O의 상황이므로 risk를 p0이라 표현하고, 양변에 지수함수를 취한다. image-20220417221220788
    3. OR(odds ratio)는 X=0일 때의 odds 분에 X=1일때의 odds (odds의 비율)이므로, 우항을 나눠보면 e^b1이 된다. image-20220417221511346
    4. 결과적으로 위험요소 노출여부 변수(X)의 회귀계수B1에다가 e^x지수함수를 취해주면 -> 질병에 대한 위험요소의 OR(odds ratio)다
암기 -> 회귀계수에다가 지수함를 취하면 [질병여부에 대한 해당 위험요소의 Odds ratio]다
multivariable logistic regression
  • 설명변수가 1-> d개로 늘어났다.
  • 설변변수 1개일 때와 똑같은 식에다가 각 설명변수 + 곱해질 회귀변수를 더해주기만 하면 된다. image-20220417223259059

-“Adjusted” logistic regression model이라고 부르기도 한다

  • 설명변수 1개일 때의 OR <-> b간의 관계가 똑같다.
    • 즉, Xj의 회귀계수에 지수함수 -> 결과변수에 대한 j번째 위험요소의 OR (나머지 독립변수들은 고정) image-20220417223634130
  • 만일 X1만 위험요소 노출여부라고 보고, 나머지 X2~Xd를 confounder들이라고 가정하면
    • e^b1 = 위험요소의 OR (단, 나머지 confounder 값들이 고정되어있을 때!)
  • 주의점: sample size에 비해 너무 많은 독립변수들을 넣으면 안된다.
    • 오버피팅: 내가 가진 데이터의 정보량에 비해, 너무 복잡한 모델(변수가 많은 모델)을 fitting하게 됨으로서, 그 데이터의 signal뿐만 아니라 데이터가 가진 noise까지 fitting에 사용되어 문제가 발생한다.
    • 기준은 없지만, 일반적으로 1개 독립변수 추가시 마다 환자 수 10~20명(데이터)가 필요하다
      • 즉, 데이터 10개마다 1개 변수, 20개마다 1개 변수를 사용할 수 있다.
      • 예를 들어, 100명의 연구대상자를 가지고 연구한다
        • 10으로 나누면 10개까지 독립변수 가능
        • 20으로 나누면 5개까지 독립변수 가능

interaction를 고려

  • 위험요소 X가 있을 때, X의 효과( = 질병과의 상관관계)를 OR(or RR)로 젤 수 있는데
    • X의 OR이, 다른 변수 A에 의해 달라질 때, X와 A간의 인터렉션이 있다고 한다.
  • X와 A간 인터렉션이 통계적으로 유의한가?도 판단해야한다.
    • A값별 X의 y에 대한 OR이 통계적으로 다른가?를 판단해야한다.
    • 로지스틱 회귀분석의 옵션 중 interaction term을 넣어서 분석하면 된다.
      • interaction term: 두 변수 X와 A를 곱한 변수(X X A)
  • 왜 두 변수를 곱한 것이 interaction term이라는 로지스틱 회귀분석의 추가 변수가 되는가
    • 앞에서 석면과 흡연간에 인터렉션이 있는가?를 본적 있다. image-20220417231152292
my) confounder 값별로, 인터렉션 유발 변수 값별로 RR, OR를 보고 검정하고 하네..
  • X석면(X1)의 효과가, A흡연(X2)에 따라 다른가?의 인터렉션에 관심있다고 가정하고 보면
    1. 로지스틱 모델 회귀식을 쓰되, 인터렉션 텀인 b3(X1 * X2)도 추가해서 작성하자
      • ln(질병O 그룹의 오즈) = b0 + b1X1(석면노출여부에 대한 term) + b2X2(흡연노출여부에 대한 term) + b3X1X2(두 변수를 곱한 인터렉션 term) image-20220417231535449
    2. 이제 interaction을 발생시키는 A(흡연, X2)변수의 값별로 X의 OR을 봐야한다.
      1. X2 = 0 (흡연여부가 비흡연인 경우)일 때, X1의 OR(X1 앞에 달려있는 회귀계수에 지수함수 취하면 됨)을 보기 image-20220417232002440
      2. X2 = 1 (흡연여부가 흡연인 경우)일 때, X1의 OR을 보기
        • X1의 회귀계수가 인터랙션 term때문에 (b1 + b3)X1으로 바뀌게 되었다. image-20220417232005888
    3. 이제 A값별 구한 OR들을 비교해보자.
      • 비흡연자의 석면OR e^b1 vs 흡연자의 석면OR e^(b1+b3)
      • 만약, b3가 0이라면, 비흡연자와 흡연자의 석면OR이 동일해진다.
      • 즉, b3 = 0이라면 인터랙션이 없다라는 뜻이다.
    4. b3 = 9인지 아닌지 검정 -> 인터랙션 없는지 있는지 검정하는 것과 같다 -> 흡연여부에 따른 석면OR이 변화가 있는지 없는지 검정하는 것과 같다

R실습(adjustment 및 interaction 고려한 )

logistric regression을 이용한 confounding과 interaction 알아보기
위험요소노출여부(X) 질병여부(Y)만 고려하기 
head(dat1)
id icd exposure dob age1st agein ageout lung.cancer nasal.cancer exposure.grp age1st.grp
3 0 5 1889.019 17.4808 45.2273 92.9808 0 0 0 0
4 162 5 1885.978 23.1864 48.2684 63.2712 1 0 0 0
6 163 10 1881.255 25.2452 52.9917 54.1644 1 0 1 1
8 527 9 1886.340 24.7206 47.9067 69.6794 0 0 0 0
9 150 0 1879.500 29.9575 54.7465 76.8442 0 0 0 1
10 163 2 1889.915 21.2877 44.3314 62.5413 1 0 0 0 
# -> logistic regression으로 보기

# 로지스틱 회귀분석은 glm()메서드를 이용한다.
# - 칼럼만 입력하며, 자동완성 안되는 메서드임.
# - family에 회귀 종류 기입
# - 회기분석은 메서드 결과를 변수에 받아놓고, summary()해서 보는게 좋다.
obj1 <- glm(nasal.cancer ~ exposure.grp, data = dat1, family = "binomial")
summary(obj1)

Call:
glm(formula = nasal.cancer ~ exposure.grp, family = "binomial", 
    data = dat1)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9263  -0.3651  -0.3651  -0.3651   2.3416  

Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept)    -2.6750     0.1615 -16.567  < 2e-16 ***
exposure.grp1   2.0508     0.3584   5.722 1.05e-08 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 386.72  on 678  degrees of freedom
Residual deviance: 359.73  on 677  degrees of freedom
AIC: 363.73

Number of Fisher Scoring iterations: 5
회귀분석 결과 분석
  • 제일 중요한 것은 Coefficents:의 회귀계수(코 에피션트) 부분이다. 
    Coefficients:
                   Estimate Std. Error z value Pr(>|z|)    
b0 <- (Intercept)    -2.6750     0.1615 -16.567  < 2e-16 ***
# 상수항
b1 <-exposure.grp1   2.0508     0.3584   5.722 1.05e-08 ***
# b1의 추정치
  • glm()의 결과물 obj -> summary()가지고는 우리가 관심있는 measure of association인 e^b의 OR가 안보인다.
    • exp( 각 회귀계수 Eistimaated 추정치)를 넣으면 직접계산도 가능하다.
  • "broom"패키지dml tidy()를 이용하여 glm() -> obj -> tidy()로
    • 회귀계수(Coefficient)들 대신 e^(b)의 OR을 구할 수 있다.


  There is a binary version available but the source version is later:
      binary source needs_compilation
broom  0.7.6  0.8.0             FALSE


installing the source package 'broom'


library(broom)

obj1 <- glm(nasal.cancer ~ exposure.grp, data = dat1, family = "binomial")
# summary(obj1)

# tidy(obj1) # 특별한 옵션을 주지 않으면, 그냥 b의 추정치가 나온다.

# exponentiate = T 옵션을 주면, estimate(변수명)은 똑같지만, OR이 찍혀있게 된다.
# tidy(obj1, exponentiate = T)

# OR에 대한 신뢰구간(conf idential int  erval)옵션까지 줄 수 있다.
tidy(obj1, exponentiate = T, conf.int = T)
term estimate std.error statistic p.value conf.low conf.high
(Intercept) 0.06890756 0.1614641 -16.567079 1.205399e-61 0.04944691 0.09326859
exposure.grp1 7.77439024 0.3584015 5.722172 1.051706e-08 3.78406082 15.55773863
  • intercept의 b0의 OR은 볼 필요 없다
  • exposure.grp1(X)의 OR을 (e^회귀계수 없이) 바로 볼 수 있다.
    • OR -> estimate:7.77
    • p.value -> 1.051706e-08
    • conf.low ~ conf.high
위험요소노출여부(X)-병여부(Y) + @의심되는 confounder(age1st.grp)까지 고려해서 로지스틱 회귀 모델 
# confounder로 생각되는 변수(첫 노출 나이)를 + 공변량으로 넣어서
# -> OR이 아니라 adjustment OR로 계산되도록 한다.
obj2 <- glm(nasal.cancer ~ exposure.grp + age1st.grp
            , data = dat1, family = "binomial")

tidy(obj2, exponentiate = T, conf.int = T)

# exposure.grp1의 OR age1st.grp1dml OR이 각각 나온다.
# 해석: age1st.grp1그룹의 값이 일정할 경우 -> exposure.grp1의 OR이 6.86으로 추정된다
#       OR의 p-value가 작고 + 신뢰구간까지 제시
term estimate std.error statistic p.value conf.low conf.high
(Intercept) 0.04012872 0.2916948 -11.024065 2.925420e-28 0.02151766 0.06810862
exposure.grp1 6.86108670 0.3635830 5.296908 1.177802e-07 3.30603509 13.86633962
age1st.grp1 2.31360264 0.3340699 2.510869 1.204342e-02 1.23213568 4.61419281
위험요소노출여부(X)-병여부(Y) + @ age1st.grp가 confounder가 아니라 effective modifier로 의심되서 interaction term까지 고려해서 로지스틱 회귀 모델 
# R에서는  a + b + a*b 대신  -> a*b만 넣어주면, 각각의 main effect가 반영된 interaction term까지 고려 모델이 된다.
obj3 <- glm(nasal.cancer ~ exposure.grp * age1st.grp
            , data = dat1, family = "binomial")

tidy(obj3, exponentiate = T, conf.int = T)

# 해석: interaction term을 의미하는 [exposure.grp1:age1st.grp1]는 OR이 중요한게 아니라 p-value가 중요하다.
# -> 3.02 x 10^-1 = 0.3 으로 유의하지가 않다.
# --> inetraction term(exposure.grp1:age1st.grp1)의 p-value가 유의하지 않으면 -> interaction이 유의하지 않다로 해석한다.
# --> 첫노출나이 그룹에 따른, exposure이 바뀌는 것이 통계적으로 유의하진 않다.
term estimate std.error statistic p.value conf.low conf.high
(Intercept) 0.04460967 0.2949937 -10.5419359 5.534766e-26 0.02364332 0.07594055
exposure.grp1 2.49074058 1.0945923 0.8337169 4.044406e-01 0.12977341 14.95214891
age1st.grp1 1.99412053 0.3529469 1.9555438 5.051892e-02 1.02220683 4.13027193
exposure.grp1:age1st.grp1 3.32556576 1.1660801 1.0304951 3.027777e-01 0.45951598 68.82990394
interaction term이 유의하게 나온 경우?! 의심변수 0일때와 1일때의 OR을 각각 확인한다. 
# -> 만약, interaction이 유의한 경우라고 가정한다면???

# interaction을 일으키는 첫노출나이의  0일때의 OR  과   1일때의 OR을 확인해야한다.
# -> 각 변수들의 회귀계수를 구해야 -> 각 OR을 알 수 있다.

# summary(obj3) # 회귀계수를 포함한 전체 정보가 나온다.
coef(obj3)# 회귀계수만 벡터정보로 얻을 수 있다.
<dl class=dl-horizontal>
(Intercept)
-3.10980466691791
exposure.grp1
0.912580089581693
age1st.grp1
0.690203115537809
exposure.grp1:age1st.grp1
1.20163981224723
</dl> 
# -> 단일 변수의 회귀계수로 OR을 구하면,  interaction변수 0일 때의 OR임
exp(coef(obj3)[2])
exposure.grp1: 2.4907405840833 
# -> 단일 변수의 회귀계수  +   interaction term의 회귀계수의 합  에 exp
exp( coef(obj3)[2] + coef(obj3)[4] )
exposure.grp1: 8.2831215970951

연구에서 Confounding과 Interaction(effect modifier)를 만났을 때 설계

실제 연구 중 confounder와 effect modifier(interactio n유발변수)가 있다고 의심이 된다?

  1. 일단은 confounder와 effect modifier 의심변수(A)을 무시하고, 관심있는 위험노출변수의 OR(crude OR)만 계산한다.
  2. A변수의 값별로 따로 따로 OR (stratum-specific OR)을 계산한다.
    • A변수=0은 a*b로 만든 회귀 모델 중 자신의 회귀계수로 OR
    • A변수=1은 a*b로 만든 회귀 모델 중 (자신의 회귀계수 + axb의 회귀계수)로 OR을 구했다.
  3. 아래와 같이 경우를 나눠 판단한다.
    1. Stratum-specific OR들이 비슷하고 crude OR과도 비슷할 경우
      • confoundin이나 interactio이 없다
    2. Stratum-specific OR들은 비슷(interaction없다)하지만 crude OR과는 차이가 클 경우
      • A가 confounder이다. M-H OR(종합 OR)이나 adjusted OR(로지스틱 회귀를 이용)을 계산한다
    3. Stratum-specific OR이 많이 다르다
      • interaction이 있다. 로지스틱 회귀분석을 이용하여 interaction이 있는지 검정하고 A변수 값 별 OR을 보고한다

2분형 결과변수의 코호트 분석 제한점

2분형 결과변수(O/X)로는 제대로 코호트 분석이 안된다? 이상적이진 않다

이분형의 발병여부 -> 분석시점에만 고려하면, FU기간 기준이 없어서 통계적으로 invalid

실제로 코호트 연구에서 2분형 변수가 ideal하진 않다. 그 이유는 많은 경우가 코호트 스터디 환자마다 follow up기간이 다르기 때문

  • 동시에 시작해서 동시에 관찰이 끝나면 문제가 없지만
  • 많은 경우, 연구 시작 -> 환자 등록은 중간 중간 이루어져서 추가됨
    • 반면, 연구가 끝나는/분석하는 시점은 딱 한 시점

  • 예를 들어, 1년마다 환자 등록이 각각된, 시작은 2010년 코호트 연구

    • Q. 4번 환자가, 연구시점/분석시점(2015)엔 발병안했지만, 2016에 발병했다면? 어떻게 해석해야할까? image-20220423122306391

    • 2015 시점 분석시, 2,3번 발병 / 1,4번 발병안함 -> 2016년엔 ??

      • 4번이 2016년에 발병했다면?, 2번과 똑같이 3년만에 발병인데, 분석시점이후라서 카운팅 안해도 되는가?
        • 공평한 분석이 아닌 것 같다. image-20220423123032051

    • FU(follow up) duration의 차이를 무시하고, 딱 분석시점에 잘라서 해도 되는 걸까?(코호트는 추가되는 환자등록이 다르고, 마지막 자르는 시점은 같음)

      • 각 환자별 동일한 기준으로 질병발병여부를 젠 것이 아니므로, 통계적으로 valid하지 않은 것이 된다.

이분형의 발병여부 -> 분석시점 기준 FU기간만 제일 짧은(제일 최근 등록)환자 기준으로 바꾸면 통계적으로 valid? but FU기준을 넘어가는 정보들을 미반영

  • Q. 발병되서 끝난 여부와는 상관없이 FU기간이 제일 짧은(분석시점에서 제일 가깝게 등록)한 환자를 기준으로 모든 환자를 똑같이 등록이후 같은 기간만 FU해서 판단한다면?
    • 젤 짧은 4번 환자 2년
      • 1번은 2년내 발병 X
      • 2번은 2년내 발병 X
      • 3번은 2년내 발병 O
      • 4번은 2년내 발병 X image-20220423123630720
    • 통계적으로는 valid하겠지만, 안타깝게도, 2번 환자가 FU기간 기준 내에선 발병X지만, 실제로 3년차에 발병O한다는 사실을 분석에 반영못하게 된다. 사실 정보를 가지고 있지만 반영 못하게 된다

생존자료분석 기법을 이용하면, 통계적 valid & 모든 정보 반영

  • 1번과 4번를 보면, 분석시점에 발병X뿐만 아니라 분석시점 이후에는 발병했을 수도/안했을 수 있고 모른다를 분석에 반영할 수 있다.

    • 코호트 분석에서 더 이상적이긴하지만, 하나의 토픽이라서 지금 다루진 않는다. image-20220423130045870

  • 만약, 생존자료분석을 이용하지 않는다면, 최근등록된 제일 짧은FU을 기준으로 삼는게 아니라 FU기간을 반영한 Person-time incidence 위주의 분석도 가능하지만 다루진 않는다

코호트 스터디 논문예제

N E J med에 발표된 천식관련 코호트

image-20220423130412452

  • 뉴질랜드 어느 마을에서 72년 4월 ~ 73년 3월까지 태어난 나이 1000명 중 90%가 참여
  • 이후 5, 7, 9, 11, 13, 15, 18, 21, 26세에 FU
    • skin-prick test: 알러지 테스트
    • 폐기능 검사
    • airway responsiveness : 천식관련 outcome
  • 천식 발병O/X를 나누기 위해 여러가지 outcome 정의(2분형 -> 여러가지로)
    • 발병이후 쭉 위징 : 펄시스턴스 위징
    • 발병이후 쭉 위징있었다가 나아지면: 리미션
    • 발병이후 쭉 위징 나았다가 다시 생기면: 릴랩스 image-20220423130654803

결과 해석

  1. 모든 FU에 다 참여한 아기 -> 613명
  2. 전체 FU에 다 참여하지 않은 아기 - > 그외 로 나눠서 분석한다. image-20220423131457473
  • 모든 FU에 다 참여한 아기 대상으로
    • 여러가지 정의했던 outcome(이분형X)에 대해 얼마나 발생했는지 본다

image-20220423131616383

  • 각 환자 특성별로 나누어 여러 outcome의 분포를 확인
    • image-20220423131709352

코호트 논문에서 관심있는 결과해석은 OR 분석표다

  • 2분형 결과변수는 아니였지만, 2가지 outcome만 OR을 계산하게 되었음.

  • 1) 위험 요소 O인 경우가 (그렇지 않은 경우에 비해) OR이 몇 배고, 2) P value 통계적으로 유의한지 보고 3) 점추정값의 신뢰구간도 보자

  • 4) 해석은 (위험요소 O)~ 한 경우가 ~하지 않는 경우에 비해 outcome O될 확률이 OR 배였다
    • Persistence하고 Relapse에 대해 OR을 계산했는데,
    • 21세에 스모킹 O인 그룹이 (그렇지 않은 그룹에 비해) Persistence가 있을 오즈가 2.05배였다 image-20220423162154170
  • 다변량 모델 = 컨파운딩 등을 보기 위해(나머지 설명변수 고정상태를 만드는 모델) OR이 유의한 설명 변수만 집어넣고 다시 모델을 만든 것.
    • 각 변수에 대해 나머지 변수는 고정일 때로 해석가능하다
    • (다변량모델에서, 나머지 변수들이 고정되어있는데도~) 21세 흡연여부가 -> OR이 1보다 높고 통계적으로 유의하다.

Cross-sectional Study

  • 딱 1시점에서 데이터를 보고 위험요소 노출 O/X 여부 및 질병 발병 O/X 여부를 보게 됨
    • 변화하는 관심 인구집단의 한 시점의 스냅샷으로
    • incidence(속도)를 볼 수 없고, prevalence만 보게 된다.
      • 위험요소가 먼저였는지, 질병이 먼저였는지 알 수 없음 image-20220423164453610

bias 생길 확률이 높다

Temporal bias

  • 템포럴 바이러스는 위험요소노출 여부 vs 질병여부 어떤 것이 먼저인지 알 수 없을 때 오는 bias다
    • 예를 들어, 운동 여부 vs CHD 발병 여부
      • CHD진단 받은 환자들이 운동을 늘렸을 가능성이 있다.
        • 그 때 Crross-section study에 잡히면 ? 상관관계가 과소추정된다.
      • 운동량 적을 수록 -> CHD 많이 발견되는 상관관계를 발견했을 지라도
        • 운동량 적 -> CHD 많? or CHD 많 -> 운동량 적어짐 선후관계를 알 수가 없다 image-20220423164550664

survival/selection bias(prevalence-incidence bias)

  • 서바이벌 바이어스 / 셀렉션 바이어스 / 프리벌런스-인시턴스 바이어스 라고도 한다.
    • 위험요소 노출후 (질병O는 물론이고) 금방 사망해버린 환자들은 연구에 포함되기가 어렵다에서 오는 바이어스다
    • 예를 들어, 흡연 -> 폐기종 환자에 있어서
      • 흡연에 의해 폐기종 발생 환자들은 비흡연 -> 폐기종 환자들에 비해, 빨리 사망함 -> Cross-sectional study에 포함되기 힘들어진다
      • 살아남아 Cross-sectional study에 포함된 폐기종 환자들의 흡연률(흡연노출여부O)는, 전체 폐기종 발병 환자의 실제 흡연률보다 낮아지게 됨 -> 흡연과 폐기종의 상관관계를 과소추정하게 된다. image-20220423164929839

한계점이 많은 연구

  • bias를 많이 가질 수 밖에 없는 study이므로, Cross-sectional Study의 결과는 다른 연구 시간의 선후관계를 가진 파악 다른 연구로 가기 전의 전 단계로서, 살펴볼 가치가 있는 가설인지 확인용으로 쓴다.

크로스-섹셔널에서의 RR & OR

  • 크로스-섹셔널 연구에서는 incidence 계산 불가 -> risk 계산 불가 -> Relative Ratio의 RR대신 Prevalence Ratio로 부른다.
  • OR에 들어가는 p1, p0가 risk가 아닌 prevalence -> prevalence OR 로 부른다. image-20220423165805421

정리하기

  • Cohort study는 연구자가 위험요소에 노출된 그룹과 노출되지 않은 그룹을 일정기간 동안 추적관찰하여, 두 그룹 간 질병 발생률 (incidence)이나 질병에 의한 사망률 등을 비교하는 연구이다
  • Cohort study에서는 RR과 OR을 추정하고 신뢰구간을 계산할 수 있다
  • Confounding이 있을 경우 Mantel-Haenszel 방법으로 OR을 계산하거나 로지스틱 회귀분석으로 Adjusted OR을 계산할 수 있다
  • Interaction이 있을 경우 로지스틱 회귀분석으로 interaction의 유의성을 검정하고 stratum-specific OR을 계산할 수 있다
  • Cross-section study에서는 prevalence ratio와 prevalence OR을 계산할 수 있다

cohort study 연습문제

Q1. 다음 중 cohort study에 대한 설명으로 적절하지 않은 것은?

  1. 관찰연구의 일종이다.
    • 연구는 실험연구 vs 관찰연구로 나뉘는데 관찰만 한다.
      • cf) 실험연구 아니면 관찰연구 -> 현재시점 기준으로 이후(전향적) 아니면 이전(후향적) 관찰 -> 이후든 이전이든 1시점만(cross-sectional) 보냐 vs 2시점 이상(Longitudinal-이후Cohort/이전Case-control)보냐
  2. 여러개의 질병을 동시에 연구할 수 없다.
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 질병발생 시점이 명확 -> incidence(질병발생 속도와 관련) 추정이 가능
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 인과관계 추론에 결정적 단서
    • 위험요소 -> 질병발생의 선후가 정해져 있어서 위험요소 노출O/X이후 나올 수 있는 여러 질병들을 동시에 연구 가능하다.(타겟질병말고도 다른 질병이 나타나면 연구가능)

  3. 희귀병은 연구하기 어렵다.

    • 코호트는 위험요소 -> 질병발생의 선후가 정해져 있어서 희귀병은 기다려도 몇 명 발생안하므로, 많은 수의 연구대상자가 확보되어야하고, 계속 follow up해야한다
    • 대안: retrospective cohort study
      • 연구시작 시점 이후 발병하는 질병을 follow 하는 것은 동일하나,
      • 연구시작 시점 이전 위험요소 노출O/X를 <기록(데이터),기억(인터뷰)>를 통해 결정하는 것
        • 질병발생률은 prospective / 위험요소 노출여부는 retrospective로서 2가지가 섞여있는 디자인이라고도 볼 수 있다.
        • 위험요소 노출O/X를 기다렸다가 나누었다가 할 필요가 없이 바로 질병발생률만 관찰하므로 전체적인 연구기간이 짧아진다.
        • 위험요소 -> 질병발생의 선후가 정해져 있어서 질병 발생할때까지 계속 기다리다가 위험요소 노출여부가 중간에 바뀌어버릴 수 있다.
          • 예> 흡연으로 위험요소 노출O였다가, 연구 중간에 금연해서 not/un exposed group으로 가버릴 수 있다
          • 대안: 위험요소 노출여부를 여러번 측정 -> longitudinal study(여러번 측정시)

  4. 위험요소 노출부터 질병 발병까지의 기간이 너무 긴 질병은 연구하기 어렵다.

  • 정답 : 2
  • 해설 : Cohort study에서는 여러개의 질병을 동시에 연구하는 것이 가능하다.

Q2. 다음은 어떤 cohort study에서 흡연과 CHD의 관계를 알아보기 위해, 흡연자 300명, 비흡연자 500명을 10년간 추적하여 CHD 발병 여부를 조사한 것이다.

image-20220423170946709

  1. RR(relative risk) = (10/300)/(10/500) = 1.67
  2. OR(odds ratio) = (290×10)/(490×10) = 0.59
  3. 이 데이터에서는 흡연자의 CHD 발병 위험이 더 높다.
  4. 흡연 여부와 CHD 발병 여부가 독립인지 검정하려면 카이제곱 검정을 시행하면 된다.

풀이

  1. row -> col 순으로 보고 위험요소2개 = interaction문제인지 / 위험요소1 질병여부1의 상관관계 문제인지 확인한다
    • 위험요소노출O/X -> 질병여부O/X 의 상관관계는 위험요소O/X에 따라 group을 나눈 -> RR/OR로 판단한다.

  1. 상관관계를 위한 RR, OR계산은 위험요소노출 O/X group별로 나눈 질병 risk인 p1/p0를 계산한다

    • p1(위험요소 노출O risk) = 10/10+290
    • p0(위험요소 노출O risk) = 10/10+490
    • rr = p1/p0 = 500/300 = 1.67

  2. OR계산은 위험요소노출 O/X group별로 나눈 질병 odds(risk / 1-risk)를 구한뒤 비율을 구해야하지만 위험요소-질병 표에서는 a(oo)d(xx) / bc로 계산해버리면 된다.

    • a(oo) = 흡연->발병 = 10
    • d(xx) = 흡연x->발병x = 490
    • bc = 반대 대각선 = 10 x 290
  1. 추가적으로 범주(위험요소 노출O/X)형 범주형(질병O/X) -> 독립성 검정으로 상관관계 없음을 검정할 수 도 있다.(카이제곱 검정)
  • 정답 : 2
  • 해설 : 위의 표는 disease positive가 왼쪽이 아닌 오른쪽 열에 배치되어있다. 이를 일반적인 형태(disease positive가 왼쪽 열, exposed가 위쪽 행)로 바꾼 후 cross-product 공식으로 OR을 구하면 (10×490)/(10×290) = 1.69이다.