ROC surfaces in R If my response variable has 2 levels, 0 and 1, I can use a ROC curve to assess the accuracy of my predictive model.  But what if my response variable has 3 levels, -1, 0, and 1?  Is there a way to plot "ROC surfaces" in R?
I started trying to code this up myself, but I couldn't figure out how to calculate and plot the surface.
library(rpart)
library(randomForest)

y <- as.factor(c(-1,-1,-1,0,0,0,1,1,1,1,1,0,0,0,1,1,1,-1,0,1,0,1,0,-1,0,-1,0,-1,1,0))
x1 <- rnorm(1:30)
x2 <- rnorm(1:30)
x3 <- rnorm(1:30)

model <- randomForest(y~x1+x2+x3)
model
y1 <- predict(model,type='prob')

DF <- data.frame(y=as.numeric(as.character(y)),'negone'=y1[,'-1'],'zero'=y1[,'0'],'one'=y1[,'1'])

DF$FPrOne <- 
DF$FPrNegOne <- 
DF$TPR <- 

library(rgl)
plot3d(DF$FPrOne,DF$FPrNegOne,DF$TPR,type=1)

 A: Some idea to stuff it into 2D in 3D is a surface defined by true positive rate of each class (He, Xin, B.D. Gallas, and E.C. Frey. 2009.), so I'll go this way. 
Let's say the g variable holds the vote proportions:
g<-model$votes;

First you need two thresholds, one for each of some two classes (let's say for 1 and -1) -- I'll call them qP and qN, which will create the surface:
qP<-qN<-seq(0,1,30)

Now, let's calculate a true positive rate for positive class (and negative, which is mostly the same):
tpP<-function(qP,qN,g) sum(g[,"1"]>qP & y=="1")/sum(y=="1")
tpN<-function(qP,qN,g) sum(g[,"-1"]>qN & y=="-1")/sum(y=="-1")

The zero class gets the rest of votes:
tpZ<-function(qP,qN,g) sum(!(q[,"-1"]>qN) & !(q[,"1"]>qP) & y=="0")/sum(y=="0")

Now we just have to apply it over the whole qPxqN, easiest way is to use outer on Vectorized functions (note that this is not the effective method, though I don't see any obvious optimizations which will not degenerate readability):
X<-outer(qP,qN,Vectorize(tpP,vec=c("qP","qN")),g)
   outer(qP,qN,Vectorize(tpN,vec=c("qP","qN")),g)->Y
Z<-outer(qP,qN,Vectorize(tpZ,vec=c("qP","qN")),g)

Now all points (X[[i]],Y[[i]],Z[[i]]) should lay on the surface; the matrix format is suitable to run rgl.surface:
rgl.surface(X,Z,Y) //Yup, the order is weired


Definitely does not seem too good, but this is what one can expect on a small random data ;-) 
