关于r语言:R语言生存分析数据分析可视化案例

64次阅读

共计 16529 个字符,预计需要花费 42 分钟才能阅读完成。

原文链接:http://tecdat.cn/?p=2858

指标

本文的目标是对如何在 R 中进行生存剖析进行简短而全面的评估。对于该主题的文献很宽泛,仅波及无限数量的(常见)问题 / 特色。
可用的 R 包数量反映了对该主题的钻研范畴。


R 包

能够应用各种 R 包来解决特定问题,并且还有代替性能来解决常见问题。以下是本次审查中用于读取,治理,剖析和显示数据的软件包。
运行以下行以装置和加载所需的包。

if (!require(pacman)) install.packages("pacman")pacman::p_load(tidyverse, survival)

数据

该评估将基于 orca 数据集,该数据集蕴含来自基于人群的回顾性队列设计的数据。
它包含 1985 年 1 月 1 日至 2005 年 12 月 31 日期间芬兰最北部省份诊断为口腔鳞状细胞癌(OSCC)的 338 名患者的一部分。患者的随访始于癌症诊断之日,并于 2008 年 12 月 31 日死亡,迁徙或随访截止日期完结。死亡起因分为两类:(1))OSCC 死亡;(2)其余起因造成的死亡。
数据集蕴含以下变量:
id= 序号,
sex= 性别,类别 1 =“女性”的因素,2 =“男性”,
age= 诊断癌症日期的年龄(年),
stage= 肿瘤的 TNM 分期(因子):1 =“I”,…,4 =“IV”,5 =“unkn”
time= 自诊断至死亡或审查的随访工夫(以年为单位),
event= 完结随访的事件(因子):1 = 活检,2 = 口腔癌死亡,3 = 其余起因造成的死亡。

 将数据从 URL 加载到 R 中。

head(orca)
  id    sex      age stage  time          event1  1   Male 65.42274  unkn 5.081          Alive2  2 Female 83.08783   III 0.419 Oral ca. death3  3   Male 52.59008    II 7.915    Other death4  4   Male 77.08630     I 2.480    Other death5  5   Male 80.33622    IV 2.500 Oral ca. death6  6 Female 82.58132    IV 0.167    Other death
summary(orca)
       id             sex           age         stage         time                   event     Min.   :  1.00   Female:152   Min.   :15.15   I   :50   Min.   : 0.085   Alive         :109   1st Qu.: 85.25   Male  :186   1st Qu.:53.24   II  :77   1st Qu.: 1.333   Oral ca. death:122   Median :169.50                Median :64.86   III :72   Median : 3.869   Other death   :107   Mean   :169.50                Mean   :63.51   IV  :68   Mean   : 5.662                        3rd Qu.:253.75                3rd Qu.:74.29   unkn:71   3rd Qu.: 8.417                        Max.   :338.00                Max.   :92.24             Max.   :23.258                       

生存数据分析

生存剖析侧重于事件数据的工夫,通常称为故障工夫,。在咱们的例子中,是诊断后的死亡工夫。ŤTŤ≥ 0T≥0ŤT

为了定义生效工夫随机变量,咱们须要:
1。工夫起源(诊断 OSCC),
2。时间尺度(诊断后的年数,年龄),
3。事件的定义。咱们将首先思考总死亡率(或全因死亡率)。

图 1:转换的框图。

         Alive Oral ca. death    Other death            109            122            107 
FALSE  TRUE   109   229 

以图形形式显示察看到的随访工夫对于生存数据的剖析十分有帮忙。

 

OSCC 死亡更有可能在诊断后晚期产生,而不是其余起因引起的死亡。审查的类型怎么样?

 'Surv' num [1:338, 1:2]  5.081+  0.419   7.915   2.480   2.500   0.167   5.925+  1.503  13.333   7.666+ ... - attr(*, "dimnames")=List of 2  ..$ : NULL  ..$ : chr [1:2] "time" "status" - attr(*, "type")= chr "right"

而后将创立的生存对象用作生存剖析的其余特定函数中的响应变量。


预计生存性能

非参数估计

咱们将首先介绍一类非参数估计(a。)。

Kaplan–Meier 

生存曲线基于每个独特死亡工夫的危险数量和事件数量的列表。包的 survfit() 性能 survival 应用不同的办法创立(预计)生存曲线。

