2x2 table에서 민감도와 특이도는 본래 McNemar’s test를 통해 \(z = \frac{(f_{21} - f_{12})}{\sqrt{f_{12} + f_{21}}}\)로 정의된 z값의 p-value를 확인하는 방법으로 진행되었습니다. 하지만 위 방식으로는 민감도와 특이도의 신뢰 구간에 대한 정보를 알 수 없기에, 신뢰 구간까지 정할 수 있는 방법인 Wald confidence interval을 적용함으로써 더 많은 정보를 얻을 수 있습니다.
개요
이 문서에서는 진단 테스트 1과 테스트 2의 민감도 차이를 기반으로 신뢰 구간(CI)을 계산하는 방법을 다룹니다. 이를 위해 두 테스트의 결과를 기반으로 2x2 테이블을 구성하고, 각 민감도를 계산합니다. 최종적으로 두 민감도의 차이에 대한 신뢰 구간을 계산하는 수식을 도출합니다.
1. 기본 개념
각 테스트의 민감도 \(p_1\)과 \(p_2\)는 아래와 같이 정의됩니다:
- \(p_1\)은 테스트 1이 진정한 양성을 탐지한 비율입니다.
- \(p_2\)은 테스트 2가 진정한 양성을 탐지한 비율입니다.
2x2 테이블은 아래와 같이 구성될 수 있습니다:
\[\begin{array}{c|c|c|c} & \text{Test 2 Positive} & \text{Test 2 Negative} & \text{Total} \\ \hline \text{Test 1 Positive} & a & b & a + b \\ \hline \text{Test 1 Negative} & c & d & c + d \\ \hline \text{Total} & a + c & b + d & n \\ \end{array}\]여기서 각 테스트의 민감도는 다음과 같이 계산됩니다.(특이도 역시 모든 데이터들이 양성인 것을 음성으로 바꾸어 생각하면 아래 과정이 동일해집니다)
\[p_1 = \frac{a + b}{n}, \quad p_2 = \frac{a + c}{n}\]
2. 신뢰 구간의 도출
두 민감도의 차이에 대한 신뢰 구간을 계산하기 위해 Wald 방법을 사용할 수 있습니다. Wald 신뢰 구간은 아래와 같이 정의됩니다:
2.1. 분산의 기본 식
두 확률 변수 \(p_1\)과 \(p_2\)의 차이의 분산은 다음과 같이 계산됩니다:
\[ \text{Var}(p_2 - p_1) = \text{Var}(p_2) + \text{Var}(p_1) - 2\text{Cov}(p_2, p_1) \]
2.2. 각 분산의 계산
우리는 \(p_1\)과 \(p_2\)를 비율로 간주할 수 있습니다. 이 비율들의 분산은 다음과 같이 계산됩니다:
\[ \text{Var}(p_1) = \frac{p_1(1 - p_1)}{n} = \frac{p_{12} + p_{11}}{n} \cdot \left(1 - \frac{p_{12} + p_{11}}{n}\right) \]
\[ \text{Var}(p_2) = \frac{p_2(1 - p_2)}{n} = \frac{p_{21} + p_{11}}{n} \cdot \left(1 - \frac{p_{21} + p_{11}}{n}\right) \]
2.3. 공분산 \(\text{Cov}(p_2, p_1)\) 의 계산
\(p_2\) 와 \(p_1\) 의 공분산은 다음과 같이 정의됩니다:
\[ \text{Cov}(p_2, p_1) = \text{Cov}(p_{11} + p_{21}, p_{11} + p_{12}) \]
이를 확장하면 다음과 같은 네 가지 항으로 분리할 수 있습니다:
\[ \text{Cov}(p_2, p_1) = \text{Cov}(p_{11}, p_{11}) + \text{Cov}(p_{11}, p_{12}) + \text{Cov}(p_{21}, p_{11}) + \text{Cov}(p_{21}, p_{12}) \]
공분산은 다음과 같이 정리할 수 있습니다.
\[ \text{Var}(p_{11}) = \frac{p_{11}(1 - p_{11})}{n} \]
\[ \text{Cov}(p_{11}, p_{12}) = -\frac{p_{11}p_{12}}{n} \]
\[ \text{Cov}(p_{11}, p_{21}) = -\frac{p_{11}p_{21}}{n} \]
\[ \text{Cov}(p_{21}, p_{12}) = -\frac{p_{21}p_{12}}{n} \]
이를 바탕으로 \(\text{Cov}(p_2, p_1)\)를 계산하면:
\[ \text{Cov}(p_2, p_1) = \frac{p_{11}(1 - p_{11} - p_{12} - p_{21}) - p_{21}p_{12}}{n} \]
여기서 \(p_{11} + p_{12} + p_{21} + p_{22} = 1\)임을 이용하여, 공분산을 다음과 같이 단순화할 수 있습니다:
\[ \text{Cov}(p_2, p_1) = \frac{p_{11}p_{22} - p_{21}p_{12}}{n} \]
2.4. 최종적으로 \(\text{Var}(p_2 - p_1)\) 의 유도
이제 분산 공식을 대입하여 최종적으로 \(p_2 - p_1\)의 분산을 계산합니다:
\[ \text{Var}(p_2 - p_1) = \text{Var}(p_2) + \text{Var}(p_1) - 2\text{Cov}(p_2, p_1) \]
따라서 \(p_2 - p_1\)의 분산은 다음과 같습니다:
\[ \text{Var}(p_2 - p_1) = \frac{(p_{12} + p_{21}) - (p_{21} - p_{12})^2/n}{n} \]
\[ \text{Var}(p_2 - p_1) = \frac{(b + c) - (b - c)^2/n}{n} \]
이를 적용한 민감도의 차이에 대한 신뢰 구간은 다음과 같이 계산됩니다. \[ \ CI_{1-\alpha/2}(\hat{\theta}) = \left[ \hat{\theta} \pm z_{1-\alpha/2} \cdot \frac{1}{n} \cdot \sqrt{b + c - \frac{(b - c)^2}{n}} \right] \ \]
만약 이항분포로 얻어진 값을 정규분포에 맞도록 조정하기 위해 continuity correction을 진행한다면 wald에서는 확률 \(p\) 하나당 \(\frac{1}{2n}\)만큼 분산을 늘려, 아래와 같은 값을 갖게 된다.
\[ CI_{1-\alpha/2}(\hat{\theta}) = \left[ \hat{\theta} \pm \left( z_{1-\alpha/2} \cdot \frac{1}{n} \sqrt{b + c - \frac{(b - c)^2}{n} }+ \frac{1}{n} \right) \right]. \]
3. Wald 이외 신뢰구간 계산 방법
3.1. Agresti 신뢰 구간(Agresti CI)
Wald 신뢰 구간을 수정한 방법 중 하나로, 샘플 크기에 특정 상수를 더한 후 위 wald와 같은 방식으로 전체 샘플 크기 𝑛에 1에서 4까지의 상수를 더한 후 이 방법들을 비교했습니다. 그 결과, n+2를 사용하는 것이 표준 Wald 신뢰 구간과 비교했을 때 표본 증가의 효과로 신뢰 구간의 성능(coverage probability)을 향상시킴을 확인했습니다. \[ \ CI_{1-\alpha/2}(\hat{\theta}) = \left[ \hat{\theta} \pm z_{1-\alpha/2} \cdot \frac{1}{n+2} \cdot \sqrt{(b+0.5) + (c+0.5) - \frac{(b - c)^2}{n+2}} \right] \ \]
3.2. Tango
Tango는 귀무가설과 함수는 두 비율 간의 차이와, 주어진 정보 행렬(\(b\), \(c\), \(n\))을 기반으로 likelihood를 계산해 신뢰구간을 계산하는 방법입니다. R에서는 scoreci.mp라는 별도의 함수를 사용해 진행합니다.
if (ci.method == "tango") {
# sensitivity
b <- tab$diseased[1,2]; c <- tab$diseased[2,1]; n <- tab$diseased[3,3]
tango <- scoreci.mp(b, c, n, conf.level=1-alpha)
sens.diff.se <- NA
sens.diff.cl <- sort(c(tango$conf.int[1], tango$conf.int[2]))
if ( (tango$conf.int[1] > sens.diff) | (tango$conf.int[2] < sens.diff))
sens.diff.cl <- sort(-1*sens.diff.cl)}4. 예시 코드
4.1. 일반적인 민감도, 특이도 검정
Loading required package: PropCIst1 <- read.tab.paired(18, 14, 0, 18,
18, 12, 2, 18)
t1Two binary diagnostic tests (paired design)
Test1: 'Noname 1'
Test2: 'Noname 2'
Diseased:
Test1 pos. Test1 neg. Total
Test2 pos. 18 14 32
Test2 neg. 0 18 18
Total 18 32 50
Non-diseased:
Test1 pos. Test1 neg. Total
Test2 pos. 18 12 30
Test2 neg. 2 18 20
Total 20 30 50sesp.diff.ci(t1, ci.method="wald", cont.corr=FALSE)$sensitivity
test1 test2 diff diff.se diff.lcl diff.ucl
0.36000000 0.64000000 0.28000000 0.06349803 0.15554615 0.40445385
$specificity
test1 test2 diff diff.se diff.lcl diff.ucl
0.60000000 0.40000000 -0.20000000 0.06928203 -0.33579029 -0.06420971
$ci.method
[1] "wald"
$alpha
[1] 0.05
$cont.corr
[1] FALSEsesp.diff.ci(t1, ci.method="wald", cont.corr=TRUE)$sensitivity
test1 test2 diff diff.se diff.lcl diff.ucl
0.36000000 0.64000000 0.28000000 0.08349803 0.11634687 0.44365313
$specificity
test1 test2 diff diff.se diff.lcl diff.ucl
0.60000000 0.40000000 -0.20000000 0.08928203 -0.37498957 -0.02501043
$ci.method
[1] "wald"
$alpha
[1] 0.05
$cont.corr
[1] TRUEsesp.diff.ci(t1, ci.method="agresti-min")$sensitivity
test1 test2 diff diff.se diff.lcl diff.ucl
0.36000000 0.64000000 0.28000000 0.06444681 0.15368658 0.40631342
$specificity
test1 test2 diff diff.se diff.lcl diff.ucl
0.60000000 0.40000000 -0.20000000 0.06954236 -0.33630053 -0.06369947
$ci.method
[1] "agresti-min"
$alpha
[1] 0.05
$cont.corr
[1] FALSEsesp.diff.ci(t1, ci.method="tango")$sensitivity
test1 test2 diff diff.se diff.lcl diff.ucl
0.3600000 0.6400000 0.2800000 NA 0.1747417 0.4166512
$specificity
test1 test2 diff diff.se diff.lcl diff.ucl
0.60000000 0.40000000 -0.20000000 NA -0.34470882 -0.06111243
$ci.method
[1] "tango"
$alpha
[1] 0.05
$cont.corr
[1] FALSE아래와 같은 방법으로 answer, test1, test2를 열로 가지고, 각 결과 데이터가 1,0(1은 diseased sample, 2는 non-diseased sample)로 표시되어 있는 경우에는 바로 paired table을 제작할 수 있습니다.
4.2. 민감도, 특이도의 비열등성 검정
어떤 테스트의 민감도와 특이도가 대조 테스트에 비해 떨어지지 않는다는 것을 확인하기 위해서는 비열등성 검정을 사용해야 합니다. 이 때는 sesp.diff.ci 함수를 통해 구한 standard error 값과 두 민감도/특이도의 차이를 통해 p value를 계산합니다. 아래 예시에서는 민감도에 대한 비열등성 마진(sens_margin)을 5%, 특이도에 대한 비열등성 마진(spec_margin)을 10%로 설정했습니다.
민감도를 예로 들자면 귀무가설은 아래와 같고, 아래 수식에 따라 p value를 계산합니다.
\[H_0 \colon \text{test1 sensitivity} - \text{test2 sensitivity} \leq -\text{ sensitivity margin}\] \[p-value = 1 - Φ\left(\frac{\text{test1 sensitivity} - \text{test2 sensitivity} + \text{sensitivity margin}}{\text{sensitivity diff.SE}}\right)\]
library(DTComPair)
t1 <- read.tab.paired(18, 14, 0, 18,
18, 12, 2, 18)
t1.wald <- sesp.diff.ci(t1, ci.method="wald", cont.corr=FALSE)
sens_margin <- 0.05
spec_margin <- 0.1
p_value_sensitivity <- pnorm((t1.wald$sensitivity['diff'] + sens_margin) / t1.wald$sensitivity['diff.se'], lower.tail = FALSE)
p_value_specificity <- pnorm((t1.wald$specificity['diff'] + spec_margin) / t1.wald$specificity['diff.se'], lower.tail = FALSE)
p_value_sensitivity diff
1.012589e-07 p_value_specificity diff
0.9255427 따라서 민감도에 대해서는 위 귀무가설이 기각 되었음으로 test1의 민감도가 test2보다 열등하지 않다고 말할 수 있지만, 특이도에 대해서는 귀무가설이 기각되지 않았기 때문에 test1의 특이도가 test2에 비해 열등하지 않다는 결론을 내릴 수 없습니다.
Reuse
Citation
@online{myung2024,
author = {Myung, Hyojong},
title = {Wald {Confidence} {Interval} for a {Difference} of {Binomial}
{Proportions} {Based} on {Paired} {Data}},
date = {2024-08-26},
url = {https://blog.zarathu.com/posts/2024-08-26-sensspec-waldinterval/},
langid = {en}
}