Как сравнить результаты регрессии Кокса Опции

[propedevt]

Уважаемые форумчане!
Для анализа выживаемости строю модель Кокса, с разными независимыми переменными. И появилась необходимость сравнить модели.
И собственно вопрос у меня к Вам - как сравнить что одна полученная модель лучше чем другая?

[p2004r]

Теперь basehazard и счет на прямую. Центрированная версия равна случаю, когда h0(t) hazard считается для случая когда z<-0. Первый вариант считает когда значения предикторов заменены их средними.

Код

> basehaz(model5, centered=FALSE)
         hazard time
1  0.0001002811   16
2  0.0002024633   45
3  0.0003061236   48
4  0.0004101635   55
5  0.0006646868   60
6  0.0007949146   62
7  0.0009255306   69
8  0.0010569939   75
9  0.0011890672   76
10 0.0013237703   83
11 0.0017833064   90
12 0.0019668331   91
13 0.0021839121   92
14 0.0024035279   99
15 0.0026280378  100
16 0.0033604084  120
17 0.0036432445  126
18 0.0039314237  140
19 0.0048821851  150
20 0.0056043508  160
21 0.0067440470  165
22 0.0092375477  180
23 0.0097008584  210
24 0.0106394178  220
25 0.0111215067  240
26 0.0116072477  250
27 0.0121144210  254
28 0.0126251301  267
29 0.0131371764  270
30 0.0136700471  278
31 0.0142122284  280
32 0.0148007218  300
33 0.0153900052  307
34 0.0159800228  356
35 0.0172082186  380
36 0.0186000766  390
37 0.0193368429  400
38 0.0201008041  420
39 0.0209019571  430
40 0.0217483690  440
41 0.0226347604  450
42 0.0236140231  475
43 0.0246231476  480
44 0.0256503391  488
45 0.0289022541  500
46 0.0300917523  505
47 0.0313090582  540
48 0.0326040644  550
49 0.0353892636  560
50 0.0368104561  570
51 0.0382620189  580
52 0.0397520960  590
53 0.0412467805  605
54 0.0427499492  620
55 0.0442892539  650
56 0.0474533644  660
57 0.0507720453  680
58 0.0524983401  690
59 0.0560595003  700
60 0.0579564065  705
61 0.0618094617  710
62 0.0658048238  720
63 0.0678481443  735
64 0.0699531004  740
65 0.0721146898  745
66 0.0766243540  750
67 0.0837802983  760
68 0.0862338454  766
69 0.0941286515  770
70 0.1033435078  780
71 0.1067425188  785
72 0.1138292088  790
73 0.1212529131  800
74 0.1212529131  805
> basehaz(model5)
         hazard time
1  0.0007172164   16
2  0.0014480302   45
3  0.0021894155   48
4  0.0029335149   55
5  0.0047538818   60
6  0.0056852791   62
7  0.0066194528   69
8  0.0075596863   75
9  0.0085042832   76
10 0.0094676883   83
11 0.0127543187   90
12 0.0140669133   91
13 0.0156194759   92
14 0.0171901823   99
15 0.0187958908  100
16 0.0240338513  120
17 0.0260567126  126
18 0.0281177885  140
19 0.0349176933  150
20 0.0400826682  160
21 0.0482338464  165
22 0.0660675194  180
23 0.0693811469  210
24 0.0760937821  220
25 0.0795417127  240
26 0.0830157622  250
27 0.0866430977  254
28 0.0902957213  267
29 0.0939579080  270
30 0.0977690327  278
31 0.1016467478  280
32 0.1058556895  300
33 0.1100702806  307
34 0.1142901239  356
35 0.1230742570  380
36 0.1330289128  390
37 0.1382983114  400
38 0.1437622096  420
39 0.1494921065  430
40 0.1555456970  440
41 0.1618852240  450
42 0.1688889721  475
43 0.1761062938  480
44 0.1834528314  488
45 0.2067107316  500
46 0.2152180967  505
47 0.2239243447  540
48 0.2331863105  550
49 0.2531062296  560
50 0.2632706871  570
51 0.2736523556  580
52 0.2843094803  590
53 0.2949995571  605
54 0.3057503141  620
55 0.3167595173  650
56 0.3393894339  660
57 0.3631248475  680
58 0.3754714163  690
59 0.4009410571  700
60 0.4145078485  705
61 0.4420651411  710
62 0.4706402210  720
63 0.4852541769  735
64 0.5003089551  740
65 0.5157687777  745
66 0.5480221780  750
67 0.5992019402  760
68 0.6167498625  766
69 0.6732140103  770
70 0.7391192393  780
71 0.7634291793  785
72 0.8141136296  790
73 0.8672084276  800
74 0.8672084276  805
> -log((survfit(model5))$surv)
[1] 0.0007172164 0.0014480302 0.0021894155 0.0029335149 0.0047538818
[6] 0.0056852791 0.0066194528 0.0075596863 0.0085042832 0.0094676883
[11] 0.0127543187 0.0140669133 0.0156194759 0.0171901823 0.0187958908
[16] 0.0240338513 0.0260567126 0.0281177885 0.0349176933 0.0400826682
[21] 0.0482338464 0.0660675194 0.0693811469 0.0760937821 0.0795417127
[26] 0.0830157622 0.0866430977 0.0902957213 0.0939579080 0.0977690327
[31] 0.1016467478 0.1058556895 0.1100702806 0.1142901239 0.1230742570
[36] 0.1330289128 0.1382983114 0.1437622096 0.1494921065 0.1555456970
[41] 0.1618852240 0.1688889721 0.1761062938 0.1834528314 0.2067107316
[46] 0.2152180967 0.2239243447 0.2331863105 0.2531062296 0.2632706871
[51] 0.2736523556 0.2843094803 0.2949995571 0.3057503141 0.3167595173
[56] 0.3393894339 0.3631248475 0.3754714163 0.4009410571 0.4145078485
[61] 0.4420651411 0.4706402210 0.4852541769 0.5003089551 0.5157687777
[66] 0.5480221780 0.5992019402 0.6167498625 0.6732140103 0.7391192393
[71] 0.7634291793 0.8141136296 0.8672084276 0.8672084276

