Определение чувствительности и специфичности комбинации маркеров Опции

[Zaycho]

Здравствуйте!

Есть задача определить диагностическую ценность ряда маркеров (допустим, А, В и С) и их комбинаций в диагностике некоторого заболевания.
Я расчитал все параметры (чувствительность, специфичность, PPV, NPV, LR) и построил в SPSS ROC-кривые отдельно для А, В и С.
Как теперь сделать то же самое для комбинаций этих маркеров, допустим, А+В и В+С? Где копать?

Заранее огромное спасибо!

[DrgLena]

А если не лес строить, а одно дерево, то при одном ветвлении, на двух конечных узлах при условии С<=16,64 классификация произойдет с точностью 65,7%
Переменная С разделяет пациентов реальной выборки, а вот В не разделяет и вообще не имеет различий в группах выживших и нет ни по каким парным критериям. И не удивительно, такие разбросы данных (скорее всего ошибки ввода) например, у умерших больных есть значение 237 и 0,14. Вряд ли В может называться маркером. Но в сгенерированных выборках значения могут повторяться или вообще не попасть.
ПО приведенным исходным данным точки разделения, и соответственно чувствительность и специфичность с моими расчетами не сходятся, но площадь ROC сходится, это очень странно. Автор должен верифицировать данные, иначе никакие многомерные методы на таких данных применять нельзя.
Автор поста проделал часть работы и задал конкретный вопрос, как объединить результаты трех тестов. Это не задача классификации. Есть понятие положительного и отрицательного теста (найдены значения для разделения на положительный и отрицательный результат теста) и нужно узнать вероятность исхода при различных сочетаниях положительных или отрицательных тестов.
Если не придираться к конкретным данным, то решить эту задачу можно следующим образом. Создать три бинарных переменных по результатам ROC, для этого использовать оптимальные точки разделения, фактор риска имеется, код 1 и нет ? 0. Эти переменные использовать в логистической регрессии. По полученным коэффициентам посчитать вероятности наступления интересующего события. При этом возможны 8 комбинаций наличия и отсутствия трех факторов риска, от отсутствия всех трех (000) до присутствия всех трех (111). Каждая из 8 групп риска будет иметь свою расчетную вероятность наступления события. При этом экспоненциальные коэффициенты уравнения регрессии для каждого предиктора будут иметь ясную интерпретацию и показывать независимый вклад каждого фактора.

[p2004r]

А если не лес строить, а одно дерево, то при одном ветвлении, на двух конечных узлах при условии С<=16,64 классификация произойдет с точностью 65,7%
Переменная С разделяет пациентов реальной выборки, а вот В не разделяет и вообще не имеет различий в группах выживших и нет ни по каким парным критериям. И не удивительно, такие разбросы данных (скорее всего ошибки ввода) например, у умерших больных есть значение 237 и 0,14. Вряд ли В может называться маркером. Но в сгенерированных выборках значения могут повторяться или вообще не попасть.
ПО приведенным исходным данным точки разделения, и соответственно чувствительность и специфичность с моими расчетами не сходятся, но площадь ROC сходится, это очень странно. Автор должен верифицировать данные, иначе никакие многомерные методы на таких данных применять нельзя.

1 по мимо того что B надо логарифмировать, имеются корреляции между предикторами

Код

> cor(cbind(data.na[,c(2,4)], log(data.na[,3])))                          A         C log(data.na[, 3])
A                 1.0000000 0.2813438         0.1289035
C                 0.2813438 1.0000000         0.2625643
log(data.na[, 3]) 0.1289035 0.2625643         1.0000000

> cor.test(data.na[,4], data.na[,2])

    Pearson's product-moment correlation

data:  data.na[, 4] and data.na[, 2]
t = 2.4177, df = 68, p-value = 0.01831
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.04965218 0.48430137
sample estimates:
      cor
0.2813438

> cor.test(data.na[,4], log(data.na[,3]))

    Pearson's product-moment correlation

data:  data.na[, 4] and log(data.na[, 3])
t = 2.2439, df = 68, p-value = 0.0281
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.02940446 0.46862610
sample estimates:
      cor
0.2625643

> cor.test(data.na[,2], log(data.na[,3]))

    Pearson's product-moment correlation

data:  data.na[, 2] and log(data.na[, 3])
t = 1.0719, df = 68, p-value = 0.2876
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
-0.1093837  0.3531800
sample estimates:
      cor
0.1289035

Группа выделена и имеет отличный состав

Код

> table(factor(cmdscale(dist(rf1$proximity), k=2)[,2]>0.2), factor(data.na$выжил.умер))
      
         0  1
  FALSE 24 21
  TRUE   9 16

во 2й группе исходы 9 к 16, в первой 24 к 21

вклад в разделение групп оказывает не один единственный предиктор, (B и С вместе)

Кружки первая группа, треугольники вторая. (На сечении С|log(B) вообще летальные исходы сгруппировались в центре и вверху. Может отработать эту нелинейность не может при разделении групп рандомфорест. Нужен какой нибудь svm или явная трансформация данных.)

