Send output to:
Browser Blue - Charts White
Browser Black/White
CSV
Data X:
0 -11.15 11.67 2.14 1.31 -1.45 -1.45 -1.67 -1.67 -3.57 -1.67 -4.81 -5.54 6.50 3.40 0 3.61 13.67 2.14 -0.33 -1.45 -1.45 -2.04 -2.04 -2.14 -1.30 -5.19 -5.54 1.33 3.40 1 18.36 13.00 3.93 0.98 -1.64 -1.64 -2.04 -2.04 -2.14 0.00 0.00 -6.00 4.67 3.96 1 19.34 6.83 3.93 0.82 -1.45 -1.45 -1.85 -1.85 3.39 0.74 -1.11 -2.31 6.00 3.77 1 -8.81 2.76 0.73 -5.18 1.43 1.43 -4.75 -4.75 4.29 -2.04 3.73 1.31 -4.91 -2.83 1 -4.41 1.72 0.73 -5.89 5.00 5.00 1.36 1.36 5.18 -2.04 2.16 -0.16 -4.91 -2.83 0 -4.07 6.03 0.73 -7.50 3.21 3.21 -1.69 -1.69 6.07 -2.04 1.18 -0.98 -2.98 -2.83 1 0.68 7.41 11.45 -4.46 3.57 3.57 -2.71 -2.71 6.25 -0.56 2.94 4.59 -2.46 -5.83 1 11.86 4.41 5.93 0.00 6.17 6.17 -2.42 -2.42 1.61 -1.48 -0.52 0.57 9.50 2.30 1 11.69 6.61 2.41 0.00 6.17 6.17 -3.06 -3.06 0.18 -2.22 -0.52 0.75 5.17 0.16 1 13.05 5.76 6.67 -0.68 4.50 4.50 -4.03 -4.03 1.43 -2.22 -1.03 1.70 5.33 -1.48 1 13.39 4.24 7.78 4.24 8.50 8.50 -5.16 -5.16 1.07 -2.22 -2.76 1.70 5.33 -0.98 0 9.83 -1.43 4.07 -1.13 -6.10 -6.10 6.15 6.15 -1.61 1.70 -3.75 -22.26 -1.17 -4.26 0 7.12 2.14 7.29 -1.13 -3.05 -3.05 -0.96 -0.96 -3.21 -0.38 -3.75 -15.28 -1.17 -4.26 0 9.83 3.93 7.97 1.13 -5.08 -5.08 -6.35 -6.35 -3.21 -0.75 -1.96 -15.09 -1.17 -3.15 0 11.36 3.93 4.07 3.06 -2.37 -2.37 -4.42 -4.42 -2.32 -1.32 -4.11 -5.85 -1.17 -1.48 0 5.42 -11.19 2.20 -2.78 -1.75 -1.75 0.18 0.18 13.75 -10.20 0.57 -2.60 0.51 -5.65 1 5.76 -12.37 2.71 -2.04 -2.81 -2.81 0.09 0.09 14.29 -7.84 1.89 -4.00 0.17 -5.65 0 7.63 -12.37 0.51 -0.56 -1.40 -1.40 -0.42 -0.42 18.39 -7.84 2.45 0.20 2.20 -2.74 0 6.44 -10.68 -0.17 -5.00 -0.70 -0.70 0.09 0.09 16.96 -0.98 3.96 -4.00 3.05 -2.42 1 -6.95 2.78 4.81 1.30 4.63 4.63 -6.96 -6.96 15.89 -4.31 -1.13 2.50 -2.54 -8.27 0 -4.92 2.78 5.37 1.48 5.56 5.56 -5.29 -5.29 15.18 -4.31 -1.13 -3.14 -1.69 -7.88 1 2.71 0.56 5.00 2.22 5.56 5.56 -5.29 -5.29 9.64 -6.08 -0.94 -1.18 -2.88 -5.96 0 -0.17 -0.19 5.37 2.78 6.48 6.48 -3.33 -3.33 10.71 -3.92 -0.38 -0.20 -3.56 -5.96
Names of X columns:
Y X1 X2 X3 X4 X5 X6 X7 X8 X9 X10 X11 X12 X13 X14
Chart options
R Code
library(brglm) roc.plot <- function (sd, sdc, newplot = TRUE, ...) { sall <- sort(c(sd, sdc)) sens <- 0 specc <- 0 for (i in length(sall):1) { sens <- c(sens, mean(sd >= sall[i], na.rm = T)) specc <- c(specc, mean(sdc >= sall[i], na.rm = T)) } if (newplot) { plot(specc, sens, xlim = c(0, 1), ylim = c(0, 1), type = 'l', xlab = '1-specificity', ylab = 'sensitivity', main = 'ROC plot', ...) abline(0, 1) } else lines(specc, sens, ...) npoints <- length(sens) area <- sum(0.5 * (sens[-1] + sens[-npoints]) * (specc[-1] - specc[-npoints])) lift <- (sens - specc)[-1] cutoff <- sall[lift == max(lift)][1] sensopt <- sens[-1][lift == max(lift)][1] specopt <- 1 - specc[-1][lift == max(lift)][1] list(area = area, cutoff = cutoff, sensopt = sensopt, specopt = specopt) } roc.analysis <- function (object, newdata = NULL, newplot = TRUE, ...) { if (is.null(newdata)) { sd <- object$fitted[object$y == 1] sdc <- object$fitted[object$y == 0] } else { sd <- predict(object, newdata, type = 'response')[newdata$y == 1] sdc <- predict(object, newdata, type = 'response')[newdata$y == 0] } roc.plot(sd, sdc, newplot, ...) } hosmerlem <- function (y, yhat, g = 10) { cutyhat <- cut(yhat, breaks = quantile(yhat, probs = seq(0, 1, 1/g)), include.lowest = T) obs <- xtabs(cbind(1 - y, y) ~ cutyhat) expect <- xtabs(cbind(1 - yhat, yhat) ~ cutyhat) chisq <- sum((obs - expect)^2/expect) P <- 1 - pchisq(chisq, g - 2) c('X^2' = chisq, Df = g - 2, 'P(>Chi)' = P) } x <- as.data.frame(t(y)) r <- brglm(x) summary(r) rc <- summary(r)$coeff try(hm <- hosmerlem(y[1,],r$fitted.values),silent=T) try(hm,silent=T) bitmap(file='test0.png') ra <- roc.analysis(r) dev.off() te <- array(0,dim=c(2,99)) for (i in 1:99) { threshold <- i / 100 numcorr1 <- 0 numfaul1 <- 0 numcorr0 <- 0 numfaul0 <- 0 for (j in 1:length(r$fitted.values)) { if (y[1,j] > 0.99) { if (r$fitted.values[j] >= threshold) numcorr1 = numcorr1 + 1 else numfaul1 = numfaul1 + 1 } else { if (r$fitted.values[j] < threshold) numcorr0 = numcorr0 + 1 else numfaul0 = numfaul0 + 1 } } te[1,i] <- numfaul1 / (numfaul1 + numcorr1) te[2,i] <- numfaul0 / (numfaul0 + numcorr0) } bitmap(file='test1.png') op <- par(mfrow=c(2,2)) plot((1:99)/100,te[1,],xlab='Threshold',ylab='Type I error', main='1 - Specificity') plot((1:99)/100,te[2,],xlab='Threshold',ylab='Type II error', main='1 - Sensitivity') plot(te[1,],te[2,],xlab='Type I error',ylab='Type II error', main='(1-Sens.) vs (1-Spec.)') plot((1:99)/100,te[1,]+te[2,],xlab='Threshold',ylab='Sum of Type I & II error', main='(1-Sens.) + (1-Spec.)') par(op) dev.off() load(file='createtable') a<-table.start() a<-table.row.start(a) a<-table.element(a,'Coefficients of Bias-Reduced Logistic Regression',5,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Variable',header=TRUE) a<-table.element(a,'Parameter',header=TRUE) a<-table.element(a,'S.E.',header=TRUE) a<-table.element(a,'t-stat',header=TRUE) a<-table.element(a,'2-sided p-value',header=TRUE) a<-table.row.end(a) for (i in 1:length(rc[,1])) { a<-table.row.start(a) a<-table.element(a,labels(rc)[[1]][i],header=TRUE) a<-table.element(a,rc[i,1]) a<-table.element(a,rc[i,2]) a<-table.element(a,rc[i,3]) a<-table.element(a,2*(1-pt(abs(rc[i,3]),r$df.residual))) a<-table.row.end(a) } a<-table.end(a) table.save(a,file='mytable.tab') a<-table.start() a<-table.row.start(a) a<-table.element(a,'Summary of Bias-Reduced Logistic Regression',2,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Deviance',1,TRUE) a<-table.element(a,r$deviance) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Penalized deviance',1,TRUE) a<-table.element(a,r$penalized.deviance) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Residual Degrees of Freedom',1,TRUE) a<-table.element(a,r$df.residual) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'ROC Area',1,TRUE) a<-table.element(a,ra$area) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Hosmer–Lemeshow test',2,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Chi-square',1,TRUE) phm <- array('NA',dim=3) for (i in 1:3) { try(phm[i] <- hm[i],silent=T) } a<-table.element(a,phm[1]) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Degrees of Freedom',1,TRUE) a<-table.element(a,phm[2]) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'P(>Chi)',1,TRUE) a<-table.element(a,phm[3]) a<-table.row.end(a) a<-table.end(a) table.save(a,file='mytable1.tab') a<-table.start() a<-table.row.start(a) a<-table.element(a,'Fit of Logistic Regression',4,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Index',1,TRUE) a<-table.element(a,'Actual',1,TRUE) a<-table.element(a,'Fitted',1,TRUE) a<-table.element(a,'Error',1,TRUE) a<-table.row.end(a) for (i in 1:length(r$fitted.values)) { a<-table.row.start(a) a<-table.element(a,i,1,TRUE) a<-table.element(a,y[1,i]) a<-table.element(a,r$fitted.values[i]) a<-table.element(a,y[1,i]-r$fitted.values[i]) a<-table.row.end(a) } a<-table.end(a) table.save(a,file='mytable2.tab') a<-table.start() a<-table.row.start(a) a<-table.element(a,'Type I & II errors for various threshold values',3,TRUE) a<-table.row.end(a) a<-table.row.start(a) a<-table.element(a,'Threshold',1,TRUE) a<-table.element(a,'Type I',1,TRUE) a<-table.element(a,'Type II',1,TRUE) a<-table.row.end(a) for (i in 1:99) { a<-table.row.start(a) a<-table.element(a,i/100,1,TRUE) a<-table.element(a,te[1,i]) a<-table.element(a,te[2,i]) a<-table.row.end(a) } a<-table.end(a) table.save(a,file='mytable3.tab')
Compute
Summary of computational transaction
Raw Input
view raw input (R code)
Raw Output
view raw output of R engine
Computing time
0 seconds
R Server
Big Analytics Cloud Computing Center
Click here to blog (archive) this computation