R语言小白绘图系列
第34弹·连续变量生存曲线
连续变量生存曲线简介
连续变量生存曲线（Continuous Variable Survival Curve）是指在生存分析中，对于一个连续变量（如年龄、血压、血糖等），通过对其不同取值的生存情况进行分析并绘制生存曲线。通常，生存分析用于研究时间到事件（如死亡、疾病复发等）的数据。
以下是连续变量生存曲线的一些关键概念和步骤：
连续变量的分类：
将连续变量划分为若干组（如分位数法、固定间隔法等），例如可以将年龄分为不同的年龄段。 或者采用平滑的方法，比如用Cox比例风险模型直接处理连续变量，而不进行分组。
生存曲线的绘制：
使用Kaplan-Meier方法绘制不同组的生存曲线。 如果使用Cox比例风险模型，则可以根据模型的预测值绘制平滑的生存曲线。
比较不同组的生存情况：
使用Log-rank检验比较不同组的生存曲线是否有显著差异。 在Cox模型中，通过检验连续变量的回归系数来判断其与生存时间的关系。
结果解释：
不同组的生存曲线可以反映出连续变量的不同取值对生存情况的影响。 在Cox模型中，连续变量的回归系数及其置信区间可以量化其影响大小。
应用实例
下图中：按照MATN3基因表达中位值将样本分为高表达组和低表达组；比较两组之间的生存率是否具有差异；横坐标代表生存时间，单位是年；纵坐标代表生存率；通过P值可以了解不同分期之间生存率是否存在差异，若p<0.05则说明不同分期的生存率存在差异；如下图，p<0.001,两组之间生存存在差异，说明该基因与生存相关。生存曲线下方表格代表不同时间点高低表达组存活的患者数目。
源文件
输入文件有四列信息
id: 样品名称
futime:生存时间，单位为年。
fustat:生存状态，0代表存活，1代表死亡。
连续变量：在这里为基因表达量
代码
环境准备
安装必要的包:
install.packages("survival")
install.packages("survminer")这两行代码用于安装
survival和survminer包。这两个包在R中用于生存分析和生存曲线的可视化。加载包:
library(survival)
library(survminer)这两行代码用于加载之前安装的
survival和survminer包。加载这些包后，可以使用其中的函数进行生存分析和绘图。设置输入输出文件和变量:
inputFile = "input.txt"
outFile = "survival.pdf"
var = "MATN3" # 用于生存分析的变量
inputFile变量设置了输入数据文件的名称为"input.txt"。outFile变量设置了输出文件的名称为"survival.pdf"，这是将保存的生存曲线图。var变量设置为"MATN3"，表示将使用输入数据中的MATN3列进行生存分析。
设置工作目录:
setwd("D:\\biowolf\\bioR\\35.survivalContinuous")
这行代码将当前工作目录设置为"D:\\biowolf\\bioR\\35.survivalContinuous"。这样，所有的文件读取和写入操作都会在这个目录下进行。
这段代码的目的是安装并加载必要的包，设置输入输出文件和变量，以及设置工作目录，为后续的生存分析做好准备。
数据准备
这段代码用于读取输入数据文件，将样本分为两组（高表达组和低表达组），然后进行生存分析并计算p值。以下是逐行解释：
rt = read.table(inputFile, header = TRUE, sep = "\t", check.names = FALSE)
这行代码读取制表符分隔的输入文件，并将其内容存储到数据框rt中。header = TRUE表示文件的第一行是列名，check.names = FALSE表示不自动修改列名。
rt = rt[, c("futime", "fustat", var)]
这行代码选择数据框rt中的三列进行分析：生存时间（futime）、生存状态（fustat）和变量var（即MATN3列）。这里假设输入数据文件中包含这些列。
# 根据中位值，把样品分为两组
group = ifelse(rt[, 3] > median(rt[, 3]), "High", "Low")
这行代码根据MATN3列的中位值将样品分为两组：高表达组（High）和低表达组（Low）。如果某样品的MATN3值大于中位值，则分为High组，否则分为Low组。
diff = survdiff(Surv(futime, fustat) ~ group, data = rt)
这行代码使用survdiff函数计算不同组间的生存差异。Surv(futime, fustat)创建一个生存对象，其中futime是生存时间，fustat是生存状态（事件发生为1，未发生为0）。~ group表示按group变量分组，data = rt指定数据来源。
pValue = 1 - pchisq(diff$chisq, df = 1)
这行代码计算卡方检验的p值。diff$chisq是survdiff计算的卡方值，df = 1是自由度（因为分为两组，自由度为1）。pchisq函数返回卡方分布的累积分布函数值，1减去该值即为p值。
if (pValue < 0.001) {
pValue = "p<0.001"
} else {
pValue = paste0("p=", sprintf("%.3f", pValue))
}
这段代码格式化p值。如果p值小于0.001，则将其表示为"p<0.001"；否则，将其格式化为小数点后三位并转换为字符串。
fit = survfit(Surv(futime, fustat) ~ group, data = rt)
这行代码使用survfit函数根据group变量分组拟合生存曲线。Surv(futime, fustat)创建生存对象，~ group表示按group变量分组，data = rt指定数据来源。fit对象包含生存曲线拟合的结果。
总结来说，这段代码的目的是读取输入数据文件，将样品分为高表达组和低表达组，计算不同组间的生存差异的p值，并拟合生存曲线。
绘制生存曲线
这段代码使用ggsurvplot函数绘制生存曲线，并将结果保存到PDF文件中。以下是逐行解释：
# 使用ggsurvplot绘制生存曲线
surPlot = ggsurvplot(fit,
data = rt,
conf.int = TRUE, # 显示置信区间
pval = pValue, # 显示p值
pval.size = 5, # p值字体大小
legend.labs = c("High", "Low"), # 图例标签
legend.title = var, # 图例标题
xlab = "Time (years)", # x轴标签
break.time.by = 1, # x轴时间间隔
risk.table.title = "", # 风险表标题为空
palette = c("red", "blue"), # 设置图形颜色
risk.table = TRUE, # 显示风险表
risk.table.height = .25) # 风险表高度
这段代码使用ggsurvplot函数创建生存曲线图。参数解释如下：
fit: 使用survfit函数拟合的生存曲线对象。data = rt: 用于绘图的数据。conf.int = TRUE: 显示置信区间。pval = pValue: 显示p值。pval.size = 5: 设置p值的字体大小。legend.labs = c("High", "Low"): 设置图例标签为High和Low。legend.title = var: 设置图例标题为变量名var（即MATN3）。xlab = "Time (years)": 设置x轴标签为“Time (years)”。break.time.by = 1: 设置x轴时间间隔为1年。risk.table.title = "": 设置风险表标题为空。palette = c("red", "blue"): 设置生存曲线的颜色，高表达组为红色，低表达组为蓝色。risk.table = TRUE: 显示风险表。risk.table.height = .25: 设置风险表的高度为0.25。
pdf(file = outFile, onefile = FALSE, width = 6, height = 5)
这行代码设置输出设备为PDF文件，文件名为"survival.pdf"，并指定PDF的宽度为6英寸，高度为5英寸。
print(surPlot)
这行代码将绘制的生存曲线图输出到PDF文件中。
dev.off()
这行代码关闭PDF设备，完成图形的保存。
总结来说，这段代码的目的是使用ggsurvplot函数绘制生存曲线，包括置信区间、p值和风险表，并将图形保存为PDF文件。
现在：
长按扫码关注：科研生信充电宝
5元赞赏本文，即喜欢作者~
即可直接解锁：
《R语言小白绘图系列·第34弹|连续变量生存曲线》哦~
您将获得示例输入文件和完整代码
看到这里你还不心动吗？
赶紧关注、转发、点赞、分享，领取你的专属福利吧~
好啦，以上就是今天推文的全部内容啦！
今天的分享就到这里啦~
另外，本公众号建立了R语言和Python学习分享群，联系后台小编，拉您入群。
版权声明：本文内容由互联网用户自发贡献，版权归作者所有，本公众号不拥有所有权，也不承担相关法律责任。
如果您发现本公众号中有涉嫌抄袭的内容，欢迎发送邮件至：kysxcdb@163.com 进行举报，一经查实，本公众号将立刻删除涉嫌侵权内容。
