동적 가능도 기반 위험률 추정: 파라메트릭·비파라메트릭의 장점을 동시에 잡다

al e" categories: [“Research”]

📝 Abstract

The best known methods for estimating hazard rate functions in survival analysis models are either purely parametric or purely nonparametric. The parametric ones are sometimes too biased while the nonparametric ones are sometimes too variable. In the present paper a certain semiparametric approach to hazard rate estimation, proposed in Hjort (1991), is developed further, aiming to combine parametric and nonparametric features. It uses a dynamic local likelihood approach to fit the locally most suitable member in a given parametric class of hazard rates, and amounts to a version of nonparametric parameter smoothing within the parametric class. Thus the parametric hazard rate estimate at time $s$ inserts a parameter estimate that also depends on $s $. We study bias and variance properties of the resulting estimator and methods for choosing the local smoothing parameter. It is shown that dynamic likelihood estimation often leads to better performance than the purely nonparametric methods, while also having capacity for not losing much to the parametric methods in cases where the model being smoothed is adequate.

💡 Analysis

**

1. 연구 배경 및 필요성

• 파라메트릭 vs. 비파라메트릭: 파라메트릭 모델(예: Weibull, Gompertz‑Makeham 등)은 추정 효율이 높지만 모델 오차에 취약하고, 비파라메트릭(Nelson‑Aalen, 커널 스무딩)은 모델 자유도가 높아 편향은 적지만 표본 변동이 커진다.
• 반파라메트릭 접근: 기존 문헌(Hjort 1991, 1992)에서는 파라메트릭 모델을 ‘베이스라인’으로 두고, orthogonal expansion이나 베이지안 비파라메트릭 사전 등을 통해 보정하는 방법을 제시했다. 그러나 실제 적용에서는 시간에 따라 파라미터가 변하는 현상을 충분히 반영하지 못했다.

2. 핵심 아이디어 – 동적 가능도(local likelihood)

• 시간‑특정 파라미터 θ(s): 각 시점 s에서 s±h/2 구간에 속한 관측치만을 이용해 부분가능도(Likelihood)를 최대화함으로써, 해당 구간에 가장 적합한 파라미터를 추정한다.
• 비파라메트릭 파라미터 스무딩: 파라메트릭 형태는 유지하되, 파라미터 자체를 스무딩(커널 가중)함으로써 “파라메트릭 안에서 비파라메트릭” 접근을 구현한다.
• 커널 선택: 대칭 커널 K(u)(예: Uniform, Epanechnikov 등)를 사용해 가중치를 부여하고, 커널의 2차 모멘트 β_K와 제곱 적분 γ_K가 편향·분산 공식에 직접 등장한다.

3. 이론적 결과

📄 Content

1. 서론 및 요약
이 논문은 수명 데이터의 위험률 함수(hazard rate function)를 추정하는 반파라메트릭(semi‑parametric) 방법들의 한 종류에 관한 것이다. 위험률을 추정하는 가장 잘 알려진 방법은 순수하게 파라메트릭(pararmetric)인 방법이나 순수하게 비파라메트릭(non‑parametric)인 방법이다. 파라메트릭 방법은 모델 자체가 보통 불완전하기 때문에 편향(bias)이 발생하기 쉽고, 비파라메트릭 방법은 추정 분산이 크게 나타나는 경우가 많다. 따라서 파라메트릭과 비파라메트릭 사이 어딘가에 위치하는 방법에 대한 여지가 있다. 진정한 위험률이 파라메트릭 모델 근처에 있을 경우 이러한 방법이 비파라메트릭 방법보다 더 나을 것이며, 파라메트릭 모델이 실제와 일치할 때는 파라메트릭 방법만큼 크게 손해 보지 않을 것이라고 기대할 수 있다.

일반적인 계수 과정(counting‑process) 모델 체계에서도 결과를 얻을 수 있지만, 여기서는 가장 단순하고 아마도 가장 중요한 검열(censoring)된 수명 데이터 모델인 ‘무작위 검열(random censorship)’ 모델을 중심으로 아이디어를 설명하고 검증한다. 이 모델은 모집단으로부터 추출된 수명 (X^{0}{1},\dots ,X^{0}{n}) 이 i.i.d.이며, 밀도 (f(\cdot )), 누적분포 (F(\cdot )), 위험률 함수 (\alpha(\cdot )) 를 갖는다고 가정한다. 위험률은

[ \alpha(s)=\frac{f(s)}{F[s,\infty)} ,\qquad \alpha(s),ds ]

