1
|
|
# pROC: Tools Receiver operating characteristic (ROC curves) with
|
2
|
|
# (partial) area under the curve, confidence intervals and comparison.
|
3
|
|
# Copyright (C) 2010-2014 Xavier Robin, Alexandre Hainard, Natacha Turck,
|
4
|
|
# Natalia Tiberti, Frédérique Lisacek, Jean-Charles Sanchez
|
5
|
|
# and Markus Müller
|
6
|
|
#
|
7
|
|
# This program is free software: you can redistribute it and/or modify
|
8
|
|
# it under the terms of the GNU General Public License as published by
|
9
|
|
# the Free Software Foundation, either version 3 of the License, or
|
10
|
|
# (at your option) any later version.
|
11
|
|
#
|
12
|
|
# This program is distributed in the hope that it will be useful,
|
13
|
|
# but WITHOUT ANY WARRANTY; without even the implied warranty of
|
14
|
|
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
15
|
|
# GNU General Public License for more details.
|
16
|
|
#
|
17
|
|
# You should have received a copy of the GNU General Public License
|
18
|
|
# along with this program. If not, see <http://www.gnu.org/licenses/>.
|
19
|
|
|
20
|
|
are.paired <- function(...) {
|
21
|
3
|
UseMethod("are.paired")
|
22
|
|
}
|
23
|
|
|
24
|
|
are.paired.auc <- function(roc1, roc2, ...) {
|
25
|
3
|
return(are.paired.roc(roc1, roc2, ...))
|
26
|
|
}
|
27
|
|
|
28
|
|
are.paired.smooth.roc <- function(roc1, roc2, ...) {
|
29
|
3
|
return(are.paired.roc(roc1, roc2, ...))
|
30
|
|
}
|
31
|
|
|
32
|
|
are.paired.roc <- function(roc1, roc2,
|
33
|
|
return.paired.rocs=FALSE,
|
34
|
|
reuse.auc = TRUE, reuse.ci = FALSE, reuse.smooth=TRUE,
|
35
|
|
...) {
|
36
|
|
# check return.paired.rocs
|
37
|
3
|
if (! is.logical(return.paired.rocs) || length(return.paired.rocs) != 1)
|
38
|
0
|
stop("'return.paired.rocs' must be either TRUE or FALSE.")
|
39
|
|
# Recover base ROC curves (not auc or smoothed)
|
40
|
3
|
if ("auc" %in% class(roc1))
|
41
|
3
|
roc1 <- attr(roc1, "roc")
|
42
|
3
|
if ("auc" %in% class(roc2))
|
43
|
3
|
roc2 <- attr(roc2, "roc")
|
44
|
3
|
if ("smooth.roc" %in% class(roc1)) {
|
45
|
3
|
oroc1 <- roc1
|
46
|
3
|
roc1 <- attr(roc1, "roc")
|
47
|
|
}
|
48
|
3
|
if ("smooth.roc" %in% class(roc2)) {
|
49
|
3
|
oroc2 <- roc2
|
50
|
3
|
roc2 <- attr(roc2, "roc")
|
51
|
|
}
|
52
|
|
# Check if the levels are the same. Otherwise it is not paired.
|
53
|
3
|
if (!identical(roc1$levels, roc2$levels))
|
54
|
3
|
return(FALSE)
|
55
|
|
# check if responses of roc 1 and 2 are identical
|
56
|
3
|
if (identical(roc1$response, roc2$response)) {
|
57
|
3
|
retval <- TRUE
|
58
|
3
|
if (return.paired.rocs) {
|
59
|
3
|
attr(retval, "roc1") <- roc1
|
60
|
3
|
attr(retval, "roc2") <- roc2
|
61
|
|
}
|
62
|
3
|
return(retval)
|
63
|
|
}
|
64
|
|
else {
|
65
|
|
# Make sure the difference is not only due to missing values
|
66
|
|
# compare original response (with NAs and response not in levels)
|
67
|
3
|
if (identical(roc1$original.response, roc2$original.response)) {
|
68
|
3
|
retval <- TRUE
|
69
|
3
|
if (! return.paired.rocs)
|
70
|
3
|
return(retval)
|
71
|
|
# else prepare paired ROCs
|
72
|
3
|
idx.exclude <- is.na(roc1$original.predictor) | is.na(roc2$original.predictor) | is.na(roc1$original.response) | ! roc1$original.response %in% roc1$levels
|
73
|
3
|
roc1.paired <- roc(roc1$original.response[!idx.exclude], roc1$original.predictor[!idx.exclude], levels=roc1$levels, percent=roc1$percent, na.rm=FALSE, direction=roc1$direction, auc=FALSE)
|
74
|
3
|
roc2.paired <- roc(roc2$original.response[!idx.exclude], roc2$original.predictor[!idx.exclude], levels=roc2$levels, percent=roc2$percent, na.rm=FALSE, direction=roc2$direction, auc=FALSE)
|
75
|
|
# Re-use auc/ci/smooth for roc1
|
76
|
3
|
if (exists("oroc1") && reuse.smooth) {
|
77
|
3
|
args <- oroc1$smoothing.args
|
78
|
3
|
args$roc <- roc1.paired
|
79
|
3
|
roc1.paired <- do.call("smooth.roc", args)
|
80
|
3
|
roc1.paired$call$roc <- as.name("roc1.paired")
|
81
|
|
}
|
82
|
3
|
if (!is.null(roc1$auc) && reuse.auc) {
|
83
|
3
|
args <- attributes(roc1$auc)
|
84
|
3
|
args$roc <- roc1.paired
|
85
|
3
|
roc1.paired$auc <- do.call("auc.roc", args)
|
86
|
|
}
|
87
|
3
|
if (!is.null(roc1$ci) && reuse.ci) {
|
88
|
0
|
args <- attributes(roc1$ci)
|
89
|
0
|
args$roc <- NULL
|
90
|
0
|
roc1.paired$ci <- do.call(class(roc1$ci)[1], c(roc=list(roc1.paired), args))
|
91
|
|
}
|
92
|
|
# Re-use auc/ci/smooth for roc2
|
93
|
3
|
if (exists("oroc2") && reuse.smooth) {
|
94
|
3
|
args <- oroc2$smoothing.args
|
95
|
3
|
args$roc <- roc2.paired
|
96
|
3
|
roc2.paired <- do.call("smooth.roc", args)
|
97
|
3
|
roc2.paired$call$roc <- as.name("roc2.paired")
|
98
|
|
}
|
99
|
3
|
if (!is.null(roc2$auc) && reuse.auc) {
|
100
|
3
|
args <- attributes(roc2$auc)
|
101
|
3
|
args$roc <- roc2.paired
|
102
|
3
|
roc2.paired$auc <- do.call("auc.roc", args)
|
103
|
|
}
|
104
|
3
|
if (!is.null(roc2$ci) && reuse.ci) {
|
105
|
0
|
args <- attributes(roc2$ci)
|
106
|
0
|
args$roc <- NULL
|
107
|
0
|
roc2.paired$ci <- do.call(class(roc2$ci)[1], c(roc=list(roc2.paired), args))
|
108
|
|
}
|
109
|
|
|
110
|
|
# Attach ROCs and return value
|
111
|
3
|
attr(retval, "roc1") <- roc1.paired
|
112
|
3
|
attr(retval, "roc2") <- roc2.paired
|
113
|
3
|
return(retval)
|
114
|
|
}
|
115
|
|
else {
|
116
|
3
|
return(FALSE)
|
117
|
|
}
|
118
|
|
}
|
119
|
|
}
|