r语言cox回归批量交互作用 r语言cox回归结果

上期咱们已经介绍了咱们绘制复杂抽样设计数据的基础图形，今天咱们来介绍一下咱们绘制复杂抽样设计cox回归生存曲线(Kaplan-Meier)。

废话不多说咱们先导入数据和R包

library(survey)
pbc<-read.csv("E:/r/test/pbc.csv",sep=',',header=TRUE)

这是一个原发性胆道胆管炎数据，公众号回复：胆管炎数据，可以获得数据，嫌麻烦的也可以在这里下载：

数据我们解释几个等下要用到的变量，age：年龄，trt：治疗方案：1D-青霉烯，2安慰剂，edema：水肿， status: 结局变量0/1/2表示审查、移植、死亡。

咱们先来一波小操作，生成一个预测值，等下好操作，不喜欢可以跳过这部分，对后面的操作没影响

pbc$randomized <- with(pbc, !is.na(trt) & trt>0)
biasmodel <- glm(randomized~age*edema,data=pbc)
pbc$randprob <- fitted(biasmodel)

生成预测值randprob后我们就可以正式分析了，我们先生成一个调查数据

dpbc <- svydesign(id=~1, prob=~randprob, strata=~edema, data=subset(pbc,randomized))

直接使用svykm函数生成预测值

s1 <- svykm(Surv(time,status>0)~sex, design=dpbc)

生成后直接绘图

plot(s1)

这样K-M生存函数图就出来了，简单把，但是这个看起来差一点意思，我们导入jskm包给它修饰一下

library(jskm)
svyjskm(s1)

嗯，这样感觉稍微好一点，我们在开始的时候开可以生成可信区间绘图

s1 <- svykm(Surv(time,status>0)~sex, design=dpbc,se=T)
svyjskm(s1)

这样好像又好一点，还有很多细节可以修改，有兴趣的可以自己研究一下。有些文章中，绘制了生存曲线之后，还会绘制一个生存人数表risk.table，就像下图这样

这个表不是很好做，费了一番功夫还是做出来了

男性组

女性组

这个数据加有权重在里面，和原数据的例数不一样是正常的。如果想比较两条曲线可以有无差异，求出P值

svylogrank(Surv(time,status==2)~sex,design=dpbc)

默认值为logrank检验，rho=1表示广义Wilcoxon检验

svylogrank(Surv(time,status==2)~sex,design=dpbc,rho=1)