가 시간 (s) 에서 아직 위험에 노출된 개체가 ([s,s+ds)) 구간에 사망할 확률을 의미한다. 실제 수명 (X^{0}{i}) 는 검열 변수 (C{i}) 로 인해 직접 관측되지 않을 수 있다; 관측되는 것은

[ X_{i}= \min (X^{0}{i},C{i}),\qquad \delta_{i}=I{X^{0}{i}\le C{i}}. ]

편의를 위해 검열 변수 (C_{i}) 가 수명과 독립이며, 누적분포 (G) 를 갖는 i.i.d.라고 가정한다. 따라서 ((X_{i},\delta_{i})) 들은 i.i.d.이며, 관측은 유한한 시간 구간 ([0,T]) ( (T) 는 알려지고 유한함) 에서만 이루어진다고 가정한다. 이는 마팅게일 수렴 이론을 적용하기에 편리하며 실제적인 제한은 아니다.

파라메트릭 접근은 위험률을

[ \alpha(s)=\alpha(s,\theta) ]

와 같이 어떤 파라미터 (\theta) (1‑차원 혹은 다차원) 로 색인되는 가족으로 가정한다. 전형적인 예로는 지수분포, Weibull, (\alpha(s)=\theta_{1}/(1+\theta_{2}s)) 형태의 단순 frailty 모델, 구간별 상수 위험률 모델, Gompertz‑Makeham 분포, 감마분포, 로그정규분포 등이 있다. 검열 데이터를 포함한 최대우도법(maximum likelihood) 의 성질은 Borgan(1984) 등에 의해 모델이 정확히 맞을 경우(즉, 실제 위험률이 (\alpha(s)=\alpha(s,\theta_{0})) 로 표현될 때) 연구되었다. 실제로는 모델이 완벽하지 않으므로, ‘최소 거짓(least false)’ 혹은 ‘가장 적합한(suitable)’ 파라미터를 찾는 문제를 다루어야 한다. 이러한 넓은 상황에서 여러 추정법의 대표본(asymptotic) 거동은 Hjort(1992) 에서 탐구되었다. 제2절에서 그 결과를 간략히 검토하고 이후 절에서 활용한다.

3절에서는 파라메트릭 위험률 추정을 위한 동적(likelihood) 접근법을 제시한다. 기본 아이디어는 임의의 파라메트릭 위험함수 (\alpha(s,\theta)) 를 잡고, 시간 (s) 에서 지역 파라미터 추정치 (\theta(s)) 를 삽입하여

[ \ell_{s}(\theta)=\int_{s-\frac{h}{2}}^{s+\frac{h}{2}} \bigl{\log \alpha(u,\theta),dN(u)-Y(u)\alpha(u,\theta),du\bigr} ]

와 같은 부분 로그우도를 최대화하는 것이다. 여기서 (N(u)) 는 사망 사건을 세는 계수 과정, (Y(u)) 는 위험에 노출된 개체 수를 나타내는 위험(at‑risk) 과정이다. 즉, (s) 이전에 살아남은 개체들만을 이용해 (\theta(s)) 를 추정한다는 의미이며, 이는 주어진 파라메트릭 가족 안에서 일종의 비파라메트릭 파라미터 스무딩이라고 볼 수 있다. 커널 함수를 이용한 보다 일반적인 스무딩 추정법도 논의한다. 1차원 파라메트릭 가족에 대해서는 3절에서, 다차원 가족에 대해서는 4절에서 편향과 분산을 각각 구한다. 주요 결과는

[ E{\hat\alpha(s)}= \alpha(s)+\frac{1}{2}\beta_{K}h^{2}b(s),\qquad \operatorname{Var}{\hat\alpha(s)}= \frac{\gamma_{K}}{n h},\alpha(s),y(s), ]

여기서 (\beta_{K},\gamma_{K}) 는 사용된 커널의 특성, (y(s)) 는 시간 (s) 에서 위험에 노출된 개체 비율의 극한값, (b(s)) 는 위험률과 파라메트릭 모델에 의존하는 편향 계수이다. 이는 경험적으로 가장 많이 쓰이는 비파라메트릭 방법(경험적 누적 위험률을 스무딩하는 방법)

[ E{\hat\alpha_{\text{np}}(s)}= \alpha(s)+\frac{1}{2}\beta_{K}h^{2}\alpha(s),\qquad \operatorname{Var}{\hat\alpha_{\text{np}}(s)}= \frac{\gamma_{K}}{n h},\alpha(s),y(s) ]

