库“TableOne”多重比较。逐行计算 p 值

Question

我收到了一位审阅者的评论，他希望在人口统计特征表（表 1）中包含每行特定变量级别的所有 p 值。尽管这个请求在我看来很奇怪（而且不准确），但我还是同意他的建议。

    library(tableone)
## Load data
library(survival); data(pbc)

# drop ID from variable list
vars <- names(pbc)[-1]

## Create Table 1 stratified by trt (can add more stratifying variables)
tableOne <- CreateTableOne(vars = vars, strata = c("trt"), data = pbc, factorVars = c("status","edema","stage"))

print(tableOne, nonnormal = c("bili","chol","copper","alk.phos","trig"), exact = c("status","stage"),  smd = TRUE)

输出：

我需要有变量status、edema 和stage 的每个级别的 p 值，并进行 Bonferroni 校正。我通过documentation 没有成功。 此外，使用卡方比较跨行的样本大小是否正确？

更新：

我不确定我的方法是否正确，但我想与您分享。我为变量 status 为每个层生成了一个虚拟变量，而不是计算 chisq 。

    library(tableone)
    ## Load data
    library(survival); data(pbc)

    d <- pbc[,c("status", "trt")]

    # Convert dummy variables
    d$status.0 <- ifelse(d$status==0, 1,0)
    d$status.1 <- ifelse(d$status==1, 1,0)
    d$status.2 <- ifelse(d$status==2, 1,0)

t <- rbind(    
    chisq.test(d$status.0, d$trt),
    # p-value = 0.7202

    chisq.test(d$status.1, d$trt),
    # p-value = 1

    chisq.test(d$status.2, d$trt)
    #p-value = 0.7818
)
t

用于多重比较的 BONFERRONI ADJ：

p <- t[,"p.value"]

p.adjust(p, method = "bonferroni")

Answer 1

这个问题是很久以前发布的，所以我想你已经回答了审稿人。

我真的不明白为什么只为三个变量计算调整后的 p 值。事实上，调整 p 值取决于进行的比较次数。如果您将p.adjust() 与 3 个 p 值的向量一起使用，则结果不会真正根据所进行的比较量“调整”（您确实做了十几个半！）

我展示了如何提取所有 p 值，以便计算调整后的值。 从包tableOne中提取pValues 有一种方法调用对象属性（先解释），还有两种快速而肮脏的方法（在底部）。

要提取它们，首先我复制您的代码来创建您的 tableOne：

library(tableone)
## Load data
library(survival); data(pbc)

# drop ID from variable list
vars <- names(pbc)[-1]

## Create Table 1 stratified by trt (can add more stratifying variables)
tableOne <- CreateTableOne(vars = vars, strata = c("trt"), data = pbc, factorVars = c("status","edema","stage"))

你可以通过attributes()查看你的“tableOne”对象有什么

attributes(tableOne)

您可以看到一个表格通常有一个用于连续变量和分类变量的表格。你也可以在其中使用attributes()

attributes(tableOne$CatTable)
# you can notice $pValues

现在您知道 pValue 的“位置”了，您可以使用 attr() 提取它们

attr(tableOne$CatTable, "pValues")

与数值变量类似的东西：

attributes(tableOne$ContTable)
# $pValues are there
attr(tableOne$ContTable, "pValues")

您有 Normal 和 NonNormal 变量的 pValue。 正如您之前设置的那样，您可以同时提取它们

mypCont <- attr(tableOne$ContTable, "pValues") # put them in an object
nonnormal = c("bili","chol","copper","alk.phos","trig") # copied from your code

mypCont[rownames(mypCont) %in% c(nonnormal), "pNonNormal"] # extract NonNormal

"%!in%" <- Negate("%in%")
mypCont[rownames(mypCont) %!in% c(nonnormal), "pNonNormal"] # extract Normal

说了这么多，你的 pValues 被提取出来，我认为有两种更方便​​的快速和肮脏的方法来完成同样的事情：

快速而肮脏的方式 A：将 dput() 与打印的 tableOne 一起使用。然后在控制台中搜索 pValues 所在的位置并将它们复制粘贴到脚本中，以将它们存储在一个对象中

快速而肮脏的方式 B：如果您查看 tableOne vignette 有一个“导出”部分，您可以使用 print(tableOne, quote = TRUE) 然后复制并粘贴到电子表格（如 LibreOffice、卓越...）。 然后我会选择带有 pValue 的列，转置它，然后将其返回到 R，以使用 p.adjust() 计算调整后的 p 值并将它们复制回电子表格以供期刊提交