Код

> pairs(cbind(data.na[,c(2,4)],
              log(data.na[,3])),
        pch=c(1,2)[factor(cmdscale(dist(rf1$proximity), k=2)[,2]>0.2)],
        col=c("red","green")[factor(data.na$выжил.умер)])

второй вариант наоборот цветом группы, значками исходы.

Эскизы прикрепленных изображений

 

[p2004r]

собственно да, svm получше выделяет

Код

> model <- svm( выжил.умер ~ A+B+C, data = data.na, type= "C-classification")
> table(model$decision.values<=0,data.na$выжил.умер)
      
         0  1
  FALSE 25 12
  TRUE   8 25
> table(model$decision.values<=-0.1,data.na$выжил.умер)
      
         0  1
  FALSE 27 13
  TRUE   6 24

## логарифмирование имеет значение!!
> model <- svm( выжил.умер ~ A+log(B)+C, data = data.na, type= "C-classification")
> table(model$decision.values<=-0.1,data.na$выжил.умер)
      
         0  1
  FALSE 27  9
  TRUE   6 28

вот собственно он их сумел обвести

Код

> pairs(cbind(data.na[,c(2,4)],
              log(data.na[,3])),
        pch=c(1,2)[factor(model$decision.values<=-0.1)],
        col=c("red","green")[factor(data.na$выжил.умер)])

более сложное ядро в svm смогло бы осилить две области выделения. а так наверное надо просто отрезать всех явно тяжелых сначала, и потом запускать svm который научен центр резать.

PS пробабилити

Код

> cbind(model$decision.values,data.na$выжил.умер)
           0/1  
1   0.98331246 0
2   0.61530576 1
3   1.00012653 0
4  -0.21684733 0
5   0.99955666 0
6   0.16078293 1
7  -1.00012726 1
8   1.14854667 1
9  -0.99985686 1
10  0.99966273 0
11  0.37811021 0
12 -0.42755238 1
13  1.00021827 0
14  1.00007588 0
15  0.76925049 0
16  0.88688175 0
17  1.01823453 0
18  0.86095764 0
19  1.01908459 0
20  1.18428922 1
21  0.45588013 1
22 -0.73663584 1
23  0.89897670 0
24  1.00026140 0
25  1.00020960 0
26  1.12620025 0
27 -0.71744634 1
28 -0.93837583 1
29 -0.99987350 1
31 -1.20454616 0
33  0.62211341 0
34 -1.37866569 1
35 -1.07486185 1
36 -1.09031422 1
37 -1.06053528 1
38 -0.50515338 1
39 -0.75995557 1
40  0.16392607 0
41 -1.00022084 1
42  0.34876086 0
43 -0.64127520 1
44 -1.15716373 1
45  0.11486255 0
46 -0.13282404 0
47 -1.04569130 0
48  0.37514348 0
49  0.56057748 0
50  0.49439228 1
51 -1.11937289 1
52 -0.89372931 1
53 -0.49529302 0
54 -0.65597539 1
55 -0.92934884 1
56  0.62306627 0
57 -1.15534228 1
58 -1.02889926 1
59 -1.00000111 1
60  0.11401070 0
61  0.68242740 1
62 -1.00003176 1
63 -0.65191720 1
64  0.74228345 1
65 -0.04860989 1
66 -1.13318704 1
67 -0.73920762 1
68  0.73537165 0
69 -1.04658105 0
70  0.16253365 0
71 -0.44523330 1
72  0.41346047 0

Сообщение отредактировал p2004r - 4.08.2012 - 11:29

Эскизы прикрепленных изображений

 

[p2004r]

линеаризирующее преобразование (!после отрезания всех явно тяжелых!):

1) берем точку среднее log(B) -- среднее C (или лучше среднее по шкалам mds) для летальных исходов

2) строим вектора с началом в этой точке и концом в точках наблюдений

3) новое пространство --- угол вектора vs его длинна

в таком пространстве должны работать линейные методы

Код

> data.leit<- data.na[data.na$C<25,]
> nrow(data.leit)
[1] 62
> model1 <- svm( выжил.умер ~ A+log(B)+C, data = data.leit, type= "C-classification")
> table(model1$decision.values<=0,data.leit$выжил.умер)
      
         0  1
  FALSE 19  5
  TRUE   6 32

вот по группе без явно безнадежных (критерий C>25) svm вырезает выживших до 32.

вот так решение графически выглядит

Код

> plot(cmdscale(dist(cbind(data.leit[,c(2,4)],
                                 log(data.leit[,3]))),
                      k=2),
             col=c("red","green")[factor(data.leit$выжил.умер)],
             pch=c(1,2)[factor(model1$decision.values<=0)])

Трансформация

Код

> tapply(data.leit[,3], data.leit$выжил.умер, mean)
       0        1