Call: survfit(formula = Surv(time, all) ~ 1, data = orca)         n     events     *rmean *se(rmean)     median    0.95LCL    0.95UCL    338.000    229.000      8.060      0.465      5.418      4.331      6.916     * restricted mean with upper limit =  23.3 

print() 函数仅返回预计的生存曲线的摘要。

   time n.risk n.event n.censor      surv     std.err     upper     lower1 0.085    338       2        0 0.9940828 0.004196498 1.0000000 0.98594012 0.162    336       2        0 0.9881657 0.005952486 0.9997618 0.97670413 0.167    334       4        0 0.9763314 0.008468952 0.9926726 0.96025924 0.170    330       2        0 0.9704142 0.009497400 0.9886472 0.95251755 0.246    328       1        0 0.9674556 0.009976176 0.9865584 0.94872286 0.249    327       1        0 0.9644970 0.010435745 0.9844277 0.9449699

该包中 ggsurvplot() 的专用性能 survminer 提供了预计的生存曲线的信息性阐明。无关 ?ggsurvplot 不同可能性(参数)的阐明。

 

默认的 KM 图表显示了生存函数。有几种代替计划 / 性能可供使用 

可升降的或精算的估算器

生命办法在精算师和人口统计学中十分广泛。它特地实用于分组数据。

为了在理论示例中显示此办法,咱们首先须要创立聚合数据,行将后续分组并在每个层中计算危险,事件和审查的人数。

基于分组的数据,咱们预计会用存活曲线 lifetab()KMsurv包。

      nsubs nlost nrisk nevent   surv    pdf hazard se.surv se.pdf se.hazard0-1     338     0 338.0     64 1.0000 0.1893 0.2092  0.0000 0.0213    0.02601-2     274     4 272.0     41 0.8107 0.1222 0.1630  0.0213 0.0179    0.02542-3     229     9 224.5     21 0.6885 0.0644 0.0981  0.0252 0.0136    0.02143-4     199    12 193.0     20 0.6241 0.0647 0.1093  0.0265 0.0140    0.02444-5     167     9 162.5     13 0.5594 0.0448 0.0833  0.0274 0.0121    0.02315-6     145    14 138.0     13 0.5146 0.0485 0.0989  0.0279 0.0131    0.02746-7     118     5 115.5      8 0.4662 0.0323 0.0717  0.0283 0.0112    0.02547-8     105     8 101.0      9 0.4339 0.0387 0.0933  0.0286 0.0126    0.03118-9      88     7  84.5      1 0.3952 0.0047 0.0119  0.0288 0.0047    0.01199-10     80     4  78.0      8 0.3905 0.0401 0.1081  0.0288 0.0137    0.038210-11    68     4  66.0      5 0.3505 0.0266 0.0787  0.0291 0.0116    0.0352

Nelson-Aalen 预计

图形比拟

能够绘制不同的生存函数估计值来评估潜在的差别。

能够从预计的存活曲线导出诸如分位数的集中趋势的度量。

      q km.quantile km.lower km.upper fh.quantile fh.lower fh.upper25 0.25       1.333    1.084    1.834       1.333    1.084    1.74750 0.50       5.418    4.331    6.916       5.418    4.244    6.91375 0.75      13.673   11.748   16.580      13.673   11.748   15.833

预计半数人的寿命超过 5。4 年。
第一个四分之一的人在 1。3 年内死亡,而前四分之三的人的寿命超过 1.3 岁。
前三分之三的人在 13。7 年内死亡,而前四分之一的人死亡工夫超过 13.7 岁。

估计量的图形示意(基于应用 KM 的生存曲线)。

 


参数估算器

咱们将思考三种常见的抉择:指数,Weibull 和 log-logistic 模型。

flexsurvreg(formula = su_obj ~ 1, data = orca, dist = "exponential")Estimates:       est      L95%     U95%     se     rate  0.11967  0.10513  0.13621  0.00791N = 338,  Events: 229,  Censored: 109Total time at risk: 1913.673Log-likelihood = -715.1802, df = 1AIC = 1432.36

同样,能够用非参数估计器图形地比拟不同的办法 

 

 


生存曲线的比拟

