导语
GUIDE ╲
前面我们介绍了一个对有害同义突变预测的方法PrDSM,可以发现,在对模型的分析中,大量的使用ROC对模型进行评估,今天我们就来介绍一下ROC的相关内容和两种ROC绘图方法:pROC、plotROC、ggROC和ROCR。
ROC介绍
ROC曲线是受试者工作特征曲线 / 接收器操作特性曲线(receiver operating characteristic curve), 是一个反映二元分类器系统在其识别阈值变化时的诊断能力的图形。
ROC曲线是通过绘制真阳性率(TPR)与假阳性率(FPR)在不同阈值设置下的曲线。在机器学习中,真阳性率也被称为灵敏度、回忆率或检出率。假阳性率也称为误报率,可以计算为(1 -特异度)。ROC曲线也可以被认为是决策规则的Type I Error 的函数(当性能仅从总体的一个样本中计算时,它可以被认为是这些量的估计值)。因此ROC曲线是敏感度或召回率作为降噪的函数。一般情况下,如果真阳性率和假阳性率分布已知,可以通过对y轴上的真阳性率和x轴上的假阳性率绘制的累积分布函数(概率分布下的面积,从-∞到判别阈值)来生成ROC曲线,因此ROC图有时被称为敏感性vs(1−特异性)图。
考虑一个两类预测问题(二元分类),其中结果被标记为正(p)或负(n)。一个二元分类器有四种可能的结果。①如果预测的结果是p,实际值也是p,则称为真正(true positive, TP)。②如果预测的结果是p,实际值为n,则称为假阳性(FP)。③当预测结果与实际值均为n时,是真阴性(TN)。④当预测结果为n而实际值为p时,是假阴性(FN)。
下图所示各个指标及计算公式:
最好的预测方法是在ROC空间的左上角或坐标(0,1)处找到一个点,表示100%的敏感性(无假阴性)和100%的特异性(无假阳性)。(0,1)点也被称为完美分类。所以ROC曲线越靠近左上角,说明该方法分类效果越好。最靠近左上角的ROC曲线上的点是分类错误最少的最好阈值,其假正例和假反例总数最少。可以对不同的学习器比较性能。将各个学习器的ROC曲线绘制到同一坐标中,直观地鉴别优劣,靠近左上角的ROC曲所代表的学习器准确性最高。
AUC是衡量学习器优劣的一种性能指标,为ROC曲线下与坐标轴围成的面积。其意义是:①因为是在1x1的方格里求面积,AUC必在0~1之间。②假设阈值以上是阳性,以下是阴性;③若随机抽取一个阳性样本和一个阴性样本,分类器正确判断阳性样本的值高于阴性样本的概率 = AUC 。④AUC值越大的分类器,正确率越高。
R包介绍
01
R包pROC
pROC是一个用于显示、平滑和比较ROC曲线的工具。(部分)曲线下面积AUC(pAUC)可以通过基于U-statistics或bootstrap的统计检验进行比较。可以计算(p)AUC或ROC曲线的置信区间。
代码语言:javascript复制install.packages("pROC")
library(pROC)
data(aSAH)
#该数据集总结了113例动脉瘤性蛛网膜下腔出血的临床和实验室变量。
1.
(1)建立ROC对象并计算AUC
代码语言:javascript复制roc1 <- roc(aSAH$outcome, aSAH$s100b)
print(roc1)
#或
roc(outcome ~ s100b, aSAH)
代码语言:javascript复制#计算部分auc
auc(roc1, partial.auc = c(1, .9))
(2)使ROC曲线平滑
代码语言:javascript复制smooth(roc1)
(3)方差计算
代码语言:javascript复制roc2 <- roc(aSAH$outcome, aSAH$wfns)
roc3 <- roc(aSAH$outcome, aSAH$ndka)
var(roc1) #计算ROC的方差
cov(roc1, roc3) #计算两组ROC的协方差
(4)绘制ROC
①plot绘制ROC
代码语言:javascript复制plot(roc1)
②绘制ROC并计算pROC
代码语言:javascript复制roc4 <- roc(aSAH$outcome,aSAH$s100b,
percent=TRUE,
#percent敏感性、特异性和AUC是否必须用百分数或分数表示
partial.auc=c(100, 90),
#计算部分auc
partial.auc.correct=TRUE,
#是否对AUC校正,以获得最大AUC为1.0和非歧视性AUC 0.5,
partial.auc.focus="sens",
#计算partial.auc根据"sens",即敏感性
ci=TRUE, boot.n=100,
#ci置信区间,boot.n是bootstrap重复数
ci.alpha=0.9, stratified=FALSE,
#stratified表示bootstrap是否分层
# arguments for plot
plot=TRUE, auc.polygon=TRUE,
#auc.polygon是否将area显示为多边形
max.auc.polygon=TRUE, grid=TRUE,
#max.auc.polygon是否将最大可能的区域显示为多边形
#grid是否添加背景网格
print.auc=TRUE, show.thres=TRUE
#print.auc,AUC的数值是否应该打印在图上
)
③添加ROC曲线
代码语言:javascript复制roc5 <- roc(aSAH$outcome, aSAH$wfns,
plot=TRUE, add=TRUE, percent=roc4$percent)
#在上述ROC绘图基础上再绘制
#add是否将其他ROC曲线将被添加到现有的plot中
2. 绘制置信区间
(1)计算置信区间
代码语言:javascript复制#ROC曲线的坐标系
coords(roc1, "best", ret=c("threshold", "specificity", "1-npv"))
coords(roc2, "local maximas", ret=c("threshold", "sens", "spec", "ppv", "npv"))
#第二个参数为要找的坐标。“all”:ROC曲线上的所有点。
#“local maximas”:ROC曲线的局部极大值。“best”:用一个方法确定最佳阈值
#ret返回坐标
ci(roc2) #计算置信区间
(2)绘制置信区间
代码语言:javascript复制sens.ci <- ci.se(roc1, specificities=seq(0, 1, .1))
#ci.se,在特定情况下计算灵敏度的置信区间
plot(sens.ci, type="shape", col="lightblue")
plot(sens.ci, type="bars")
#type,plot的形状
(3)加入ROC曲线
代码语言:javascript复制plot(roc2, add=TRUE) #加入roc2
plot(ci.thresholds(roc2)) #计算特异性和敏感性CI阈值
3. 比较ROC曲线
(1)Bootstrap比较ROC曲线
代码语言:javascript复制roc.test(roc1, roc2, reuse.auc=FALSE, partial.auc=c(1, .8),
#reuse.auc=FALSE(默认值),“roc”对象包含一个“auc”字段
#请在测试中重用这些specifications。
allow.invalid.partial.auc.correct = TRUE,
#当试图校正对角线以下的pAUC时,指示校正是否必须返回NA(带有警告)。
#设置为TRUE则返回一个(可能无效)已更正的AUC。
#这对于在bootstrap操作中避免对低pAUCs造成偏倚尤其有用
partial.auc.focus="se", partial.auc.correct=TRUE)
(2)修改部分参数
代码语言:javascript复制roc.test(roc1, roc2, reuse.auc=FALSE, partial.auc=c(1, .8),
allow.invalid.partial.auc.correct = TRUE,
partial.auc.correct=TRUE,
boot.n=1000, boot.stratified=FALSE)
4. ROC曲线的样本量power计算
计算ROC曲线的样本量、power、显著性水平或最小曲线下面积
(1)一条曲线
(2)两条曲线
(3)限定参数
代码语言:javascript复制power.roc.test(ncases=41, ncontrols=72, sig.level=0.05, power=0.95)
#ncases,观察样本数
#ncontrols,对照样本数
#sig.level,期望的显著性水平(第一类错误的probability)
#power,测试的期望power(第二类错误的1 -probability)
02
R包plotROC
大多数ROC曲线绘图模糊了cutoff 值,限制了多条曲线的解释和比较。plotROC试图通过提供绘图和交互式工具来解决这些缺点。提供可以生成用于web使用的交互式ROC曲线图,以及打印版本的功能。plotROC是基于ggplot2绘图的。
代码语言:javascript复制install.packages("plotROC")
library(plotROC)
1.生成展示绘图的ROC数据
以下一段语句是用于生成展示绘图的ROC数据(在我们自己绘图的时候可用自己的数据)
代码语言:javascript复制D.ex <- rbinom(50, 1, .5)
#rbinom生成一组二项分布随机数,参数为(生成的随机数数量,进行随机试验的次数,二项分布概率)
rocdata <- data.frame(D = c(D.ex, D.ex),
#D是分类标签,必须为0和1。
M = c(rnorm(50, mean = D.ex, sd = .4),
#rnorm是生成50个服从正态分布的随机数,均值是D.ex,方差是0.4
rnorm(50, mean = D.ex, sd = 1)),
#M是连续marker值或预测的类标签
Z = c(rep("A", 50), rep("B", 50))
#z是分组参数,指的是A、B两种分类
)
rocdata
2. ggplot绘制ROC曲线
代码语言:javascript复制ggroc2 <- ggplot(rocdata, aes(m = M, d = D, color = Z)) geom_roc()
3. calc_auc计算ROC曲线下面积
代码语言:javascript复制calc_auc(ggroc2)
4. direct_label添加标签
代码语言:javascript复制direct_label(ggroc2, labels = "Biomarker",
label.angle = 45,
nudge_x=0.05,nudge_y = -0.1)
#labels标签向量直接添加到图中
#label.angle调整标签角度
#nudge_x, nudge_y水平和垂直的调整,以推动标签。
#Size标签大小
5. 更改cutoffs
代码语言:javascript复制ggplot(rocdata, aes(m = M, d = D)) geom_roc(n.cuts = 20)
代码语言:javascript复制ggplot(rocdata, aes(m = M, d = D)) geom_roc(cutoffs.at = c(1.5, 1, .5, 0, -.5))
6. geom_rocci添加ROC曲线的置信区间
(1)置信区间
代码语言:javascript复制ggplot(rocdata, aes(m = M, d = D)) geom_roc() geom_rocci()
(2)置信区间的显著性水平
代码语言:javascript复制ggplot(rocdata, aes(m = M, d = D, color = Z)) geom_roc() geom_rocci(sig.level = .01)
#sig.level置信区间的显著性水平
(3)biomarker处展示置信区间
代码语言:javascript复制ggplot(rocdata, aes(m = M, d = D)) geom_roc(n.cuts = 0)
geom_rocci(ci.at = quantile(rocdata$M, c(.1, .25, .5, .75, .9)),linetype = 1)
#在biomarker处展示置信区间
(4)stat_rocci()展示置信区间
代码语言:javascript复制ggplot(rocdata, aes(m = M, d = D)) geom_roc()
stat_rocci(ci.at = quantile(rocdata$M, c(.1, .3, .5, .7, .9)))
7. plot_interactive_roc生成交互式的ROC曲线
代码语言:javascript复制rocdata <- data.frame(D = c(D.ex, D.ex),
M = c(rnorm(50, mean = D.ex, sd = .4), rnorm(50, mean = D.ex, sd = 1)),
Z = c(rep("A", 50), rep("B", 50)))
rocplot <- ggplot(rocdata, aes(m = M, d = D)) geom_roc()
plot_interactive_roc(rocplot aes(color = Z))
8. stat_roc()计算经验ROC曲线
代码语言:javascript复制D.ex <- rbinom(50, 1, .5)
rocdata <- data.frame(D = c(D.ex, D.ex),
M = c(rnorm(50, mean = D.ex, sd = .4), rnorm(50, mean = D.ex, sd = 1)),
Z = c(rep("A", 50), rep("B", 50)))
ggplot(rocdata, aes(m = M, d = D)) stat_roc()
9. style_roc()在ROC图添加指南和注释
代码语言:javascript复制D.ex <- rbinom(50, 1, .5)
fakedata <- data.frame(M1 = rnorm(50, mean = D.ex),
D = D.ex)
ggplot(fakedata, aes(m = M1, d = D)) geom_roc()
style_roc(xlab = "1 - Specificity",theme = theme_grey)
03
R包ggROC
代码语言:javascript复制install.packages("ggROC")
library(ggROC)
data(roc)
直接输出PDF文件
代码语言:javascript复制ggroc(roc, 0.01,"green",sp =19,output="F:/roc.pdf")
04
R包ROCR
代码语言:javascript复制install.packages("ROCR")
library(ROCR)
data(ROCR.simple)
pred <- prediction( ROCR.simple$predictions, ROCR.simple$labels)
perf <- performance(pred,"tpr","fpr")
plot(perf)
小编总结
ROC曲线在生信分析中是一个非常常见的分析工具,一般用在我们对自己构建的方法模型的进行验证分析的时候,有一点要注意的是,前提需要有金标准做对照。今天介绍了四种ROC曲线的绘制方法,大家可以试试哦!