2. AUC 구하는 방법.
요약: 진단 방법의 효율성을 판단하는 방법 중 널리 사용되는 것이 ROC curve 이다. 민감도(sensitivity)와 특이도(specificity)가 어떤 관계를 갖고 변하는지를 이차원 평면 상에 표현한 것이 ROC curve인데, ROC curve 아래의 면적(AUC, area under curve) 이 넓을수록 좋은 진단 방법이라 할 수 있다. 이 글은 ROC curve의 AUC를 구하는 간단한 방법1을 설명한다.
1. ROC Curve란.
ROC Curve는 Receiver-Operating Characteristic curve의 줄임말로, 특정 진단 방법의 민감도와 특이도가 어떤 관계를 갖고 있는지를 표현한 그래프이다.
진단 방법에 의한 환자의 분류는, 점수에 따라 사람을 환자/정상 으로 분류하는 것으로 결국 분류(classification)의 한 예이고, 따라서 아래 나오는 진단 방법은 classifier로 바꾸어 이해해도 큰 무리가 없겠다.
우선 민감도와 특이도가 무엇인지 간단히 알아 보자. 보다 자세한 내용은 별도의 글을 참고하자. 민감도는 진짜 환자 중 검사 방법이 환자를 얼마나 잘 골라 내는가, 하는 것이고, 특이도는 진짜 정상 중 검사 방법이 정상을 얼마나 잘 골라 내는가를 의미한다. 예를 들어 보자. 내가 혈압으로 심근경색을 진단할 수 있다고 가정해 보자. '혈압이 얼마 이상이면 심근경색이다'라는 진단은 민감도 100%, 특이도 100%를 갖게 혈압의 기준을 설정할 수 있다. 가령, '혈압이 5 이상이면 심근경색' 이러면 전체 인구가 모두 환자로 판별되어, 누구나 환자로 판정이 나고, 따라서 환자를 100% 환자라고 진단할 수 있다. 문제는, 환자가 아닌 사람도 몽땅 환자라고 하기 때문에
위의 예에서, 혈압을 얼마로 하는 것이 민감도/특이도를 높이게 될까? 어차피 민감도를 높이면 특이도는 줄어 들고, 특이도를 높이면 민감도는 감소하기 때문에 (즉, 하나가 다른 하나에 영향을 받기 때문에), 혈압에 대한 기준을 변경시켜 가면서 각각의 값을 그래프에 그려 보는 것이 좋겠다. (혈압, 심근경색여부)에 대한 데이터가 다음과 같다고 하자.
80/0
92/0
93/1
98/0
102/1
110/1
112/0
119/1
이 때 심근경색 여부가 1 이면 심근경색이 있는 것이다. 위의 경우 '혈압이 N이상이면 심근경색이다'라 한다면, N 을 80 으로 잡으면 민감도가 1 이 되나, 특이도가 0 이 된다. 따라서 N 을 변화시키면서 각각의 값을 구해 보면 다음과 같다.
| 혈압 | class | TP | TN | FP | FN | P(=TP+FN) | N(=TN+FP) | 특이도 | 1-특이도 | 민감도 |
| 80 | 0 | 4 | 0 | 4 | 0 | 4 | 4 | 0 | 1 | 1 |
| 92 | 0 | 4 | 1 | 3 | 0 | 4 | 4 | 0.25 | 0.75 | 1 |
| 93 | 1 | 4 | 2 | 2 | 0 | 4 | 4 | 0.5 | 0.5 | 1 |
| 98 | 0 | 3 | 2 | 2 | 1 | 4 | 4 | 0.5 | 0.5 | 0.75 |
| 102 | 1 | 3 | 3 | 1 | 1 | 4 | 4 | 0.75 | 0.25 | 0.75 |
| 110 | 1 | 2 | 3 | 1 | 2 | 4 | 4 | 0.75 | 0.25 | 0.5 |
| 112 | 0 | 1 | 3 | 1 | 3 | 4 | 4 | 0.75 | 0.25 | 0.25 |
| 119 | 1 | 1 | 4 | 0 | 3 | 4 | 4 | 1 | 0 | 0.25 |
위의 테이블을 그래프로 그려 보면 다음과 같다.
x
축을 1-특이도로 한 이유는 저렇게 해야 우리가 익히 보아 오던 그래프 형태와 비슷하기 때문이다. 위 표와 그래프에서 알 수
있듯이, 진단 기준을 '혈압 80 이상이면 모두 환자'라 하면 민감도는 1 이지만 특이도가 0 이 되서 그래프 상에는 제일
오른쪽 상위 점으로 찍히게 된다.
위와 같은 그래프의 아래 면적이 1/2에서 멀먼 멀수록 좋은 진단 기준이라 할
수 있겠다. 즉 (민감도,1-특이도) 가 딱 두 점 (0,0), (1,1)에 찍히게 하는 진단 방법은 정확히 AUC가 1 이어서
가장 완벽한 방법이라 할 수 있다.
2. AUC 구하는 방법.
진단값을 오름차순으로 정렬한 후, 1-특이도가 변하는 지점에서의 sensitivity 합계를 전체 정상(negative) 수로 나누어 주면 된다. 이 때 낮은 진단값이 정상이라 가정한다. 반대일 경우 그에 맞게 적당히 변경해서 사용하면 된다. 작은 진단값에서부터 시작하여 큰 값으로 이동하면서 negative (즉, 이 경우 정상(0)) 을 만나는 점에서의 민감도값만 남겨 준다.
| 혈압 | class | TP | TN | FP | FN | P(=TP+FN) | N(=TN+FP) | 특이도 | 1-특이도 | 민감도 |
| 80 | 0 | 4 | 0 | 4 | 0 | 4 | 4 | 0 | 1 | 1 |
| 92 | 0 | 4 | 1 | 3 | 0 | 4 | 4 | 0.25 | 0.75 | 1 |
| 93 | 1 | 4 | 2 | 2 | 0 | 4 | 4 | 0.5 | 0.5 | 1 |
| 98 | 0 | 3 | 2 | 2 | 1 | 4 | 4 | 0.5 | 0.5 | 0.75 |
| 102 | 1 | 3 | 3 | 1 | 1 | 4 | 4 | 0.75 | 0.25 | 0.75 |
| 110 | 1 | 2 | 3 | 1 | 2 | 4 | 4 | 0.75 | 0.25 | 0.5 |
| 112 | 0 | 1 | 3 | 1 | 3 | 4 | 4 | 0.75 | 0.25 | 0.25 |
| 119 | 1 | 1 | 4 | 0 | 3 | 4 | 4 | 1 | 0 | 0.25 |
위의 경우 민감도 중 노랑색으로 칠한 값들의 합인 3 이 되겠다. 이 값을 전체 negative 수 (정상 수, 이 경우 4명) 로 나누어 준다. 그러면 AUC가 0.75 가 나온다.
여기까지가 기본적인 방법이고, 실제로 구현을 할 때는 계산의 용이함과 약간(?)의 속도를 위하여 살짝 계산을 줄여서 다음과 같이 계산한다.
where
Pi = the number of passed positive sample at ith negative point,
P = total positive number,
N = total negative number.
실제 C++ 코드로 보면 다음과 같다. 중요: 기존의 코드는 값에 같은 값(80, 80 처럼)이 있는 경우 오류가 있음이 밝혀 졌다.
bool sort_by_class(const std::pair<double,bool>& op1, const std::pair<double,bool>& op2) {
return op1.first < op2.first;
}
// yul code.
double get_AUC(std::vector<std::pair<double,double> >& xy_pair) {
double AUC = 0;
std::sort(xy_pair.begin(), xy_pair.end(), sort_by_class);
std::vector<std::pair<double, double> >::const_iterator roc_pos = xy_pair.begin();
for (roc_pos; roc_pos != xy_pair.end(); roc_pos++){
if (roc_pos + 1 == xy_pair.end()){
continue;
}
if (roc_pos->first < (roc_pos+1)->first){
double temp_area = 0;
double del_x = (roc_pos+1)->first - roc_pos->first;
double del_y = (roc_pos+1)->second - roc_pos->second;
temp_area = (del_x * roc_pos->second) + (del_x * del_y * 0.5);
AUC = AUC + temp_area;
}
else continue;
}
return AUC;
}
double get_ROC_AUC(std::vector<std::pair<double,bool> >& value)
{
std::sort(value.begin(), value.end(), sort_by_class);
std::vector<std::pair<double, double> > roc_pair;
roc_pair.push_back(std::pair<double, double>(0, 0));
int negative_total = 0;
int positive_total = 0;
std::vector<std::pair<double,bool> >::const_iterator pos;
for (pos = value.begin(); pos != value.end(); pos++){
if(pos->second == true) positive_total++;
if(pos->second == false) negative_total++;
}
for (pos = value.begin(); pos != value.end(); pos++){
double classifier = -1;
if (pos + 1 == value.end()) classifier = (pos->first) + 1;
else classifier = (pos + 1)->first;
if (classifier - (pos->first) < 0.00000001) continue;
double positive_count = 0;
double negative_count = 0;
std::vector<std::pair<double, bool> >::const_iterator snsp_pos = value.begin();
for (snsp_pos; snsp_pos <= pos; snsp_pos++){
if (snsp_pos->second == true) positive_count++;
else if (snsp_pos->second == false) negative_count++;
}
roc_pair.push_back(std::pair<double, double>((negative_count/static_cast<double>(negative_total)), (positive_count/static_cast<double>(positive_total))));
}
double AUC = get_AUC(roc_pair);
return AUC;
}
(위 코드는 연구실 후배가 작성한 코드임)
문제가 있던 기존 코드는 다음에 있다.
위에서는 일부러 처리하지 않았는데, AUC는 0.5 이하가 될 수 없다. 그런데 만약 위의 계산 결과가 0.5 보다 작다면 질병과 정상군을 반대로 입력한 것이다. 따라서 반환값을 보고 잘 판단하면 된다.
위의 코드를 갖고 실제로 구한 값은 0.823508 이 나오고 GraphPad Prism 으로 구한면 0.8235 가 나온다. 입력 파일은 다음에 있다.
data.txt이에 대한 ROC 를 그려 보면 다음과 같다.
- 원래 복잡한 것을 간단히 설명한다기보다는, 방법 자체가 간단하다는 뜻. ㅋ [본문으로]
'연구관련 > Bioinfo류' 카테고리의 다른 글
| 정규화(normalization) (18) | 2010.04.19 |
|---|---|
| drug 관련 싸이트 (0) | 2010.04.06 |
| 생물정보학(bioinformatics)에 대하여 (3) | 2010.01.12 |
| outlier 빼고 상관계수 구하기 : Mahalanobis 거리 (1) | 2010.01.04 |
| cytoscape: network 연구 프로그램 (0) | 2009.03.06 |
