在对一组数据进行ROC分析之后,如何计算p值?使用相同的统计信息,我看到可以在SPSS中输出p值。示例代码如下:
library(pROC)
data(aSAH)
head(aSAH)
# gos6 outcome gender age wfns s100b ndka
# 29 5 Good Female 42 1 0.13 3.01
# 30 5 Good Female 37 1 0.14 8.54
# 31 5 Good Female 42 1 0.10 8.09
# 32 5 Good Female 27 1 0.04 10.42
# 33 1 Poor Female 42 3 0.13 17.40
# 34 1 Poor Male 48 2 0.10 12.75
(rr <- roc(aSAH$outcome, aSAH$s100b, plot=T))
# Setting levels: control = Good, case = Poor
# Setting direction: controls < cases
#
# Call:
# roc.default(response = aSAH$outcome, predictor = aSAH$s100b, plot = F)
#
# Data: aSAH$s100b in 72 controls (aSAH$outcome Good) < 41 cases (aSAH$outcome Poor).
# Area under the curve: 0.7314
在SPSS中计算出的p值为0.000007,但通过计算出的p值为verification::roc.area()
0.000022546,是否roc.area()
与SPSS的计算方法不一致?
levels(aSAH$outcome) <- c(0, 1)
library(verification)
ra <- roc.area(as.numeric(as.vector(aSAH$outcome)), rr$predictor)
ra$p.value
# [1] 0.00002254601
没有选择来获取p值pROC::roc
,您可以设置选项ci=TRUE
来获取置信区间。pROC::roc
产生不可见的输出,您可以通过将其分配给对象来获取。
library(pROC)
data(aSAH)
rr <- pROC::roc(aSAH$outcome, aSAH$s100b, ci=TRUE)
使用str(rr)
揭示了如何访问ci
:
rr$ci
# 95% CI: 0.6301-0.8326 (DeLong)
因此,您已经有了一个置信区间。
此外,您还可以使用pROC::var
*获得方差,从中可以手动计算标准误差。
(v <- var(rr))
# [1] 0.002668682
b <- rr$auc - .5
se <- sqrt(v)
(se <- sqrt(v))
# [1] 0.05165929
*注意,还有一个bootstrap选项pROC::var(rr, method="bootstrap")
。
这与Stata计算的结果相同,
# . roctab outcome_num s100b, summary
#
# ROC -Asymptotic Normal--
# Obs Area Std. Err. [95% Conf. Interval]
# ------------------------------------------------------------
# 113 0.7314 0.0517 0.63012 0.83262
# .
# . display r(se)
# .05165929
其中的Stata基本参考手册14 -roctab
状态(第2329):
默认情况下,
roctab
使用DeLong,DeLong和Clarke-Pearson(1988)建议的算法和渐近正态置信区间来计算曲线下面积的标准误差。
一旦有了标准误差,我们还可以基于z分布(参考)来计算p值。
z <- (b / se)
2 * pt(-abs(z), df=Inf) ## two-sided test
# [1] 0.000007508474
该p值接近您的SPSS值,因此很可能是使用类似于Stata的算法计算的(比较:IBM SPSS Statistics 24 Algorithms,第888:889页)。
但是, ROC分析的p值的计算可能会引起争议。例如,您在编辑中显示的方法(另请参见下面的第一个链接)是基于Mann-Whitney U统计量的。
在决定哪种方法最适合您的分析之前,您可能需要对主题进行更深入的研究。我在这里为您提供一些阅读建议:
本文收集自互联网,转载请注明来源。
如有侵权,请联系[email protected] 删除。
我来说两句