临床研究中常需要绘制两组或多组患者(如非AKI组和AKI组)的基线特征表。
下图就是临床中常见的基线特征表。
那么在R中怎么快速绘制绘制临床论文中的基线特征表1?
今天介绍一个新的绘制基线表的包——compareGroups
包。
目 录
compareGroups
包可以通过分组变量来创建单变量分析结果的基线特征表,在创建出表格后可以导出各种格式用于报告。
在使用之前先安装和加载R包。
install.packages("compareGroups") # 安装包
library(compareGroups) # 加载包
PREDIMED研究是一项随机、多中心队列研究,共7000余名研究对象,选取其中部分数据进行演示说明。
研究人群在纳入研究前时没有心血管疾病,但是有心血管风险。
将研究人群随机分为3组,每组采用不同的饮食(对照组+低脂饮食、橄榄油+地中海饮食、坚果+地中海饮食),然后随访观察主要不良心血管事件的发生率。
data(predimed) # 加载数据集
View(predimed) # 预览数据集
数据集中变量介绍:
group # 分组变量,不同的饮食方式(Control; MedDiet + Nuts; MedDiet + VOO)
sex # 性别,男性和女性
age # 年龄
smoke # 吸烟,三个水平Never、Current、Former
bmi # 体重指数
waist # 腰围
wth # 腰高比
htn # 因子,是否为高血压,No和Yes
diab # 因子,是否为糖尿病,No和Yes
hyperchol # 因子,是否为高血脂,No和Yes
famhist # 因子,是否有冠心病家族史,No和Yes
hormo # 因子,是否使用激素替代疗法,No和Yes
p14 # MeDiet坚持得分
toevent # 主要结局的随访时间(年)。
event # 因子,是否发生感兴趣结局,No和Yes
看下数据集各变量信息。
str(predimed) # 查看数据集结构
从上面我们可以看到,数据集中的分类变量都显示为因子,并且都添加了标签。
在使用compareGroups
包前需要注意下:
今天用来绘制基线特征表的主要是compareGroups
包的descrTable()
函数。
先不分组,描述下总样本人群。
descrTable( ~ ., data = predimed)
.
,则默认数据集的全部变量进行统计计算。从上图可以看出,基线特征表的结果显示的很清楚,虽然大部分变量都没有缺失值,但是hormo变量存在缺失值,只有5661例患者。
在上面我们简单的统计描述了下总样本人群的基线特征,下面可以添加分组变量分析看看。
数据集中group
为分类变量,表示不同的饮食方式,分为三组。
descrTable(group ~ ., data = predimed)
如上图所示,我们就可以看到基线特征表显示了3组患者的数据,并且自动输出了总体比较的P值,很强大有木有。
上面我们简单统计描述了下总研究人群以及添加分组变量后研究人群的基线特征,但是我们纳入的是该数据集中的所有变量,有时候我们不需要纳入这么多的变量进行统计分析。
下面我们只纳入数据集中的部分变量进行统计分析。
比如我只纳入五个变量进行分析。
descrTable(group ~ age + sex + smoke + waist + hormo, # 左边为分组变量,右边为基线表行变量
data = predimed) # 数据集
如果基线表中纳入的变量较多,不想这么麻烦,也可以选择用移除变量的形式来绘制基线特征表。
比如说可以通过-
号的形式移除下面这四个变量。
descrTable(group ~ . - toevent - event - diab - p14,
data = predimed)
除了选择部分变量进行统计分析外,我们还可以选择亚组人群进行统计分析。
比如说只选取女性进行分析。
descrTable(group ~ age + smoke + waist + hormo + toevent + event + diab + p14,
data = predimed,
subset = sex == "Female") # 只选取女性患者绘制基线表
除了选择亚组人群外,还可以在亚组人群基础上选取特定变量进行研究。
descrTable(group ~ age + sex + smoke + waist + hormo,
data = predimed,
selec = list(hormo = sex == "Female", waist = waist > 20))
# 这里选择女性患者且waist变量中waist>20的患者进行分析
另外基线特征表中的变量可以在公式中出现两次,比如说bmi变量:
descrTable(group ~ age + sex + bmi + bmi + waist + hormo,
data = predimed,
selec = list(bmi.1 = !is.na(hormo)))
输出的基线特征表中会报告两次bmi的统计结果,第一个bmi表示所有患者的bmi结果,第二个bmi是输出hormo变量中无缺失值时研究者的bmi结果。
前面输出的基线表并没有涉及到统计检验的计算,下面来介绍下基线表的统计检验。
在这个包中,默认情况下,连续变量认为是正态分布变量,在生成基线特征表时,将使用均值+标准差
描述连续变量。
如果要指定某一连续变量为非正态分布变量,比如说指定waist
为非正态分布变量,则:
descrTable(group ~ age + smoke + waist + hormo,
data = predimed,
method = c(waist = 2))
如上图所示,在上面的结果中waist
变量被指定为非正态分布的连续变量,数据被描述为中位数+四分位数。
还可以指定多个连续变量为非正态分布变量;
descrTable(group ~ age + smoke + waist + hormo + bmi + toevent,
data = predimed,
method = c(waist = 2, bmi = 2, toevent = 2))
除了上面两种方法外,我们还可以在参数method中将变量设置为NA,表示该变量会自动执行Shapiro-Wilks检验来确定变量是正态分布还是非正态分布。
descrTable(group ~ age + smoke + waist + hormo + bmi + toevent,
data = predimed,
method = c(waist = NA, bmi = NA, toevent = NA))
method中的数字解释:1表示指定连续变量为正态分布;2表示指定连续变量为非正态分布;3表示将连续变量指定为分类变量;NA表示变量自动执行Shapiro-Wilks检验来确定是正态分布还是非正态分布。
如果参数method使用的是NA,还可以使用alpha
参数来指定统计学意义的阈值,默认为0.05。
当分组变量是二分类变量时,可以计算OR值;当分组变量是time-to-event变量时,则计算HR值。
对于大多数分类变量来说,类别水平编码一般为1、2、3等数字,因此我们可以使用ref
参数来指定参考类别水平,设置show.ratio为TRUE表示在基线表中显示OR/HR值。
descrTable(htn ~ age + sex + bmi + smoke,
data = predimed,
ref = 1,
show.ratio = TRUE)
正常来说,一般分类变量的第一水平默认为参考类别水平。
如果想指定某一变量的参考类别水平,比如说性别sex
指定第二水平为参考类别,吸烟smoke
仍指定第一水平为参考类别,则:
descrTable(htn ~ age + sex + bmi + smoke,
data = predimed,
ref = c(smoke = 1, sex = 2),
show.ratio = TRUE)
分类变量除了编码为数字123外,可能类别水平还会编码为yes/no
,这时指定参考水平的参数为ref.no
,默认情况下指定no
类别为参考类别水平。
descrTable(htn ~ age + sex + bmi + hormo + hyperchol,
data = predimed,
ref.no = "NO",
show.ratio = TRUE)
这里的编码不区分大小写,no/No/NO结果是一样的。
连续变量也是可以计算OR或HR值的,默认情况下,连续变量每增加一个单位,计算OR/HR
。
descrTable(htn ~ age + bmi + waist + p14,
data = predimed,
show.ratio = TRUE)
如果想修改每增加单位,比如说连续变量每增加10个单位,计算OR/HR,可以通过fact.ratio
参数来修改:
descrTable(htn ~ age + bmi + waist + p14,
data = predimed,
fact.ratio = c(age = 10,
bmi = 2),
show.ratio = TRUE)
如上所示,年龄修改为每增加10年,bmi修改为每增加2个单位计算OR/HR,其余两个变量还是每增加1个单位。
在计算OR/HR时,默认情况下是选定响应变量(分组变量)的第一水平作为参考类别。
descrTable(htn ~ age + sex + bmi + hyperchol,
data = predimed,
show.ratio = TRUE)
当我们想调整响应变量也就是分组变量的参考类别水平时,比如说指定第二水平作为参考类别,可以修改ref.y
参数来调整。
descrTable(htn ~ age + sex + bmi + hyperchol,
data = predimed,
ref.y = 2, # 指定第二水平为参考类别
show.ratio = TRUE)
在上面的输出的基线特征表中,默认二分类变量、多分类变量的各类别水平的结果都输出来。
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed)
但是对于二分类变量,有时我们只需要显示一个类别结果就可以了,比如说性别,显示女性结果即可,因此,我们可以使用hide参数隐藏某一类别的结果。
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed,
hide = c(sex = "Male"))
而对于二分类变量编码为yes/no的,如果需要隐藏某一类别结果,可以修改hide.no
参数:
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed,
hide = c(sex = "Male"),
hide.no = "no")
这里的编码不区分大小写,no/No/NO结果是一样的。
在输出的基线特征表中,如果需要调整结果中的有效数字位数,可以修改digits
参数。
在前面的表格中,年龄的有效数字位数为2位,性别为1位,想分别修改为4位、3位。
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed,
digits = c(age = 4, sex = 3))
在基线特征表中,分类变量显示结果默认使用频率+百分比
形式显示,如果需要修改显示形式可调整type
参数。
type参数的取值有3个:1表示百分比;3表示病例数;2或NA则两个都显示(默认)。
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed,
type = 1)
如上所示,性别等分类变量只显示百分比等结果。
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed,
type = 3)
如上所示,性别等分类变量只显示病例数等结果
对于分组变量是三分类或多分类变量时,可以修改show.p.mul = TRUE
来计算组间两两比较的p值。
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed,
show.p.mul = TRUE)
在基线表中,有时候需要加入总研究人群,也就是overall
列的统计描述,可以通过修改show.all=TRUE
来显示。
descrTable(group ~ age + sex + bmi + waist + hormo,
data = predimed,
show.all = TRUE)
除了以上输出结果调整外,还可以调整p值、OR/HR值的小数有效位数、显示置信区间、修改表头、按行合并基线表等等,有需要的自行查阅帮助文件。
有时我们需要绘制分层后的基线特征表。
我们可以先绘制一个基线特征表,然后再使用strataTable()函数来添加分层变量,比如说这里我们将性别sex
变量分层。
restab <- descrTable(group ~ age + smoke + bmi + waist + hormo,
data = predimed)
strataTable(restab, "sex")
在绘制好基线特征表后,我们就需要将基线特征表输出来。
R包支持输出的格式有很多(CSV/HTML/LaTeX/PDF/Markdown/Word/Excel
),我们一般常用csv/xls/word/PDF
四种格式,所以就演示怎么输出这几种格式文件。
## 先绘制一个需要输出的基线特征表,并储存在restab中。
restab <- descrTable(group ~ age + smoke + bmi + waist + hormo,
data = predimed)
restab
export2csv(restab, file='table1.csv')
export2xls(restab, file='table1.xlsx')
export2word(restab, file='table1.docx')
export2pdf(restab, file='table1.pdf')
好了,这个包学到这里,基本就OK了。
◆ ◆ ◆ ◆ ◆