生存分析时时P值不显著？生存曲线交叉怎么办？ggscitable包完美复现顶刊（柳叶刀）landmark分析

我们在临床研究中，经常要研究疾病与生存率的关系，cox回归是用得比较常见的模型之一。Cox 比例风险模型依赖于风险随时间变化的假设（PH假设），意思是协变量对结局的影响随着时间变化是固定的。
然而然而实际操作中经常有P值大于0.05，生存曲线交叉，表明，两组生存率没有明显区别，那咱们要不要放弃了呢？
咱们来看看顶刊NEJM新格兰杂志是如何处理，在文章Differential clinical outcomes after 1 year versus 5 years in a randomised comparison of zotarolimus-eluting and sirolimus-eluting coronary stents (the SORT OUT III study): a multicentre, open-label, randomised superiority trial（在随机比较雷帕霉素洗脱支架和西罗莫司洗脱支架1年后和5年后的临床差异结局（SORT OUT III研究）：一项多中心、开放标签、随机 superiority 试验） 中，

这篇文章主要时研究两种支架放置后，观察它们不良事件的发生率， （下图）

作者使用了k-m曲线，发现两个支架的不良事件无明显区别，P=0.40, 但是作者发现在12个的时候，两条曲线还是分得比较开的，有可能在0-12个月这段时间，两条曲线时有区别的。因此作者做了个landmark分析，简单来说就是把时间点在12个月这里分成2段，观察这两段时间不良发生率的区别

结果表明在0-12月这段时间，蓝色先不良事件发生率明显高于红色，OR：2.13,大于12个月后2个支架的不良事件发生率无明显区别，作者的结论也是这样的。

今天咱们通过ggscitable包来复现这篇文章的B图的landmark分析。

先导入R包和数据，使用的时我既往的体检数据，

library(ggscitable)
bc<-read.csv("E:/r/test/demo.csv",sep=',',header=TRUE)
bc$SEX<-as.factor(bc$SEX)
bc$OCCU.NEW<-as.factor(bc$OCCU.NEW)
names(bc)
str(bc)

数据变量很多，我解释几个我等下要用的，HBP：是否发生高血压，结局指标，AGE：年龄，是我们的协变量，SEX：性别，OCCU.NEW这个我也不知道时什么，反正是个2分类变量。公众号回复：体检数据，可以获得数据。

假设我想知道不同OCCU.NEW的分类结局的高血压结局有没有区别，先来个常规的

###常规绘图
scikm(data = bc,x= "OCCU.NEW",y="HBP",time = "AGE",
      legendposition=c(0.25,0.8))

我们可以看到两条曲线胶合在一起，P值大于0.05，看下累积发病率，也是一样的

scikm(data = bc,x= "OCCU.NEW",y="HBP",time = "AGE",cumhaz=T,
      legendposition=c(0.25,0.8))

我们可以看到在30岁这个位置，两条曲线还是分得比较开的，有可能在0-30岁这个区间两个支架的不良事件结局时有区别的，所以我们可以把节点设置为30

####landmark
scikm(data = bc,x= "OCCU.NEW",y="HBP",time = "AGE",cumhaz=T,
      legendposition=c(0.25,0.8),landmark.point = 30,HR=T)

上面HR时通过COX回归算出来的，和下面的基本一样，我们把下面的P值关掉

###关掉P值显示
scikm(data = bc,x= "OCCU.NEW",y="HBP",time = "AGE",cumhaz=T,
      legendposition=c(0.25,0.8),landmark.point = 30,HR=T,pval = F)

右边这个P值和曲线重合了，咱们调整一下

###调整位置
scikm(data = bc,x= "OCCU.NEW",y="HBP",time = "AGE",cumhaz=T,
      legendposition=c(0.25,0.8),landmark.point = 30,HR=T,pval = F,
      hrtxt2.x = 65)

上图这个可以用于发表的图片就基本画好了，scikm函数提供了2种landmark的算法，都是可靠的权威的算法，放心用，使用landmark.type这个参数来改变，默认的时A算法，如果你想换另一种

###方法2
scikm(data = bc,x= "OCCU.NEW",y="HBP",time = "AGE",cumhaz=T,
      legendposition=c(0.25,0.8),landmark.point = 30,HR=T,landmark.type = "B",pval = F)

如果你研究HR的变化，还可以加入协变量

###加入协变量
scikm(data = bc,x= "OCCU.NEW",y="HBP",time = "AGE",cumhaz=T,cov = "SEX",
      legendposition=c(0.25,0.8),landmark.point = 30,HR=T,pval = F,
      hrtxt2.x = 65)

我们可以看到，加入性别这个变量后，P值明显变化，说明这个时对结局影响明显的变量。

未完待续。下面还有视频介绍

生存分析时时P值不显著？生存曲线交叉怎么办？ggscitable包复现顶刊（柳叶刀）landmark分析