11
11
12
12
library(survival) # survfit
13
13
library(ggplot2)
14
- library(ggfortify) # autoplot
15
14
library(glmnet) # Cox Models
16
- library(VGAM) # R >= 4.4.0
17
- library(INLA)
15
+ library(ggsurvfit)
16
+ # library(VGAM) # R >= 4.4.0
17
+ # library(INLA)
18
18
```
19
19
20
- 生存分析可以用于用户流失分析,注册、激活、活跃。 分析次日留存、7日留存、15日留存。有学生来上体验课,多久来付费上课。 有一个人医院看病之后,多久办理住院。 最早,生存分析用于研究飞机出去之后,有多少返回的。还是要回归到原始文献去了解基本概念,及其背后的思考和应用
20
+ 生存分析可以用于用户流失分析,注册、激活、活跃。 分析次日留存、7日留存、15日留存。有学生来上体验课,多久来付费上课。 有一个人医院看病之后,多久办理住院。 最早,生存分析用于研究飞机出去之后,有多少还能返回的。生存分析的学习还是要回归到原始文献去了解基本概念,及其背后的思考和应用。
21
+
22
+ 以一个生存问题引出本章主题,讲述和展示一个数据集,先探索和分析数据,之后建立和拟合模型,结果解释。
21
23
22
- 以一个问题提出本章主题,讲述和展示一个数据集。建立模型,拟合模型,结果解释。
24
+ lung 数据集 ** survival ** 包(模型)和 ** ggsurvfit ** 包(可视化)
23
25
24
26
## 问题背景 {#sec-aml}
25
27
@@ -28,7 +30,7 @@ library(INLA)
28
30
``` {r}
29
31
library(survival)
30
32
data(cancer, package = "survival")
31
- str(aml )
33
+ str(lung )
32
34
```
33
35
34
36
数据的分布情况如下
@@ -40,7 +42,7 @@ str(aml)
40
42
#| fig-width: 4.5
41
43
#| fig-height: 3
42
44
43
- ggplot(data = aml , aes(x = time, y = status, color = x )) +
45
+ ggplot(data = lung , aes(x = time, y = status, color = factor(sex) )) +
44
46
geom_jitter(height = 0.2) +
45
47
theme_minimal()
46
48
```
@@ -55,46 +57,46 @@ Cox 比例风险回归模型与 Box-Cox 变换 [@Box1964]
55
57
- ` MASS::boxcox() ` Box-Cox 变换
56
58
- ` glmnet::glmnet(family = "cox") `
57
59
- INLA 包的函数 ` inla() ` 与 ` inla.surv() ` 一起拟合,[ 链接] ( https://becarioprecario.bitbucket.io/inla-gitbook/ch-survival.html )
58
- - [ survstan] ( https://github.com/fndemarqui/survstan ) Stan 与生存分析
59
- - rstanarm 包的函数 ` stan_jm() ` 使用说明 Estimating Joint Models for Longitudinal and Time-to-Event Data with rstanarm [ 链接] ( https://cran.r-project.org/web/packages/rstanarm/vignettes/jm.html )
60
- - rstanarm 包的[ 生存分析分支] ( https://github.com/stan-dev/rstanarm/pull/323 )
61
60
62
61
### survival
63
62
64
63
R 软件内置了 [ survival] ( https://github.com/therneau/survival ) 包,它是实现生存分析的核心 R 包 [ @Terry2000 ] ,其函数 ` survfit() ` 拟合模型。
65
64
66
65
``` {r}
67
- aml_survival <- survfit(Surv(time, status) ~ x , data = aml )
68
- summary(aml_survival)
66
+ lung_surv <- survfit(Surv(time, status) ~ sex , data = lung )
67
+ lung_surv
69
68
```
70
69
71
70
拟合 Cox 比例风险回归模型(Cox Proportional Hazards Regression Model)
72
71
73
72
``` {r}
74
- aml_coxph <- coxph(Surv(time, status) ~ 1 + x , data = aml )
75
- summary(aml_coxph )
73
+ lung_coxph <- coxph(Surv(time, status) ~ 1 + sex , data = lung )
74
+ summary(lung_coxph )
76
75
```
77
76
78
- 展示拟合结果。可以绘制生存分析的图的 R 包有很多,比如 ggfortify 包、[ ggsurvfit] ( https://github.com/ddsjoberg /ggsurvfit/ ) 包和 [ survminer] ( https://github.com/kassambara/survminer ) 包等。ggfortify 包可以直接针对函数 ` survfit() ` 的返回对象绘图,[ ggsurvfit] ( https://github.com/ddsjoberg/ggsurvfit/ ) 包提供新函数 ` survfit2() ` 拟合模型、函数 ` ggsurvfit() ` 绘制图形,画面内容更加丰富,而 [ survminer] ( https://github.com/kassambara/survminer ) 包依赖很多。
77
+ 展示拟合结果。可以绘制生存分析的图的 R 包有很多,比如 ggfortify 包、[ ggsurvfit] ( https://github.com/pharmaverse /ggsurvfit ) 包和 [ survminer] ( https://github.com/kassambara/survminer ) 包等。ggfortify 包可以直接针对函数 ` survfit() ` 的返回对象绘图,ggsurvfit 包提供新函数 ` survfit2() ` 拟合模型、函数 ` ggsurvfit() ` 绘制图形,画面内容更加丰富,而 survminer 包依赖很多。
79
78
80
79
``` {r}
81
80
#| label: fig-leukemia-surv
82
81
#| fig-cap: "急性粒细胞白血病生存数据"
83
82
#| fig-showtext: true
84
83
#| fig-width: 6
85
- #| fig-height: 3
86
-
87
- library(ggplot2)
88
- library(ggfortify)
89
- autoplot(aml_survival, data = aml) +
90
- theme_minimal()
91
- ```
92
-
93
- 参数化的生存分析模型(参数模型,相对于非参数模型而言)
94
-
95
- ``` {r}
96
- aml_surv_reg <- survreg(Surv(time, status) ~ x, data = aml, dist = "weibull")
97
- summary(aml_surv_reg)
84
+ #| fig-height: 5
85
+
86
+ p <- survfit2(Surv(time, status) ~ sex, data = lung) |>
87
+ ggsurvfit(linewidth = 1) +
88
+ add_confidence_interval() +
89
+ add_risktable() +
90
+ add_quantile(y_value = 0.6, color = "gray50", linewidth = 0.75) +
91
+ scale_ggsurvfit()
92
+ p +
93
+ # limit plot to show 8 years and less
94
+ coord_cartesian(xlim = c(0, 1000)) +
95
+ # update figure labels/titles
96
+ labs(
97
+ y = "Percentage Survival",
98
+ title = "Recurrence by Time From Surgery to Randomization",
99
+ )
98
100
```
99
101
100
102
### glmnet
@@ -106,19 +108,22 @@ glmnet 包拟合 Cox 比例风险回归模型 [@simon2011] 适合需要多变量
106
108
107
109
library(glmnet)
108
110
# alpha = 1 lasso
109
- aml_glmnet <- glmnet(x = aml$x , y = Surv(aml $time, aml $status), family = "cox", alpha = 1)
110
- aml_glmnet_cv <- cv.glmnet(x = aml$x , y = Surv(aml $time, aml $status), family = "cox", alpha = 1)
111
+ lung_glmnet <- glmnet(x = lung$sex , y = Surv(lung $time, lung $status), family = "cox", alpha = 1)
112
+ lung_glmnet_cv <- cv.glmnet(x = lung$sex , y = Surv(lung $time, lung $status), family = "cox", alpha = 1)
111
113
```
112
114
113
115
### INLA
114
116
115
117
INLA 包拟合 Cox 比例风险回归模型 [ @Virgilio2020 ] 采用近似贝叶斯推断。
116
118
117
119
``` {r}
120
+ #| eval: false
121
+
118
122
library(INLA)
119
123
inla.setOption(short.summary = TRUE)
120
- aml_inla <- inla(inla.surv(time, status) ~ x, data = aml, family = "exponential.surv", num.threads = "1:1")
121
- summary(aml_inla)
124
+ lung_inla <- inla(inla.surv(time, status) ~ sex, data = lung,
125
+ family = "exponential.surv", num.threads = "1:1")
126
+ summary(lung_inla)
122
127
```
123
128
124
129
## Tobit 回归 {#sec-tobit-regression}
@@ -135,17 +140,7 @@ Tobit (Tobin's Probit) regression 起源于计量经济学中的 Tobit 模型,
135
140
136
141
library(VGAM) # Vector Generalized Linear and Additive Models
137
142
# VGAM::vglm(family = tobit(Upper = 800)) # Tobit regression
138
- ```
139
-
140
- ``` {r}
141
- library(VGAM)
142
- with(aml, SurvS4(time, status))
143
- ```
144
-
145
- ``` {r}
146
- #| eval: false
147
- #| echo: false
148
-
149
- aml_vglm <- vglm(SurvS4(time, status) ~ x, data = aml, family = cens.poisson)
150
- summary(aml_vglm)
143
+ with(lung, SurvS4(time, status))
144
+ lung_vglm <- vglm(SurvS4(time, status) ~ sex, data = lung, family = cens.poisson)
145
+ summary(lung_vglm)
151
146
```
0 commit comments