hazard у одиночного случая по времени
h(t)=h0(t)*exp(sum(z_случая*coef(fit)))

считаем формулу по шагам для первого случая выборки

Код

> coef(model5)
        efmss           gal           cys
-3.7569726150  0.0944770843  0.0001896144
> data[1,c(7,4,6)]
  efmss gal  cys
1   0.3   6 1800
> data[1,c(7,4,6)]*coef(model5)
      efmss       gal       cys
1 -1.127092 0.5668625 0.3413059
> sum(data[1,c(7,4,6)]*coef(model5))
[1] -0.2189233
> exp(sum(data[1,c(7,4,6)]*coef(model5)))
[1] 0.8033833
>  -log((survfit(model5))$surv)*exp(sum(data[1,c(7,4,6)]*coef(model5)))
[1] 0.0005761997 0.0011633233 0.0017589398 0.0023567369 0.0038191892
[6] 0.0045674583 0.0053179578 0.0060733257 0.0068321990 0.0076061827
[11] 0.0102466066 0.0113011232 0.0125484260 0.0138103053 0.0151003048
[16] 0.0193083947 0.0209335277 0.0225893616 0.0280522916 0.0322017461
[21] 0.0387502666 0.0530775416 0.0557396545 0.0611324736 0.0639024834
[26] 0.0666934767 0.0696076175 0.0725420743 0.0754842140 0.0785460079
[31] 0.0816612994 0.0850426929 0.0884286250 0.0918187766 0.0988758024
[36] 0.1068732067 0.1111065535 0.1154961580 0.1200994615 0.1249628150
[41] 0.1300558851 0.1356825793 0.1414808550 0.1473829406 0.1660679492
[46] 0.1729026242 0.1798970785 0.1873379871 0.2033413173 0.2115072727
[51] 0.2198477318 0.2284094878 0.2369977169 0.2456346956 0.2544793055
[56] 0.2726598026 0.2917284374 0.3016474645 0.3221093486 0.3330086822
[61] 0.3551477507 0.3781044927 0.3898451007 0.4019398581 0.4143600213
[66] 0.4402718645 0.4813888306 0.4954865382 0.5408488915 0.5937960517
[71] 0.6133262514 0.6540452922 0.6967007661 0.6967007661

[p2004r]

Код

> exp(model5$linear.predictors)
  [1]   0.11232890   0.04297787   0.02450795   0.06653547   0.44724814
  [6]   0.10442370   0.04559316   0.08998441   0.08979395   1.14709744
