我在前面的推文中详细的介绍了使用ggplot2包绘制生存曲线,虽然R中ggplot2包绘制出来的图形很优美,但是基础绘图函数也不差。
今天使用基础绘图包来绘制KM生存曲线。
需要的包为:survival包和基础绘图包。
使用survival包的lung数据集,运行以下代码加载数据:
install.packages("survival") # 安装包library(survival) # 加载包data(lung) # 加载数据集View(lung) # 预览数据集names(lung) # 输出数据框变量名称
数据集中变量介绍见常用内置数据集介绍。
R中使用Surv()函数创建生存对象。
attach(lung)Surv(time, status)# 或 with(lung, Surv(time, status)) 等效
上图为输出的生存对象;带"+"号表示为右删失数据。
is.Surv(x) # 可判断x是否为生存对象
library(survival)fit <- survfit(Surv(time, status) ~ sex, data = lung)fit # 可输出部分信息
summary(fit) # 查看拟合对象的详细信息# 结果太长,贴部分数据。
图中的survival为生存函数在生存时间点处的KM估计值。同时输出结果还给出了估计的标准误差和95%置信区间。
使用survival包中plot.survfit()
函数用来绘制生存曲线
plot(fit)
plot(fit, conf.int = T) # 实线为生存曲线,虚线为置信区间,默认绘制95%。# plot(fit, conf.int = 0.99) # 为数字,表示绘制99%置信区间
plot(fit, conf.int = T, mark.time = TRUE)
plot(fit, conf.int = T, mark.time = TRUE, pch = c(1,2), cex=2)
plot(fit, conf.int = T, col = c("#e42c64","#2C2C2C"), lwd = 2, xlab = "follow-up time", # 设置x轴标签 ylab = "Survival probability", # 设置y轴标签 main = "KM curve") # 设置图形标题
legend(800, 0.8, c("Male", "Female"), lwd = 2, lty = 1, col = c("#e42c64","#2C2C2C")) # 设置图例坐标和标签
生存率的比较见这篇生存率的比较。
这里就不计算了,假设p值为0.001。
text(200,0.2,"p-value = 0.001")
至此,图形基本绘制完成了,还有一些细节方面可以参考基础绘图函数的图形参数来设置。
Surv(time, time2, event, type=c('right', 'left', 'interval', 'counting', 'interval2', 'mstate'), origin=0)is.Surv(x) # 判断x是否为生存对象
time:# 对于右删失数据,time表示随访时间;对于区间数据,time表示事件起始时间。event:# 生存状态,多以0/1、1/2或FALSE/TRUE来表示删失/死亡;对于区间删失数据,一般用0/1/2/3表示右删失/事件时间/左删失/区间删失;# 对于多个终点数据,event变量是一个因子变量,第一个水平认为是删失。 time2:# 区间的结束时间仅用于区间删失或计数过程数据,区间为左开右闭区间(start, end),对于区间删失数据,time2表示事件结束时间;# 对计数过程,event指标用来表示区间右端事件是否发生。type:# 字符串变量,用来指定删失的类型,可能的取值有:right/left表示右/左删失;counting表示计数过程;# interval/interval2表示区间删失;还有一个mstate类型。origin:# 适用于计数过程,用以给定风险函数的初值。
当type参数没有设定时,将按以下规则定义type参数:
如果有两个未命名的参数,将按time、event顺序匹配;
如果有三个未命名的参数,将按time、time2和event的顺序进行匹配;
如果event变量是一个因子变量,设定type类型为mstate;如果没有time2参数,则type类型为右删失
attach(heart)Surv(start, stop, event)
输出的生存对象为区间形式。
删失可细分为左删失(left censored)、区间删失(interval censored)和右删失(right censored)3类。
如果只知道感兴趣终点事件会在目前知晓时间(如截至时间、失访时间、死于其他疾病时间)之前发生,则称为左删失;
如果只知道感兴趣终点事件会在某一区间内发生,则称为区间删失;
如果只知道感兴趣终点事件会在知晓时间之后发生,则称为右删失。
右删失在实际工作中最常见,即大多数情况下获得的删失个体生存时间应该比知晓的时间更长。采用“+”号表示个体为右删失。
plot(x, conf.int=, mark.time=FALSE, pch=3, col=1, lty=1, lwd=1, cex=1, log=FALSE, xscale=1, yscale=1, xlim, ylim, xmax, fun, xlab="", ylab="", xaxs="r", conf.times, conf.cap=.005, conf.offset=.012, conf.type = c("log", "log-log", "plain", "logit", "arcsin"), mark, noplot="(s0)", cumhaz=FALSE, firstx, ymin, ...)
x:survfit函数拟合的生存对象; conf.int:确定是否添加置信区间,输入TRUE,则默认绘制 95% 置信区间;# 输入数字 0.99,绘制 99% 置信区间## 置信条conf.times:一个数字向量,指定在x轴上某时间点显示置信条,如果指定了向量,则将使用向量中的值代替置信区间。 conf.cap:置信条顶部的水平帽的宽度;仅在conf.times有效时使用;1表示绘图区域的宽度 conf.offset:当图上有多条曲线时,置信条的偏移度;1表示图形区域的宽度;如果为单个数字,表示每个曲线条从先前曲线条的偏移度;如果为向量,则设置每条曲线的偏移度。conf.type:默认log,不能使用缩写,要使用全称。包括以下四种:# plain表示不给出置信区间# log表示基于累计风险或log(survival)计算区间;# log-log基于log风险或log(-log(survival))计算区间# logit基于log(survival/(1-survival)))计算区间;
mark.time:是否添加删失点标签;为FALSE,不绘制标签;为TRUE,在每一个删失时间点绘制标签。 如果mark.time为指定时间点的数字向量,则在曲线上指定时间点添加删失点标签。## 删失点的符号pch:为单个数字,则所有曲线删失点符号类型一样; 为一字符向量,则根据向量中的值按顺序更改曲线上删失点的符号类型; 单个字符串如abcd可视为c("a","b","c","d");如果向量的长度小于曲线的条数,则循环使用向量直至补齐。mark:pch的曾用名。noplot:对于多状态模型,曲线的删失点不会绘制。 ## 删失点的大小cex:为单个数字,则指定删失点的大小,所有标记的大小相同; 表示为向量,则按向量中的值分别修改删失点的大小。
## 设置线条颜色、类型、线宽col:一个整数向量,用来指定曲线的颜色,默认为1。 lty:一个整数向量,用于指定每条曲线的线型。默认为1。 lwd:一个数字向量,用来指定线条宽度,默认为1。## 设置坐标轴范围、标签xlim,ylim:设置x轴和y轴坐标轴范围;ymin:y轴最小值,通常作为ylim参数的一部分。firstx:x轴起始点,通常作为xlim参数的一部分。xmax:设置水平坐标轴的最大值;在绘制曲线之前缩短曲线,因此和xlim参数不同,不会产生超出坐标轴范围的警告信息.xlab:设置x轴标签ylab:设置y轴标签## 设置坐标轴显示标签xscale:与yscale类似,给定数字365.25将以年为单位而不是原来的日期。yscale:一个数字,用来乘以y轴上的标签,用来改变y轴标签显示方式,但不改变实际绘图坐标的意义;比如说y轴上的标签为0%-100%,将标签乘以100,得到标签为0-100的显示方式。log:逻辑词,为TRUE,则y轴将变换为对数坐标轴;或者指定"x","y",或"xy"将分别变换x轴、y轴或xy轴为对数坐标轴## 设置置信区间
fun:定义生存曲线转换的函数;使用字符参数来指定5个常见的转换函数:# "S"绘制常见的生存曲线;# "log"绘制对数生存曲线,但是坐标轴标签为log(S)值,与log=T等效;# "event"或"F"绘制empirical CDF F(t)= 1-S(t)(f(y) = 1-y);# "cloglog"创建一个互补的log-log生存图。# "cumhaz"绘制累积风险函数。cumhaz:逻辑词,为TRUE,则绘制累积风险曲线。
联系客服