NHANES数据(复杂调查数据)COX回归亚组交互函数1.8尝鲜版（P for interaction）发布

admin管理员组

文章数量:1589781

写在前面的话，本函数仅适用于NHANES数据的COX回归亚组交互分析，不适用于其他情况，请注意甄别。电子产品，售出不能退换。

在SCI文章中，交互效应表格（通常是表五）能为文章锦上添花，增加文章的信服力，增加结果的可信程度，还能进行数据挖掘。什么是亚组，通常就是特殊类型人群，比如男女，种族等，就是说你的数据放入特殊人群中结果还可靠吗？如果在各个特殊人群中，你的结果很稳定，说明你的结论很可靠。如果亚组的结论和你的数据数据结论相反，你可以拿来做个新论题。还可以比较不同亚组之间有无区别，比如做了心脏支架和没做支架的区别，可以发现很多新思路，易于数据挖掘。

在既往文章《代码+视频–NHANES数据(复杂调查数据)亚组交互函数1.6尝鲜版（P for interaction）发布—用于一键生成交互效应表》中，咱们发布了svy.scitb5函数，反响还不错，基本没啥大问题， 但是svy.scitb5函数只能进行线性回归和cox回归。本次发布能用于cox回归的svyscitb5.coxph函数，就是在之前的svy.scitb5函数加个后缀名，还是很好记的。
下面我来介绍一下，咱们先导入数据

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

我介绍一下数据，SEQN：序列号，RIAGENDR, # 性别， RIDAGEYR, # 年龄，RIDRETH1, # 种族，DMDMARTL, # 婚姻状况，WTINT2YR,WTMEC2YR, # 权重，SDMVPSU, # psu，SDMVSTRA，# strata，LBDGLUSI, #血糖mmol表示，LBDINSI, #胰岛素( pmmol/L)，PHAFSTHR #餐后血糖，LBXGH #糖化血红蛋白，SPXNFEV1, #FEV1：第一秒用力呼气量，SPXNFVC #FVC：用力肺活量，ml（估计肺容量），LBDGLTSI #餐后2小时血糖，factor.FVC是我把肺活量分为了2分类，方便用于测试。

把分类变量转成因子

bc$DMDMARTL<-ifelse(bc$DMDMARTL==1,1,0)
bc$RIAGENDR<-as.factor(bc$RIAGENDR)
bc$RIDRETH1<-as.factor(bc$RIDRETH1)
bc$DMDMARTL<-as.factor(bc$DMDMARTL)
bc$oGTT2<-as.factor(bc$oGTT2)

建立抽样调查函数

bcSvy2<- svydesign(ids = ~ SDMVPSU, strata = ~ SDMVSTRA, weights = ~ WTMEC2YR,
                   nest=TRUE,data = bc)

设定协变量和交互变量，这里要注意一下cov1指协变量，Interaction为交互变量，交互变量

cov1<-c("DMDMARTL","RIAGENDR","LBDINSI","RIDRETH1")	
Interaction<-c("DMDMARTL","RIAGENDR","RIDRETH1")

导入我写的函数

source("E:/r/test/svyscitb5cox18.R")

成功导入后,右上界面应该出现10个函数

函数的主体为

svy.scitb5 (data,x,y,Interaction,cov=NULL,time=NULL,family=NULL,method=NULL,
                     ids=NULL,strata=NULL,weights=NULL,svydstr=NULL)

我来解释一下data是数据，必须数据框形式，x是你研究的目标变量，y是你的结局变量，Interaction是你的分层变量，这个必须是分类变量并转成因子，cov是你的协变量。统一使用family=“svycoxph”，svydstr这里填入调查函数，这里只支持COX回归。如果你使用过我写的svy.scitb5函数，这个函数使用起来完全没有压力，就是多个时间参数而已。

out<-svyscitb5.coxph(data=bc,x="RIDAGEYR",y="factor.FVC",time="time",Interaction=Interaction,
                     cov = cov1,family="svycoxph",svydstr=bcSvy2)

一键生成表格

绘图

plotforest(out)

下面咱们来验证一下算得对不对，就以DMDMARTL婚姻状态这个分层为例子

