load(url("http://myweb.uiowa.edu/pbreheny/data/bcTCGA.RData"))
require(ncvreg)

# Simulation --------------------------------------------------------------

# Simulation code (14-22) is complicated and not particularly portable;
# I still need to work on how best to share it


# TCGA data: Lasso --------------------------------------------------------

# Slide 24
cvfit <- cv.ncvreg(X, y, penalty="lasso")
fit <- cvfit$fit
plot(fir(fit))
plot(fir(fit), type="EF")
plot(fir(fit), log.l=TRUE)
plot(cvfit)

# Slide 25
FIR <- fir(fit)
max(FIR$S[FIR$FIR < .05])
fit$lambda[which.max(FIR$S[FIR$FIR < .05])]
cvfit$lambda.min
FIR$S[cvfit$min]
FIR$FIR[cvfit$min]


# TCGA data: MCP ----------------------------------------------------------

# Slide 26
cvfit <- cv.ncvreg(X, y, lambda.min=0.1)
fit <- cvfit$fit
plot(fir(fit))
plot(fir(fit), type="EF")
plot(fir(fit), log.l=TRUE)
plot(cvfit)

# Slide 27
FIR <- fir(fit)
max(FIR$S[FIR$FIR < .1])
fit$lambda[which.max(FIR$S[FIR$FIR < .1])]
cvfit$lambda.min
FIR$S[cvfit$min]
FIR$FIR[cvfit$min]