例如,肿瘤阶段是癌症存活钻研中的重要预后因素。咱们能够预计和绘制不同色彩的不同组(阶段)的生存曲线。

  stage  D       Y  x      pt       rate      lower     upper conf.level1     I 25 336.776 25 336.776 0.07423332 0.04513439 0.1033322       0.952    II 51 556.700 51 556.700 0.09161128 0.06646858 0.1167540       0.953   III 51 464.836 51 464.836 0.10971611 0.07960454 0.1398277       0.954    IV 57 262.552 57 262.552 0.21709985 0.16073995 0.2734597       0.955  unkn 45 292.809 45 292.809 0.15368380 0.10878136 0.1985862       0.95

通常,与具备高阶段肿瘤的患者相比,具备较低阶段肿瘤的诊断患者具备较低的(死亡率)。能够应用该 survfit() 函数执行生存函数的整体比拟。

Call: survfit(formula = su_obj ~ stage, data = orca)            n events median 0.95LCL 0.95UCLstage=I    50     25  10.56    6.17      NAstage=II   77     51   7.92    4.92   13.34stage=III  72     51   7.41    3.92    9.90stage=IV   68     57   2.00    1.08    4.82stage=unkn 71     45   3.67    2.83    8.17

因为低肿瘤阶段的发病率较低,因而肿瘤分期减少的中位生存工夫也会缩小。能够察看到雷同的行为,别离针对不同的肿瘤阶段绘制 KM 存活曲线。

也能够为每个阶段级别构建整个生存表。这里是每个肿瘤阶段生存表的前 3 行。

# Groups:   strata [5]
    time n.risk n.event n.censor  surv std.err upper lower strata   <dbl>  <dbl>   <dbl>    <dbl> <dbl>   <dbl> <dbl> <dbl> <fct>  1 0.17      50       1        0 0.98   0.0202 1     0.942 I      2 0.498     49       1        0 0.96   0.0289 1     0.907 I      3 0.665     48       1        0 0.94   0.0357 1     0.876 I      4 0.419     77       1        0 0.987  0.0131 1     0.962 II     5 0.498     76       1        0 0.974  0.0186 1     0.939 II     6 0.665     75       1        0 0.961  0.0229 1     0.919 II     7 0.167     72       1        0 0.986  0.0140 1     0.959 III    8 0.249     71       1        0 0.972  0.0199 1     0.935 III    9 0.413     70       1        0 0.958  0.0246 1     0.913 III   10 0.085     68       2        0 0.971  0.0211 1     0.931 IV    11 0.162     66       1        0 0.956  0.0261 1     0.908 IV    12 0.167     65       1        0 0.941  0.0303 0.999 0.887 IV    13 0.162     71       1        0 0.986  0.0142 1     0.959 unkn  14 0.167     70       2        0 0.958  0.0249 1     0.912 unkn  15 0.17      68       1        0 0.944  0.0290 0.999 0.892 unkn  
 arrange_ggsurvplots(glist, print = TRUE, ncol = 2, nrow = 1)

Mantel-Haenszel logrank 测试

默认参数 rho = 0 实现 log-rank 或 Mantel-Haenszel 测试。

Call:
survdiff(formula = su_obj ~ stage, data = orca)            N Observed Expected (O-E)^2/E (O-E)^2/Vstage=I    50       25     39.9     5.573     6.813stage=II   77       51     63.9     2.606     3.662stage=III  72       51     54.1     0.174     0.231stage=IV   68       57     33.2    16.966    20.103stage=unkn 71       45     37.9     1.346     1.642 Chisq= 27.2  on 4 degrees of freedom, p= 2e-05 

Peto&Peto Gehan-Wilcoxon 测试

survdiff(formula = su_obj ~ stage, data = orca, rho = 1)            N Observed Expected (O-E)^2/E (O-E)^2/Vstage=I    50     14.5     25.2     4.500     7.653stage=II   77     29.3     39.3     2.549     4.954stage=III  72     30.7     33.8     0.284     0.521stage=IV   68     40.3     22.7    13.738    21.887stage=unkn 71     32.0     25.9     1.438     2.359 Chisq= 30.9  on 4 degrees of freedom, p= 3e-06 

依据失败工夫,不同的测试应用不同的权重来比拟生存函数。在理论例子中,他们给出了可比拟的后果,表明不同肿瘤阶段的生存性能是不同的。


建模生存数据

当比拟因子程度的生存函数时,非参数检验特地可行。它们十分弱小,高效,通常简略 / 直观。
然而,随着感兴趣因素的数量减少,非参数测试变得难以进行和解释。相同,回归模型对于摸索生存与预测因子之间的关系更为灵便。