12.67121 20.19569
> tapply(data.leit[,4], data.leit$выжил.умер, mean)
       0        1
13.48000 11.40541
> tapply(log(data.leit[,3]), data.leit$выжил.умер, mean)
       0        1
1.775605 1.189165
> ((data.leit$C-12.67121)^2+(log(data.leit$B)-1.775605)^2)^0.5
[1]  5.3311377  5.1087745  2.1709325  6.0049577  6.3416986  2.5971127
[7]  7.6488885  9.4283268  2.9613371  7.3290071  2.0874944  2.6590571
[13]  0.7077666  1.7968509  1.6400452  4.3481172  1.3288040  5.5428207
[19]  2.6856595  9.7595860  3.7072333  4.6724333  4.3954970  7.6958694
[25]  7.8144577  5.1120372  3.4003874  6.7942236  0.7232215  2.2810260
[31]  0.6200533  3.0480755  5.3432953  1.8203928 12.6712147  7.3712229
[37]  2.4614863  6.8202322  4.3332743  0.6811193 10.3318129  7.2930514
[43]  9.8974529  4.0868944  2.6899526  6.9645308  2.7174036  3.7413277
[49]  2.2007804 12.0539703  3.0041289 11.3308822  8.6715569  3.7608580
[55]  7.3350855  3.3930430  4.4240841  7.5608856  6.8246717  4.3312637
[61] 10.4496927  6.0334436
> l<-((data.leit$C-12.67121)^2+(log(data.leit$B)-1.775605)^2)^0.5
> atan2(data.leit$C-12.67121, log(data.leit$B)-1.775605)
[1]  1.60047511  1.01092706 -0.87854661 -1.90576287  1.63461177  2.02937132
[7] -2.08195084  1.42536019 -2.91294791  1.57849395 -2.21332718  2.61831456
[13] -1.24799125  0.18401803  0.20184394  1.66511774  1.56620199  1.29199482
[19] -1.67457571 -1.43611893  0.08880546 -1.54791341 -0.65318030 -1.65087067
[25] -1.37902877 -1.98913076 -2.23799124 -1.76137710 -1.95235967 -2.31937471
[31]  2.58268440  2.69050581  1.49709553 -1.16313148 -1.57165444  1.67814707
[37]  1.90064623 -1.78022440  1.52529843  0.50375801  1.54660568 -1.98674087
[43] -1.78502201 -2.02573395 -1.68891279 -1.27953913  1.02940171 -2.34643526
[49]  2.99163425 -1.82347545 -1.09554895  1.55157912 -1.56185158 -1.78957673
[55]  1.52936211  1.76571555 -2.16280092  1.81921224 -1.78326322  1.60459531
[61]  1.72306160  1.08266063
>  a<-atan2(data.leit$C-12.67121, log(data.leit$B)-1.775605)
> plot(a,
       l,
       col=c("red","green")[factor(data.leit$выжил.умер)],
       pch=c(1,2)[factor(model1$decision.values<=0)])

и в пространстве трансформированных шкал mds

Код

> data.mds<-cmdscale(dist(cbind(data.leit[,c(2,4)],log(data.leit[,3]))))
> tapply(data.mds[,1], data.leit$выжил.умер, mean)
         0          1
-1.2605935  0.8517523
> tapply(data.mds[,2], data.leit$выжил.умер, mean)
         0          1
-0.2579902  0.1743177
> l.mds<-((data.mds[,1]+1.2605935)^2+(data.mds[,2]+0.2579902)^2)^0.5
> a.mds<-atan2(data.leit$C-12.67121, log(data.leit$B)-1.775605)
> plot(a.mds,
       l.mds,
       col=c("red","green")[factor(data.leit$выжил.умер)],
       pch=c(1,2)[factor(model1$decision.values<=0)])

похожую трансформацию, судя по внешнему сходству, рандом_форест и пытался сделать.

PS забыл отнормировать данные прежде чем считать длину вектора и угол

Код

> data.mds.scale<-scale(data.mds, center=c(-1.2605935, -0.2579902))
> l.mds.scale<-((data.mds.scale[,1])^2+(data.mds.scale[,2])^2)^0.5
> a.mds.scale<-atan2(data.mds.scale[,1], data.mds.scale[,2])
> plot(a.mds.scale,l.mds.scale, col=c("red","green")[factor(data.leit$выжил.умер)], , pch=c(1,2)[factor(model1$decision.values<=0)])

на картинке я зону svm приблизительно обвел, горизонтальная ось замкнута (это радианы угла вектора)

поскольку все читать утомительно я резюмирую:

решатель состоящий из правила -- все кто тяжелее 25 считать безнадежными и svm-решателя для остальных случаев дает
11 ошибок на 70 ~ 0.16

Сообщение отредактировал p2004r - 5.08.2012 - 12:47

Эскизы прикрепленных изображений

 

[Zaycho]

Коллеги, а если вернуться к посту #14, то как определить оптимальные точки разделения (cut-off) каждого из маркеров в комбинации?