와 형태가 거의 동일함을 보여준다.

5절에서는 새로운 방법이 전통적인 비파라메트릭 방법보다 우수한 상황을 규정한다. 6절에서는 지역 스무딩 폭 (h) 를 선택하는 여러 방법을 논의한다. 특히 각 (s) 에 대해 (s\pm h/2) 구간을 확대해 가며 적합도 검정이 모델을 기각할 때까지 진행하는 방법을 제시한다. 전반적으로 동적 로그우도 추정량은 비파라메트릭 방법보다 더 나은 성능을 보이는 경우가 많으며, 파라메트릭 모델이 실제 위험률에 가깝다면 파라메트릭 방법과 큰 차이를 보이지 않는다. 마지막으로 7절에서는 보충 결과와 몇 가지 논평을 제공한다.

본 논문은 Hjort(1991) 에서 제시된 기본 결과들을 여러 방면에서 확장한다. 그 논문에서는 또한 두 가지 추가적인 반파라메트릭 추정 스킴을 제안했는데, 하나는 초기 파라메트릭 추정치를 정규 직교 전개(orthogonal expansion)로 보정하는 방법이고, 다른 하나는 주어진 파라메트릭 위험 모델 주변에 비파라메트릭 사전분포를 두는 베이지안 절차이다.

순수 비파라메트릭 및 순수 파라메트릭 추정

이 절에서는 기본 기호들을 소개하고, 비파라메트릭 경우의 누적 위험률(Nelson‑Aalen) 추정량과 파라메트릭 경우의 최대우도 및 가중 최대우도(maximum weighted likelihood) 추정량의 성질을 검토한다. 이러한 결과들은 이후 절에서 활용될 것이며, 파라메트릭 모델이 실제와 일치하지 않을 때의 거동도 함께 고려한다.

계수 과정과 위험 과정은

[ N(t)=\sum_{i=1}^{n} I{X_{i}\le t,;\delta_{i}=1},\qquad Y(s)=\sum_{i=1}^{n} I{X_{i}\ge s} ]

로 정의한다. 누적 위험률 (A(t)=\int_{0}^{t}\alpha(s),ds) 의 Nelson‑Aalen 추정량은

[ \hat A(t)=\int_{0}^{t}\frac{dN(s)}{Y(s)} . ]

이 추정량의 성질은 마팅게일

[ B(t)=N(t)-\int_{0}^{t}Y(s)\alpha(s),ds ]

을 이용해 설명한다. (y(s)=\lim_{n\to\infty}Y(s)/n) 를 위험에 노출된 개체 비율의 극한값이라 하면, (B(t)/\sqrt{n}) 은 독립적 증분을 갖는 가우시안 마팅게일 (V(t)) 로 수렴하고

[ \operatorname{Var}{dV(s)}=y(s)\alpha(s),ds . ]

이로부터

[ \sqrt{n}\bigl{\hat A(t)-A(t)\bigr};\Longrightarrow; \int_{0}^{t}y(s)^{-1},dV(s) ]

가 얻어진다. 자세한 내용은 Andersen, Borgan, Gill & Keiding(1993, Chapter II)를 참고한다. 위험률 자체를 추정하려면 Nelson‑Aalen 추정량을 스무딩하고 미분하면 된다(식 5.1).

파라메트릭 모델은 (\alpha(t)=\alpha(t,\theta)) 형태이며, (\theta\in\mathbb{R}^{p}) 를 파라미터라 한다. 관측된 데이터에 대한 로그우도는

[ L_{n}(\theta)=\int_{0}^{T}\bigl{\log\alpha(t,\theta),dN(t)-Y(t)\alpha(t,\theta),dt\bigr} ]

이며, 이를 최대화한 (\hat\theta) 가 최대우도 추정량이다. 대표본 거동을 설명하기 위해

[ U_{n}(\theta)=\int_{0}^{T}\psi(t,\theta)\bigl{dN(t)-Y(t)\alpha(t,\theta),dt\bigr}, \qquad \psi(t,\theta)=\frac{\partial\log\alpha(t,\theta)}{\partial\theta}, ]

를 정의한다. (\hat\theta) 가 0이 되는 점을 만족하면 (\hat\theta) 는 (\theta_{0}) (‘least‑false’ 파라미터) 로 수렴한다. 이때 (\theta_{0}) 는

[ u(\theta)=\int_{0}^

이 글은 AI가 자동 번역