咱们将介绍两种不同的宽泛模型:半参数(即比例危险)和参数(减速生效工夫)模型。

CoxPH 模型

在咱们的例子中,咱们将思考将死亡工夫建模为性能性别,年龄和肿瘤阶段。
能够应用包装中的 coxph() 性能来装置 Cox 比例危险模型survival

 summary(m1)
Call:coxph(formula = su_obj ~ sex + I((age - 65)/10) + stage, data = orca)  n= 338, number of events= 229                     coef exp(coef) se(coef)     z Pr(>|z|)sexMale          0.35139   1.42104  0.14139 2.485 0.012947I((age - 65)/10) 0.41603   1.51593  0.05641 7.375 1.65e-13stageII          0.03492   1.03554  0.24667 0.142 0.887421stageIII         0.34545   1.41262  0.24568 1.406 0.159708stageIV          0.88542   2.42399  0.24273 3.648 0.000265stageunkn        0.58441   1.79393  0.25125 2.326 0.020016                 exp(coef) exp(-coef) lower .95 upper .95sexMale              1.421     0.7037    1.0771     1.875I((age - 65)/10)     1.516     0.6597    1.3573     1.693stageII              1.036     0.9657    0.6386     1.679stageIII             1.413     0.7079    0.8728     2.286stageIV              2.424     0.4125    1.5063     3.901stageunkn            1.794     0.5574    1.0963     2.935Concordance= 0.674  (se = 0.02)Rsquare= 0.226   (max possible= 0.999)Likelihood ratio test= 86.76  on 6 df,   p=<2e-16Wald test            = 80.5  on 6 df,   p=3e-15Score (logrank) test = 82.86  on 6 df,   p=9e-16

咱们能够应用函数检查数据是否与每个变量的比例危险假如别离和全局统一cox.zph()

                      rho    chisq     p
sexMale          -0.00137 0.000439 0.983I((age - 65)/10)  0.07539 1.393597 0.238stageII          -0.04208 0.411652 0.521stageIII         -0.06915 1.083755 0.298stageIV          -0.10044 2.301780 0.129stageunkn        -0.09663 2.082042 0.149GLOBAL                 NA 4.895492 0.557

显然没有找到拥护比例假如的证据。

Cox 模型的结果表明性别,年龄和阶段的显着影响。特地是,每减少 10 年,死亡率就会减少 50%。与男性和女性相比,全因死亡率的 HR 为 1.42。此外,估计数中第一阶段和第二阶段之间未发现任何差别。另一方面,阶段未知的群体是来自不同实在阶段的患者的简单混合物。因而,审慎的做法是将这些主题从数据中排除,并将前两个阶段组合为一个。

round(ci.exp(m2), 4)
                 exp(Est.)   2.5%  97.5%sexMale             1.3284 0.9763 1.8074I((age - 65)/10)    1.4624 1.2947 1.6519st3III              1.3620 0.9521 1.9482st3IV               2.3828 1.6789 3.3818

显示和图形化比拟多变量 Cox 模型的后果的便捷形式是通过森林图。

 

让咱们逐渐绘制预测的生存曲线,依据拟合的模型确定性别和年龄的值 

newd
      sex age  st3 id1    Male  40 I+II  12  Female  40 I+II  23    Male  80 I+II  34  Female  80 I+II  45    Male  40  III  56  Female  40  III  67    Male  80  III  78  Female  80  III  89    Male  40   IV  910 Female  40   IV 1011   Male  80   IV 1112 Female  80   IV 12


AFT 模型

参数模型假如存活工夫的散布。

Call:flexsurvreg(formula = Surv(time, all) ~ sex + I((age - 65)/10) +     st3, data = orca2, dist = "weibull")Estimates:                   data mean  est       L95%      U95%      se        exp(est)  L95%    shape                   NA    0.93268   0.82957   1.04861   0.05575        NA        NAscale                   NA   13.53151   9.97582  18.35456   2.10472        NA        NAsexMale            0.53184   -0.33905  -0.66858  -0.00951   0.16813   0.71245   0.51243I((age - 65)/10)  -0.15979   -0.41836  -0.54898  -0.28773   0.06665   0.65813   0.57754st3III             0.26966   -0.32567  -0.70973   0.05839   0.19595   0.72204   0.49178st3IV              0.25468   -0.95656  -1.33281  -0.58030   0.19197   0.38421   0.26374                  U95%    shape                   NAscale                   NAsexMale            0.99053I((age - 65)/10)   0.74996st3III             1.06012st3IV              0.55973N = 267,  Events: 184,  Censored: 83Total time at risk: 1620.864Log-likelihood = -545.858, df = 6AIC = 1103.716

