Robins [5]는 임상시험에서 여러 가지 이유로, 예를 들면 대조군과 처리군의 비교성이 깨어지는 경우 인과관계를 밝힐 수 있는 Rank Preserving Structural Failure Time (RPSFT) 모델을 제안하였다. 만약 대조군에 배정된 환자가 악화 또는 재발(progression)된 시점에서 실험군 약으로 치료를 변경한 경우(single crossover), 최종 관측된 대조군의 OS는 실험약으로 인한 효과가 포함되어 있다. 따라서 만약 이 환자가 실험약으로 치료받지 않았다는 가정(counterfactual)에서의 잠재적 OS를 생각할 수 있다. 먼저 관측된 OS는 다음과 같이 표현할 수 있다.

여기서 T_i는 i번째 대상의 관측된 OS로서 대조약을 치료받으면서 악화되어 치료가 변경되기까지의 기간(Tioff)과 치료가 변경된 시점부터 사망 또는 중도절단 시점까지의 시간(Tion)의 합으로 정의된다. 잠재적 전체 생존기간(U_i)은 다음과 같이 모형화할 수 있다.

여기서 ψ는 가속화 요인(acceleration factor)이다. 따라서 가속화 요인 ψ가 0보다 작다는 것은 새로운 치료약 즉, 새로운 항암제의 효과가 더 우수하다는 것을 의미한다. 가속화 요인 ψ의 추정은, 만약 두 군의 배정이 무작위하게 이루어졌다면 grid searching을 통해 할 수 있다. 일정 범위에서 선택한 ψ값을 ψ^*라고 하면, 잠재적 변수(latent variable) U_i의 값을 얻을 수 있고, 두 군의 로그-순위(log-rank) 검정을 이용할 수 있다. 이 때 얻어진 p-value가 최대가 될 수 있는 ψ^*를 가속화 요인으로 추정한다. 여기서 ψ^*가 음수일 경우 exp(ψ^*)는 1보다 작아져 U_i의 시간은 줄어들게 된다. 따라서 기존 중도절단 시간 C_i를 재설정하는 재 중도절단(re-censoring) 과정이 필요하다. 재 중도절단 시간 C_i(ψ^*)는 다음과 같이 정의할 수 있다.

Jonsson [6]에 따르면 연구 중인 치료제는 사망할 때까지의 시간을 지속적으로 감소시키고 재발 전과 후 시점에서 모든 환자들에게 동일하다는 것이 RPSFT 모델의 중요 가정이며, 사망이나 악화 또는 재발같은 사건이 발생할 때까지의 시간이 비례적으로 감소하는 경우인 연구에 해당 모델을 적용하는 것이 타당하다고 설명하였다.

본 연구에서는 RPSFT 모델을 바탕으로 한 잠재적 변수로 GMI를 새롭게 추정하는 방법을 제안하였다. 먼저 다음과 같은 두 가지 형태로 잠재적 변수를 생각할 수 있다. 편의상 1차로 투여받는 항암제를 A, 2차로 투여받는 항암제를 B라고 지정하자. 첫 번째 형태는 연구시작 시점에서 1차 항암제로 A를 투여받고 첫 번째 재발하기까지 관찰된 생존시간 PFS₁에 해당하는 잠재적 변수로 만약 연구시작 시점에서 항암제 B를 투여 받았다고 가정하였을 때 첫 번째 재발하기까지의 생존시간은 PFS₁ × e^β 로 생각할 수 있다. 또 다른 형태의 잠재적 변수는 1차 항암치료 후 첫 번째 재발시점에서 2차 항암제 B를 투여받고 두 번째 재발하기까지의 생존시간 PFS₂에 해당하는 잠재적 변수로 만약 첫 번째 재발시점에서 항암제 B 대신 항암제 A를 지속적으로 투여 받았다고 가정하였을 때 두 번째 재발하기까지의 생존시간은 PFS₂ ×e^β 로 생각할 수 있다. 이렇게 잠재적 변수를 정의하게 되면 각 대상자로부터 동일 시점에서 서로 다른 항암제를 투여 받았을 때 생존시간을 짝으로 정의할 수 있고 앞에서 설명한 RPSFT 모형을 활용하여 추론할 수 있다. 본 연구에서는 2차 항암제의 유효성 평가에 초점을 둔 유효성 평가변수로 GMI를 고려하므로 두 번째 방법으로 잠재적 변수를 정의하였다. 생존시간 T를 PFS₁과 PFS₂를 더한 전체 관측된 시간으로 정의할 때, U는 PFS₁과 PFS₁^*을 더한 잠재적 생존시간으로 정의하며 i번째 관측대상의 관측시간 T와 U는 다음과 같이 표현할 수 있다.

