1

我正在尝试使用序数逻辑回归(特别是 MASS 包中的 polr)计算一组生物标志物的优势比和相关 p 值。我一直在使用 epiDisplay 包中的 'ordinal.or.display()' 函数来查看回归结果,但请注意显示的 p 值与我手动计算的值之间存在差异……大约是两倍当我使用正态分布计算它时很大。我是否遗漏了特定于序数逻辑回归的某些内容,或者这是 epiDisplay 函数的问题?

我尝试查看 epiDisplay 包的文档(https://cran.r-project.org/web/packages/epiDisplay/epiDisplay.pdf),但没有找到任何解释 p 值是如何计算的. 非常感谢任何帮助或其他知识!

#the model using polr from MASS
#generating artifical outcome var from mt cars

mtcars <- mtcars
mtcars$outcome <- round(runif(nrow(mtcars))*5)
myMod <-polr(ordered(outcome) ~ factor(am)+ factor(carb)+ wt,
             data = mtcars,
             Hess = TRUE)

summary <- summary(myMod)
(ctable <- coef(summary(myMod)))
p <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2

## combined table: p = 0.624
(ctable <- cbind(ctable, "p value" = p))


ordinal.or.display(myMod) #p value = 0.312

我希望 p 值是相同的——epiDisplay 可能不会乘以 2 吗?

4

1 回答 1

0

模拟数据时请设置随机数种子。否则,我们无法重现您的结果。

set.seed(123)

当您应该使用 t 分布时,您已经使用正态分布来获得 p 值。所以你计算p值的代码应该是:

df <- summary(myMod)$df.residual
p <- pt(abs(ctable[, "t value"]), df=df, lower.tail=FALSE) * 2

这也是ordinal.or.displayepiDisplay 包中的函数进行计算的方式。您可以通过查看大约第 9 行的代码来检查。

library(epiDisplay)
ordinal.or.display

现在检查:

data(mtcars)
mtcars$outcome <- round(runif(nrow(mtcars))*5)

myMod <- polr(ordered(outcome) ~ factor(am)+ factor(carb)+ wt,
              data = mtcars,
              Hess = TRUE)

ctable <- coef(summary(myMod))

df <- summary(myMod)$df.residual
p <- pt(abs(ctable[, "t value"]), df=df, lower.tail=FALSE) * 2

options(digits=2)

(ctable <- cbind(ctable, "p value:" = p))

       Value Std. Error t value p value:
am1   -2.335       1.19  -1.962    0.064
carb2  0.045       0.92   0.049    0.961
carb3  1.840       1.40   1.319    0.202
carb4 -0.618       1.16  -0.534    0.599
carb6  1.721       1.85   0.933    0.362
carb8  2.418       2.10   1.152    0.263
wt    -0.875       0.68  -1.295    0.210

ordinal.or.display(myMod)
#Waiting for profiling to be done...
      Ordinal OR lower95ci upper95ci P value
am1   0.097      0.009     0.971     0.064  
carb2 1.046      0.166     6.477     0.961  
carb3 6.298      0.428     112.827   0.202  
carb4 0.539      0.051     5.065     0.599  
carb6 5.59       0.121     273.826   0.362  
carb8 11.218     0.162     864.225   0.263  
wt    0.417      0.106     1.573     0.210

在我看来是一样的。

于 2020-06-27T01:44:38.310 回答