能够证实,假如指数或威布尔散布的 AFT 模型能够从新参数化为比例危险模型(具备来自指数 / 威布尔散布族的基线危险)。
这能够应用包中的 weibreg() 性能显示eha

Call:weibreg(formula = Surv(time, all) ~ sex + I((age - 65)/10) +     st3, data = orca2)Covariate           Mean       Coef Exp(Coef)  se(Coef)    Wald psex           Female    0.490     0         1           (reference)            Male    0.510     0.316     1.372     0.156     0.043 I((age - 65)/10)   -0.522     0.390     1.477     0.062     0.000 st3             I+II    0.551     0         1           (reference)             III    0.287     0.304     1.355     0.182     0.095               IV    0.162     0.892     2.440     0.178     0.000 log(scale)                    2.605    13.532     0.156     0.000 log(shape)                   -0.070     0.933     0.060     0.244 Events                    184 Total time at risk        1620.9 Max. log. likelihood      -545.86 LR test statistic         68.7 Degrees of freedom        4 Overall p-value           4.30767e-14

系数的(指数)具备与 Cox 比例模型的系数的等效解释(预计也相似)。

通过将参数提供 fnsummaryplot 办法,能够汇总或绘制拟合模型的参数的任何函数。例如,Weibull 模型下的中位存活率能够概括为

newd <- data.frame(sex = c("Male", "Female"), age = 65, st3 = "I+II")summary(m2w, newdata = newd, fn = median.weibull, t = 1, B = 10000)
sex=Male, I((age - 65)/10)=0, st3=I+II   time      est      lcl      ucl1    1 6.507834 4.898889 8.631952sex=Female, I((age - 65)/10)=0, st3=I+II   time      est      lcl      ucl1    1 9.134466 6.801322 12.33771

将后果与 Cox 模型的后果进行比拟。

survfit(m2, newdata = newd)
Call: survfit(formula = m2, newdata = newd)    n events median 0.95LCL 0.95UCL1 267    184   7.00    5.25    10.62 267    184   9.92    7.33    13.8

泊松回归

能够证实,Cox 模型在数学上等效于对数据的特定变换的泊松回归模型。
这个想法是每次在每个工夫距离仅蕴含一个事件时以这种形式察看事件时宰割后续工夫。在这个加强的数据集中,能够屡次示意主题(多行)。

咱们首先定义察看事件(all == 1)的惟一工夫,并应用包中的 survSplit() 函数 survival 来创立加强或宰割数据。

head(orca_splitted, 15)

包中的 gnm() 函数 gnm 适宜决裂数据上的条件泊松,其中工夫的影响(作为因子变量)能够被边缘化(不预计以进步计算效率)。

mod_poi <- gnm(all ~ sex + I((age-65)/10) + st3, data = orca_splitted,                family = poisson, eliminate = factor(time))summary(mod_poi)

将从条件 Poisson 取得的估计值与 cox 比例危险模型进行比拟。

round(data.frame(cox = ci.exp(m2), poisson = ci.exp(mod_poi)), 4)
                 cox.exp.Est.. cox.2.5. cox.97.5. poisson.exp.Est.. poisson.2.5. poisson.97.5.sexMale                 1.3284   0.9763    1.8074            1.3284       0.9763        1.8074I((age - 65)/10)        1.4624   1.2947    1.6519            1.4624       1.2947        1.6519st3III                  1.3620   0.9521    1.9482            1.3620       0.9521        1.9482st3IV                   2.3828   1.6789    3.3818            2.3828       1.6789        3.3818

如果咱们想要预计基线危险,咱们还须要预计泊松模型中工夫的影响(OBS 咱们还须要包含工夫距离的(对数)长度作为偏移)。

orca_splitted$dur <- with(orca_splitted, time - tstart)mod_poi2 <- glm(all ~ -1 + factor(time) + sex + I((age-65)/10) + st3,                 data = orca_splitted, family = poisson, offset = log(dur))

基线危险包含阶梯函数,其中速率在每个工夫距离内是恒定的。

