A Bayesian location-scale joint model for time-to-event and multivariate longitudinal data with association based on within-individual variability

기술적 요약: 시점-사건 및 다변량 종단 데이터에 대한 베이지안 위치-척도 결합 모델

문제 정의
생물 의학적 종단 연구에서, 건강 지표의 개인 내 변동성(within-individual variability, WIV)은 질병 진행 및 사망과 같은 사건 발생 시간(time-to-event) 결과의 결정적인 예측 인자로 점점 더 중요하게 인식되고 있습니다. 표준 편차와 같은 단순 요약 통계량이 WIV를 정량화하기 위해 자주 사용되지만, 이들은 심각한 한계를 가집니다. 즉, 시간에 따른 변화를 설명하지 못하고, 방문 횟수에 민감하며, "불멸 시간(immortal time)"과 같은 편향을 유발할 수 있습니다. 또한, 이러한 요약 통계량이 생존 모델의 공변량으로 사용될 때, 생성되는 해저드 비율(hazard ratio) 추정치는 회귀 희석 편향(regression dilution bias)을 겪게 되어 부정확한 연관성을 초래합니다. 기존의 결합 모델링 접근 방식은 단일 바이오마커에 대한 WIV를 다루거나 빈도주의 방법을 사용해 왔으나, 변동성 자체가 사건 결과에 영향을 미치는 베이지안 결합 모델 프레임워크 내에서 다변량(상관관계가 있는) 종단 바이오마커의 WIV를 성공적으로 통합한 사례는 거의 없습니다.

방법론
저자들은 **WIV 기반 연관성을 가진 베이지안 다변량 결합 모델(JM-WIV)**을 제안합니다. 이 모델은 두 개의 상호 연결된 하위 모델로 구성됩니다:

1. 종단 하위 모델 (혼합 효과 위치-척도 모델 - MELSM):

   • 분산을 일정하다고 가정하는 대신, 모델은 $K$개의 종단 바이오마커의 평균($\mu$)과 표준 편차($\sigma$)를 모두 시간, 공변량 및 피험자 특이적 랜덤 효과의 함수로 취급합니다.
   • 형식은 다음과 같습니다:
     $$y_{ijk} = \eta^{\mu}_{ik}(t_{ij}) + \varepsilon_{ijk}, \quad g(\sigma_{ijk}) = \eta^{\sigma}_{ik}(t_{ij})$$
     여기서 $g(\cdot) = \log(\cdot)$는 비음수 분산을 보장합니다.
   • 공유된 랜덤 효과 벡터 $\mathbf{b}_i$는 각 바이오마커의 평균과 변동성 프로세스를 연결하며, 서로 다른 바이오마커 간의 상관관계(예: 바이오마커 A의 평균과 바이오마커 B의 변동성 사이의 상관관계)를 모델링할 수 있게 합니다.
   • 공분산 구조는 비구조적(unstructured)이며, 이를 통해 여러 바이오마커에 걸친 평균 및 표준 편차의 랜덤 효과 간의 복잡한 상관관계를 허용합니다.

2. 사건 하위 모델 (비례 해저드 모델):

   • 사건 발생 시간은 B-스플라인(B-splines)을 통해 유연한 기저 해저드(baseline hazard)를 지정한 콕스 비례 해저드 모델을 사용하여 모델링됩니다.
   • 결정적으로, 연관 항 $\eta^{\alpha}_i$는 평균의 현재 값뿐만 아니라 표준 편차(변동성)의 현재 값으로부터 종단 하위 모델의 정보를 통합합니다.
   • 저자들은 세 가지 잠재적 연관 구조를 정의합니다: 평균과 SD의 현재 값(CV), 선형 예측자(LP), 또는 랜덤 효과(RE)를 직접 사용하는 방식입니다.

3. 구현:

   • 이 모델은 효율적인 파라미터 탐색을 위해 해밀토니안 몬테카를로(HMC)를 사용하는 Stan을 통해 구현되었습니다.
   • 우도는 랜덤 효과에 대해 사건 및 종단 프로세스로 인수분해됩니다.
   • 사전 분포에는 고정 효과를 위한 정규 분포, 랜덤 효과 표준 편차를 위한 half-t 분포, 그리고 상관 행렬을 위한 LKJ 사전 분포가 포함됩니다.

주요 기여