[11]   0.04227987   0.08024934   0.17445771   0.55891127   0.09342941
[16]   0.09313293   0.08613578   0.09323165   0.35756428   0.25403493
[21]   0.14375222   0.11315116   0.11241344   0.03704815   0.12805265
[26]   0.05922301   0.14368160   0.06504510   0.44795261   0.08632247
[31]   0.07419565   0.35569425   0.08513071   0.12207266   0.04828226
[36]   0.07304558   0.07435988   0.06792582   0.05937370   0.10827094
[41]   0.21430708   1.28573887   0.02740159   0.10619491   1.38923493
[46]   0.26800790   0.27805913   0.27806554   0.07879766   0.14102941
[51]   0.39052389   0.08982980   0.15872279   0.17735826   0.12601166
[56]   0.08164121   6.89159980   0.43838287   0.46470958   0.89579887
[61]   0.36878196   1.52229755   0.71595552   0.57227080   0.49934867
[66]   2.00364842   0.29531268   7.31765007   1.90864101   0.09352414
[71]   1.28573887   2.53417668   1.43656366   0.24760979   0.41366511
[76]   1.46681930   1.80377799   1.04091624   0.05380667   0.17995576
[81]   0.13074055   0.18007876   1.33455272   0.23941348   1.31278772
[86]   0.15193422   0.40575200   0.52799318   0.61701292   2.68516850
[91]   1.58009249   2.05840242   0.41321756   0.16692950   0.11265189
[96]   5.85199422   0.13072850   1.74029800   0.89956176   0.53228300
[101]   1.49479436   0.51734466   4.92009247   1.21235927   0.42130538
[106]   0.79911269   0.59143095   0.66516120   1.58054137   1.14709744
[111]   0.27818298   0.28924135   0.81981082   0.42204997   0.84835882
[116]   2.18361003   2.48969131   0.15490086   1.94586396   0.61742043
[121]   1.84224835   2.22216563   1.16831089   0.25570773   1.06363906
[126]   0.11017759   0.31847704   1.73906924 236.80047306   1.04552952
[131]   0.71491192   3.31452879   0.54903263  15.56796569 117.75340276
[136]   6.89922294  25.94175269   5.61212163   1.92017850   3.86127465
[141]   9.33199040   3.19081830  13.87818390   6.75546260   6.88705121
[146]   2.02107380  27.90562252   7.44022727   9.32848045   3.06868455
[151] 104.88115464  43.13367980   4.55163006 124.35055369  31.28870201
[156]  25.94055751  78.88807278  20.66781100   7.45074491  20.29440783
[161]  18.45010077   3.76962888  14.67530276  19.51116438   2.43700853
[166]  10.67116195  14.71094846   4.22513952   1.79979410   1.67202854
[171]  12.16406643   4.47296628   7.43251999   1.80377799   8.65323018
[176]   5.50514631   8.19937467   2.18052759   4.91292080   8.03215988
[181]   1.49201973   9.86753570   1.35859638   2.38734662   8.49630599
[186]  27.90433687  14.95978216   3.00707267   4.72588644  35.08515480
[191]   1.04520851  10.86800471  27.94507043   9.16578236   4.55163006
[196]   1.25886070   4.50591302
>

риск считает как экспоненту от линейного предиктора

[p2004r]

Линейный предиктор он считает в момент фита


names(coef) <- dimnames(x)[[2]]
lp <- c(x %*% coef) + offset - sum(coef*coxfit$means)

Код

# напоминаю что было вот так
195  0.636968961  0.45032478  0.428191680
196  0.411550604 -0.11653772 -0.064805777
197  0.561829509  0.82823312  0.115327909
attr(,"constant")
[1] 1.967401
# константа получается вот так
> model5$mean
[1]    0.2795431   25.2335025 3341.7766497
> model5$mean*model5$coefficients
    efmss       gal       cys
-1.050236  2.383988  0.633649
> sum(model5$mean*model5$coefficients)
[1] 1.967401

теперь линейный предиктор для 197го случая

Код

# так он в принципе и все сразу посчитать может, а не только один случай
>  c(as.numeric(data[197,c(7,4,6)])%*%as.numeric(model5$coefficients)) - sum(model5$mean*model5$coefficients)
[1] 1.505391

# вот оно же "по простому", для одиночного случая

>  sum(data[197,c(7,4,6)]*model5$coefficients) - sum(model5$mean*model5$coefficients)
[1] 1.505391

Все рассчитывается за вычетом средней по выборке.

[p2004r]

Конкретно в этой модели первая переменная похоже не выполняет условия пропорциональности риска.

[propedevt]

Конкретно в этой модели первая переменная похоже не выполняет условия пропорциональности риска.

Хм, при проверке cox.zph она немного больше 0,05 выдала

[p2004r]

Хм, при проверке cox.zph она немного больше 0,05 выдала

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