newd <- data.frame(time = cuts, dur = 1,                   sex = "Female", age = 65, st3 = "I+II")blhaz <- data.frame(ci.pred(mod_poi2, newdata = newd))ggplot(blhaz, aes(x = c(0, cuts[-138]), y = Estimate, xend = cuts, yend = Estimate)) + geom_segment() +  scale_y_continuous(trans = "log", limits = c(.05, 5), breaks = pretty_breaks()) +  theme_classic() + labs(x = "Time (years)", y = "Baseline hazard")

 

更好的办法是通过应用例如具备节点 \(k \)的样条来灵便地模仿基线危险。

                     exp(Est.)  2.5%  97.5%(Intercept)              0.074 0.040  0.135ns(time, knots = k)1     0.402 0.177  0.912ns(time, knots = k)2     1.280 0.477  3.432ns(time, knots = k)3     0.576 0.220  1.509ns(time, knots = k)4     1.038 0.321  3.358ns(time, knots = k)5     4.076 0.854 19.452ns(time, knots = k)6     1.040 0.171  6.314sexMale                  1.325 0.975  1.801I((age - 65)/10)         1.469 1.300  1.659st3III                   1.360 0.952  1.942st3IV                    2.361 1.665  3.347

比拟不同的策略

咱们能够依据特定协变量模式的预测生存曲线比拟之前的策略,称 65 岁的女性患有肿瘤 I 期或 II 期。

newd <- data.frame(sex = "Female", age = 65, st3 = "I+II")

生存函数的图形示意便于比拟。

 


其余剖析

非线性

咱们假如年龄对(log)死亡率的影响是线性的。放宽这一假如的可能策略是适宜 Cox 模型,其中年龄用二次效应建模。

Call:coxph(formula = Surv(time, all) ~ sex + I(age - 65) + I((age -     65)^2) + st3, data = orca2)  n= 267, number of events= 184                      coef exp(coef)  se(coef)     z Pr(>|z|)sexMale         2.903e-01 1.337e+00 1.591e-01 1.825   0.0681I(age - 65)     3.868e-02 1.039e+00 6.554e-03 5.902 3.59e-09I((age - 65)^2) 9.443e-05 1.000e+00 3.576e-04 0.264   0.7917st3III          3.168e-01 1.373e+00 1.838e-01 1.724   0.0847st3IV           8.691e-01 2.385e+00 1.787e-01 4.863 1.16e-06                exp(coef) exp(-coef) lower .95 upper .95sexMale             1.337     0.7481    0.9787     1.826I(age - 65)         1.039     0.9621    1.0262     1.053I((age - 65)^2)     1.000     0.9999    0.9994     1.001st3III              1.373     0.7284    0.9576     1.968st3IV               2.385     0.4193    1.6801     3.385Concordance= 0.674  (se = 0.022)Rsquare= 0.216   (max possible= 0.999)Likelihood ratio test= 64.89  on 5 df,   p=1e-12Wald test            = 63.11  on 5 df,   p=3e-12Score (logrank) test = 67.64  on 5 df,   p=3e-13

非线性(即二次项)的值很高,因而没有证据能够回绝零假如(即线性假如是适合的)。

如果关系是非线性的,则年龄系数不再能够间接解释。咱们能够将 HR 作为年龄的函数以图形形式出现。咱们须要指定一个批示值; 咱们抉择 65 岁的中位年龄值。

age <- seq(20, 80, 1) - 65  geom_vline(xintercept = 65, lty = 2) + geom_hline(yintercept = 1, lty = 2)

 

工夫依赖系数

cox.zph() 函数可用于绘制个体预测因子随工夫的影响,因而可用于诊断和了解非比例危险。

 

咱们能够通过拟合的阶梯函数来放宽比例危险假如,这意味着在不同的工夫距离内有不同的 