RPSFT 모델을 바탕으로 하는 방법을 적용하여 추정한 i번째 관측대상의 새로운 GMI를 GMI₁^*라고 정의하며 다음과 같이 정의한다.

여기서 β를 추정하기 위해, 먼저 PFS_1i와 PFS_1i^* 간에는 감마 프레일티(γ_i) 공유로 인한 연관관계가 있다고 가정한다. PFS_1i와 PFS_1i^*를 구분하는 그룹변수를 만들고, Cox 회귀분석 시 그룹변수와 감마 프레일티를 모형에 포함하여 분석한다. 두 그룹 간의 차이를 나타내는 그룹변수의 p-value가 가장 큰 값을 가질 때의 특정 범위 내 가속화 요인 β값을 grid searching하여 GMI₁^*를 추정한다. 추정한 β는 모든 관측대상에서 동일한 값을 가지므로, 추정된 GMI₁^* 역시 모든 관측대상에서 동일한 상수의 값을 가지게 된다. β가 양수인 경우 e^β 는 1보다 큰 값을 가지므로 기존 중도절단 시간 Ci를 그대로 사용한다. 하지만 β가 음수인 경우 e^β 는 1보다 작은 값을 가지므로 PFS_1i^*는 중도절단 시간 C_i- PFS_1i보다 작아질 수 있어, 재 중도절단에 대한 고려가 필요하다. 본 연구에서 재 중도절단 조정 방법은 다음과 같다. 먼저 PFS_1i^*의 중도절단 여부는 PFS_2i의 중도절단 여부를 그대로 사용하며 PFS_1i이 중도절단된 경우 PFS_1i^*는 결측치(NA)로 처리한다. 하지만 PFS_1i이 중도절단 되지 않고 PFS_2i가 중도 절단되었다면 PFS_1i^*의 재 중도절단 된 시간은 다음과 같이 표현한다.

콕스 비례위험(Cox proportional hazards, Cox PH) 모형은 생존 시간에 영향을 미치는 공변량(covariate)을 찾아내는 분석 방법으로서, 반응변수인 생존 시간과 공변량의 관계가 위험함수(hazard function)를 통해 표현된다는 특징을 가지고 있다. Cox PH 모형은 생존 시간에 대한 기저분포(baseline distribution) 가정을 필요로 하지 않는 준모수(semi-parametric) 방법이다. 하지만 위험률의 비는 시간에 의존하지 않고 항상 일정하다는 비례위험에 대한 가정이 만족하지 않으면 문제가 발생할 수 있으며, 추정된 값에 대한 해석이 어렵다는 단점이 존재한다. 이와 같은 경우, 생존시간을 직접 모형화하는 가속실패 시간모형(accelerated failure time, AFT)을 고려할 수 있다[7]. AFT 모델은 Cox PH 모형에 비해 직관적인 해석이 가능하다는 이유로 임상에서 유용하게 사용할 수 있으나 모수적 분포 가정을 해야 하는 한계점이 존재한다[8].

본 연구에서 두 번째 방법으로 GMI를 추정하는 방법으로 먼저 다음과 같은 모형을 구축한다[9].

