N <- 250
        x1 <- runif(N,0,1)
        x2 <- sample(c(0,1),N,prob=c(.5,.5),replace=T)
        pr <- exp(2-0.2*x1+x2)/(1+exp(2-0.2*x1+x2))
        Y <- rep(0,N)
        for (i in 1:N) {
            Y[i] <- sample(c(0,1),1,prob=c(1-pr[i],pr[i]))
        }
        glm1 <- glm(Y ~ x1+x2,family=binomial)
        phat <- fitted(glm1)
        pc <- sort(unique(phat))
        n2 <- length(pc)
        m2 <- matrix(0,nrow=n2,ncol=2)
        for (i in 1:n2) {
            t1 <- table(factor(phat>pc[i],levels=c(T,F)),Y)
            m2[i,] = c((t1[1,2])/sum(t1[,2]),1-t1[2,1]/sum(t1[,1]))
        }
        plot(m2[,2],m2[,1], type='l',xlim=c(0,1),ylim=c(0,1), xlab='1-specificity',
            ylab='sensitivity')
        abline(0,1)
        AUC <- sum((m2[1:(n2-1),1]+m2[2:n2,1])/2*(-m2[2:n2,2]+m2[1:(n2-1),2]))
        mytitle <- paste("ROC CURVE (AUC=",paste(round(AUC,3)),")",sep="")
        title(mytitle)