包中的 survSplit() 函数 survival 将数据集合成为工夫相干的局部。

  id    sex      age stage          event  st3 tstart  time all tgroup1  2 Female 83.08783   III Oral ca. death  III      0 0.419   1      12  3   Male 52.59008    II    Other death I+II      0 5.000   0      13  3   Male 52.59008    II    Other death I+II      5 7.915   1      24  4   Male 77.08630     I    Other death I+II      0 2.480   1      15  5   Male 80.33622    IV Oral ca. death   IV      0 2.500   1      16  6 Female 82.58132    IV    Other death   IV      0 0.167   1      1
    I((age - 65)/10) + st3, data = orca3)                                                 coef exp(coef) se(coef)      z        pI((age - 65)/10)                              0.38184   1.46498  0.06255  6.104 1.03e-09st3III                                        0.28857   1.33451  0.18393  1.569   0.1167st3IV                                         0.87579   2.40076  0.17963  4.876 1.09e-06relevel(sex, 2)Male:strata(tgroup)tgroup=1    0.42076   1.52312  0.19052  2.209   0.0272relevel(sex, 2)Female:strata(tgroup)tgroup=1       NA        NA  0.00000     NA       NArelevel(sex, 2)Male:strata(tgroup)tgroup=2   -0.10270   0.90240  0.28120 -0.365   0.7149relevel(sex, 2)Female:strata(tgroup)tgroup=2       NA        NA  0.00000     NA       NArelevel(sex, 2)Male:strata(tgroup)tgroup=3    1.13186   3.10142  1.09435  1.034   0.3010relevel(sex, 2)Female:strata(tgroup)tgroup=3       NA        NA  0.00000     NA       NALikelihood ratio test=68.06  on 6 df, p=1.023e-12n= 416, number of events= 184 

尽管不显着,但男女比拟的危险比在第二期间(5 至 15 年)低于 1,而在其余两个期间高于 1。该 cox.zph() 函数可用于查看在宰割剖析中是否依然偏离比例假如。

模仿生存百分位数

一个不同但乏味的视角包含模仿生存工夫的百分位数。

Call:ctqr(formula = Surv(time, all) ~ st3, data = orca2, p = p)Coefficients:             p = 0.25(Intercept)   2.665  st3III       -1.369  st3IV        -1.877  Degrees of freedom: 267 total; 225 residuals

β0= 2.665β0=2.665 是参考组中死亡概率等于 0.25 的工夫。另一个被解释为绝对度量。

该信息能够直观地比拟在肿瘤阶段的程度上别离预计的存活曲线。

  p = c(p, p - .005, p + .005))[-1,]  = 1 - p, xend = time_ref,                                  yend = 1 - p))

 

对 Cox 模型中示意的那个进行补充剖析_m2_评估生存工夫百分位数的可能差别,作为诊断,性别和肿瘤阶段年龄的函数

ctqr(formula = Surv(time, all) ~ sex + I((age - 65)/10) + st3,     data = orca2, p = seq(0.1, 0.7, 0.1))Coefficients:                  p = 0.1   p = 0.2   p = 0.3   p = 0.4   p = 0.5   p = 0.6   p = 0.7 (Intercept)        1.44467   2.44379   4.65302   7.73909  10.81386  12.18348  15.19359sexMale           -0.09218  -0.27385  -0.85720  -2.49580  -3.27962  -2.81428  -4.01656I((age - 65)/10)  -0.19026  -0.39819  -1.20278  -1.93144  -2.39229  -3.03915  -3.52711st3III            -0.60994  -1.08534  -1.89357  -2.23741  -3.10478  -2.00037  -1.59213st3IV             -1.07679  -1.59566  -2.92700  -3.16652  -4.74759  -4.80838  -5.25810Degrees of freedom: 267 total; 220 residuals

后果包含不同百分位数下每种协变量的存活工夫差别。图形示意能够促成后果的解释。

coef_q <- data.frame(coef(fit_q)) %>%    .96 * se    )

 

或者,能够针对一组特定的协方差模式预测存活工夫的百分位数。


CIF 累积发生率函数

对于因某一特定起因而死亡的危险或危险,通常是次要的趣味。因为竞争起因阻止受试者倒退事件,可能无奈察看到起因特定事件。竞争事件不仅产生在特定起因的死亡率上,而且每次事件都会阻止并发事件产生时更多。
在咱们的例子中,咱们感兴趣的是模仿口腔癌死亡危险,并且死于其余起因将被视为竞争事件。

在竞争危险情景中,事件特定的生存取得审查其余事件(也就是天真的 Kaplan-Meier 对特定起因生存的预计)通常是不适合的。
咱们将思考事件的累积发生率函数(CIF)

CIF