여기서 w_i1은 i번째 관측대상의 첫 번째 관측된 사건 발생까지의 생존 시간을, w_i2는 i번째 관측대상의 첫 번째 사건 발생부터 두 번째 사건 발생까지의 관측된 갭타임(gap time)으로 정의한다. 즉 w_ij는 i번째 관측대상의 j번째 갭타임을 의미한다. I는 indicator 함수로써 I(j =2)는 i번째 관측대상의 두 번째 관측시간이 관측되었을 경우 1의 값을 가지며, 관측되지 않았을 경우 0의 값을 가지게 된다. σ_j는 척도모수이고 ε_ij는 오차항에 해당한다.

본 연구에서는 로그노말 분포를 따르는 AFT 모델을 이용하여 GMI₂^*라는 새로운 GMI를 추정한다. i번째 관측대상의 PFS_1i와 PFS_2i는 w_i1, w_i2와 동일하다. 위의 AFT 모델은 다음과 같이 표현할 수 있다.

따라서 log(w_i1)와 log(w_i2)의 기댓값은 다음과 같다.

따라서 기하평균(geometric mean) 방법으로 구한 새로운 GMI₂^*는 다음과 같이 정의할 수 있다.

위의 AFT 모델을 이용한 방법으로 구한 GMI₂^*는 잠재적 변수를 이용한 GMI₁^*와는 달리 모든 관측대상이 서로 다른 값을 가지고 있다.

본 논문의 모의실험에서는 다음과 같은 단계로 PFS_1i와 PFS_2i 그리고 중도절단 자료를 생성하였다.

① i번째 대상의 프레일티를 평균이 1이고 분산이 α인 감마분포로부터 발생하였다.

② PFS_1i과 PFS_2i는 다음과 같이 발생하였다. 즉,

여기서 U_1i, U_2i,는 각각 균등분포(uniform distribution)를 따르는 0과 1 사이의 임의의 실수이며, λ₁/λ₂는 hazard ratio로 참 GMI이며 본 연구에서는 1.0과 1.33으로 설정하였다. 두 PFS 간의 상관관계를 나타내는 Kendall's τ는 0.1과 0.4를 가정하였다.

③ PFS₁의 중도절단 시간 C₁은 모의실험 시 고려하는 모든 모수에 대하여 10만 개의 생존 시간을 생성한 후 중도절단 비율이 5%, 10%, 20%가 되도록 설정하였다. PFS₂의 중도절단을 결정하는 C₂는 생성된 PFS₁과 PFS₂를 더한 기간이 C₁의 k배가 되도록 하였으며 k는 무한대(∞)와 1.5로 하였다. 정의한 PFS₁의 중도절단 비율에 따른 경험적으로 추정된 중도절단 시간은 다음과 같다(Table 1).

④ 기존 GMI 추정[2]은 PFS₁의 값이 중도절단된 경우는 제외하였다. 잠재적 변수를 이용한 방법, 그리고 AFT 모형을 이용한 방법은 PFS₁의 값이 중도절단 되어도 모두 분석에 포함하였다. 세 가지 방법으로 추정된 추정치의 편의와 평균제곱오차(mean square error, MSE)를 1,000번 반복 수행을 통해 그 특성을 비교하였다. 잠재적 변수를 이용한 방법의 경우 모든 관측대상에서 동일한 값을 가지므로, 1,000번 추정된 GMI₁^*의 평균값을 구하였다. 기존의 GMI와 AFT 모델을 이용한 GMI₂^*는 각각 Kaplan-Meier 방법을 통해 추정된 생존함수의 중앙값(median)을 구하고, 1,000번 반복하여 추정한 중앙값의 평균값으로 구하였다. 실제 암 2상 임상시험에서는 많은 환자 자료가 얻어지는 경우가 적기 때문에 표본의 크기는 100으로 제한하였다. 또한 모의실험에서 사용한 모든 분석은 R 프로그램 3.2.3 버전을 사용하여 분석하였고, survival과 reshape 패키지를 사용하였다.