• 새로운 프레임워크: 본 논문은 다변량 종단 데이터를 처리하면서, 변동성이 사건 결과와의 연관성을 명시적으로 주도하는 베이지안 결합 모델을 도입했습니다.
• 회귀 희석 해결: 두 단계 요약 통계 방식을 사용하는 대신 결합 프레임워크 내에서 WIV를 직접 모델링함으로써, 회귀 희석 편향을 완화하는 것을 목표로 합니다.
• 유연한 연관성: 모델은 해저드 함수가 바이오마커의 변동성에 의존할 수 있게 하며, 이는 표준 결합 모델이나 WIV를 생존 하위 모델에 포함하지 않은 이전의 위치-척도 확장 모델에는 없는 특징입니다.
• 소프트웨어 가용성: 저자들은 기존의 결합 모델링 도구를 이변량(및 잠재적 다변량) 위치-척도 설정으로 확장하는 Stan 기반의 오픈 소스 구현체(rstanjmwiv 패키지)를 제공합니다.

결과

• 시뮬레이션 연구:

  • 시뮬레이션 연구에서는 제안된 JM-WIV 모델을 일정한 분산을 가정하는 표준 결합 모델(JMbayes2 구현)과 비교했습니다.
  • 참값 데이터 생성 과정에 WIV와 사건 간의 연관성이 포함되어 있을 때, 표준 결합 모델은 생존 계수에 대해 편향된 추정치를 생성하였고 신뢰 구간의 포함 확률(coverage probability)이 명목 수준인 95%를 훨씬 밑돌았습니다.
  • JM-WIV 모델은 편향이 거의 없는 추정치를 제공하였으며 신뢰 구간 포함 확률이 명목 수준과 일치하여, WIV가 유의미한 예측 인자일 때 파라미터를 정확하게 복구할 수 있음을 입증했습니다.

• 낭성 섬유증(Cystic Fibrosis, CF) 데이터 적용:

  • 모델은 영국 낭성 섬유증 레지스트리 데이터(여성 3,282명, 방문 29,094회)를 활용하여 폐 기능(FEV1)과 영양 상태(BMI)의 연관성을 평가했습니다.
  • 종단적 발견: 모델은 연령에 따른 FEV1의 감소와 BMI의 증가를 성공적으로 포착했습니다. 또한 FEV1의 변동성은 연령에 따라 감소하는 반면, BMI의 변동성은 F508 동형접합성(homozygosity)과 특정한 연관성을 보임을 정량화했습니다.
  • 생존 관련 발견:
    • 낮은 평균 FEV1과 낮은 평균 BMI는 높은 사망 위험(Hazard Ratio < 1)과 관련이 있었습니다.
    • FEV1의 개인 내 변동성(WIV)이 높을수록 사망 위험이 유의미하게 증가했습니다 (HR = 1.374, 95% CI: 1.268–1.499).
    • BMI 변동성과 사망 사이의 명확한 연관성은 발견되지 않았습니다.
  • 이 결과는 평균적인 수준과는 독립적으로, 폐 기능의 변동(fluctuations)이 사망의 강력한 예측 인임을 확인시켜 줍니다.

의의 및 주장
본 논문은 JM-WIV 모델이 다변량 환경에서 개인 내 변동성을 적절히 정량화하고, 요약 통계 접근 방식의 편향 없이 이를 생존 결과와 직접 연결함으로써 방법론적 진보를 이루었다고 주장합니다. 저자들은 WIV의 연관성을 무시할 경우 생존 계수의 편향과 신뢰 구간 포함 능력의 상실이 발생한다는 시뮬레이션 결과를 제시했습니다. 임상 적용 측면에서, 이 모델은 폐 기능의 변동성이 평균적인 추세와는 별개로 사망의 뚜렷한 위험 인자임을 확인시켜 주며, 평균적인 경향뿐만 아니라 변동성을 모니터링하는 것의 임상적 유용성을 강화합니다.

저자들은 위치-척도 모델에서 WIV가 측정 오차(measurement error)와 구별되지 않을 수 있다는 점과, 바이오마커나 랜덤 효과의 수가 많아질 경우 계산 부하가 크게 증가한다는 한계점을 인정했습니다. 향에서는 새로운 치료법의 영향을 모델링하기 위해 비정규 분포, 경쟁 위험(예: 이식), 또는 변화점(change-point) 구조를 탐구할 수 있을 것이라고 제언했습니다.