包中的 Cuminc() 函数 mstate 计算竞争事件的非参数 CIF(也称为 Aalen-Johansen 预计)和相干的标准误差。

 head(cif)
   time      Surv CI.Oral ca. death CI.Other death      seSurv seCI.Oral ca. death1 0.085 0.9925094       0.007490637    0.000000000 0.005276805         0.0052768052 0.162 0.9887640       0.011235955    0.000000000 0.006450534         0.0064505343 0.167 0.9812734       0.011235955    0.007490637 0.008296000         0.0064505344 0.170 0.9775281       0.011235955    0.011235955 0.009070453         0.0064505345 0.249 0.9737828       0.011235955    0.014981273 0.009778423         0.0064505346 0.252 0.9662921       0.014981273    0.018726592 0.011044962         0.007434315  seCI.Other death1      0.0000000002      0.0000000003      0.0052768054      0.0064505345      0.0074343156      0.008296000

咱们能够绘制 CIF(一个重叠的另一个)以及派生的无事件生存函数。

 

曾经施行了扩大以通过因子变量的程度来预计累积发生率函数。

grid.arrange(
  ncol = 2)

 

咱们能够看到,IV 期口腔癌死亡的 CIF 高于 III,甚至更高于 I + II。相同,对于其余起因死亡率,曲线仿佛不随肿瘤阶段而变动。

当咱们想要在竞争危险设置中对生存数据进行建模时,有两种常见的策略能够解决不同的问题:
– 针对事件特定危害的 Cox 模型,例如,趣味在于预测因素对死亡率的生物效应十分疾病,往往导致相干的后果。
– 当咱们想要评估因子对事件总体累积发生率的影响时,细分模型的细分模型为。Cc

CIF Cox 模型

 round(ci.exp(m2haz2), 4)
                 exp(Est.)   2.5%  97.5%sexMale             1.8103 1.1528 2.8431I((age - 65)/10)    1.4876 1.2491 1.7715st3III              1.2300 0.7488 2.0206st3IV               1.6407 0.9522 2.8270

起因特异性 Cox 模型的后果与起因特异性 CIF 的图形示意统一,即肿瘤 IV 期仅是口腔癌死亡率的重要危险因素。年龄减少与两种起因的死亡率减少相干(口腔癌死亡率 HR = 1.42,其余起因死亡率 HR = 1.48)。仅依据其余起因死亡率察看到性别差异(HR = 1.8)。

CRR 模型

crr()cmprsk 竞争危险的状况下,包中的函数可用于子散布函数的回归建模。咱们应用与上述雷同的协变量,给出了针对口腔癌死亡和其余起因死亡的散布危险的细灰模型的后果。

Call:crr(ftime = time, fstatus = event, cov1 = model.matrix(m2), failcode = "Oral ca. death")                    coef exp(coef) se(coef)      z p-valuesexMale          -0.0953     0.909    0.213 -0.447 6.5e-01I((age - 65)/10)  0.2814     1.325    0.093  3.024 2.5e-03st3III            0.3924     1.481    0.258  1.519 1.3e-01st3IV             1.0208     2.775    0.233  4.374 1.2e-05                 exp(coef) exp(-coef)  2.5% 97.5%sexMale              0.909      1.100 0.599  1.38I((age - 65)/10)     1.325      0.755 1.104  1.59st3III               1.481      0.675 0.892  2.46st3IV                2.775      0.360 1.757  4.39Num. cases = 267Pseudo Log-likelihood = -501 Pseudo likelihood ratio test = 31.4  on 4 df,
m2fg2 <- with(orca2, crr(time, event, cov1 = model.matrix(m2), failcode = "Other death"))summary(m2fg2, Exp = T)
Competing Risks RegressionCall:crr(ftime = time, fstatus = event, cov1 = model.matrix(m2), failcode = "Other death")                   coef exp(coef) se(coef)      z p-valuesexMale           0.544     1.723   0.2342  2.324   0.020I((age - 65)/10)  0.197     1.218   0.0807  2.444   0.015st3III            0.130     1.139   0.2502  0.521   0.600st3IV            -0.212     0.809   0.2839 -0.748   0.450                 exp(coef) exp(-coef)  2.5% 97.5%sexMale              1.723      0.580 1.089  2.73I((age - 65)/10)     1.218      0.821 1.040  1.43st3III               1.139      0.878 0.698  1.86st3IV                0.809      1.237 0.464  1.41Num. cases = 267Pseudo Log-likelihood = -471 Pseudo likelihood ratio test = 9.43  on 4 df, 还有问题,分割咱们!

正文完
 0