load(url('http://myweb.uiowa.edu/pbreheny/data/Golub1999.RData'))
source('http://myweb.uiowa.edu/pbreheny/7600/s16/notes/fun.R')
require(glmnet)
require(ncvreg)

# Note: I used ncvreg for lasso models here just to keep the fold assignments
# the same between MCP and lasso; the general code for fitting lasso-penalized
# logistic models is:

cvfit <- cv.glmnet(X, y, family="binomial")


# Deviance ----------------------------------------------------------------

cvfit.lasso <- cv.ncvreg(X, y, family="binomial", penalty="lasso", nfolds=25, seed=9)
cvfit.mcp <- cv.ncvreg(X, y, family="binomial", nfolds=25, seed=9)

plot(cvfit.lasso)
plot(cvfit.mcp)


# Rsq ---------------------------------------------------------------------

plot(cvfit.lasso, type="rsq")
plot(cvfit.mcp, type="rsq")
summary(cvfit.lasso)
summary(cvfit.mcp)

# Misclassification error -------------------------------------------------

plot(cvfit.lasso, type="pred")
plot(cvfit.mcp, type="pred")

cvfit <- cv.glmnet(X, y, family="binomial", type.measure="class") # glmnet version
plot(cvfit, xlim=rev(range(log(cvfit$lambda))), las=1)


# AUC ---------------------------------------------------------------------

cvfit.auc <- cv.glmnet(X, y, family="binomial", type.measure="auc", nfolds=5)
plot(cvfit.auc, xlim=rev(range(log(cvfit.auc$lambda))), las=1, xlab=expression(log(lambda)))
max(cvfit.auc$cvm)