####亚组
dat1<-subset(bc,bc$DMDMARTL==0)
dat2<-subset(bc,bc$DMDMARTL==1)
bcSvy21<- svydesign(ids = ~ SDMVPSU, strata = ~ SDMVSTRA, weights = ~ WTMEC2YR,
                    nest=TRUE,data = dat1)
bcSvy22<- svydesign(ids = ~ SDMVPSU, strata = ~ SDMVSTRA, weights = ~ WTMEC2YR,
                    nest=TRUE,data = dat2)
fit1<-svycoxph(Surv(time, factor.FVC) ~ RIDAGEYR+RIAGENDR+LBDINSI+RIDRETH1, design = bcSvy21)
fit2<-svycoxph(Surv(time, factor.FVC) ~ RIDAGEYR+RIAGENDR+LBDINSI+RIDRETH1, design = bcSvy22)
summary(fit1)
summary(fit2)

可以看到DMDMARTL=0的时候，HR是0.99，P值是0.203，DMDMARTL=1的时候HR是1.01，P值是0.141，和咱们算出来是非常接近的，所以可靠性是没有问题的。

在既往函数中，我要求协变量是要包含有交互变量的，现在也可以通过contain参数关掉这个设置
咱们重新设置一下交互变量和协变量

cov1<-c("DMDMARTL","RIAGENDR","LBDINSI")	
Interaction<-c("DMDMARTL","RIAGENDR","RIDRETH1")

我们可以看到协变量没有"RIDRETH1"这个变量，没关掉前是做不出来的

out<-svyscitb5.coxph(data=bc,x="RIDAGEYR",y="factor.FVC",time="time",Interaction=Interaction,
                     cov = cov1,family="svycoxph",svydstr=bcSvy2)

关掉后就可以生成结果了

out<-svyscitb5.coxph(data=bc,x="RIDAGEYR",y="factor.FVC",time="time",Interaction=Interaction,
                     cov = cov1,family="svycoxph",svydstr=bcSvy2,contain=F)

接下来介绍一下X是分类变量的情况
假设我研究的变量x是婚姻状态，y变量是肺活量，想了解在这个分层的关系

cov1<-c("RIDAGEYR","RIAGENDR","LBDINSI","RIDRETH1")	
Interaction<-c("RIAGENDR","RIDRETH1")

生成方法都一样

out<-svyscitb5.coxph(data=bc,x="DMDMARTL",y="factor.FVC",time="time",Interaction=Interaction,
                     cov = cov1,family="svycoxph",svydstr=bcSvy2)

上图我要解释一下，做分类变量的时候需要设定一个参考，在DMDMARTL这个分层比较的时候，RIAGENDR.1_DMDMARTL.0是被认定做参考的，什么意思呢？就是当RIAGENDR等于1和DMDMARTL等于0这个亚组的人群是被默认为做参考比较的，其他组都是和它进行比较，分类变量进行亚组交互的时候，分类最好不要太多，要不数据会很大，而且有些层分不到数据就会显示数据缺失NA。

一键生成森林图

plotforest(out)

另外提一下，svyscitb5.coxph 函数method参数，这个参数是用来计算P for interaction值的，它有(“LRT”, “Wald”)两个选项，这个参数平时你不用理也不用填。为什么我要提这个东西呢？因为有时候你的模型数据极端值太小，使用LRT方法算不出来，然后咱们需要改一下，平时一般不用理。
如果有如下报错

需要对函数添加一个method="Wald。
参看文章《NHANES数据(复杂调查数据)亚组交互函数1.4尝鲜版（P for interaction）发布—用于一键生成交互效应表》

如果可以你还不会, 参看我既往的svy.scitb5函数的视频，差不多的

代码+视频--NHANES数据(复杂调查数据)亚组交互函数1.4尝鲜版（P for interaction）发布---用于一键生成交互效应表

需要获取svyscitb5.coxph函数，请看下面这篇文章

NHANES数据(复杂调查数据)COX回归亚组交互函数1.8尝鲜版（P for interaction）发布

本文标签： 数据函数NHANESCOXInteraction