library(hdrm)
library(magrittr)
library(kableExtra)

cvfit <- vector('list', 8)

# # Elastic net
# cvfit[[1]] <- cv.glmnet(X, y)
# cvfit[[2]] <- cv.glmnet(X, y, alpha=0.75)
# cvfit[[3]] <- cv.glmnet(X, y, alpha=0.5)
# cvfit[[4]] <- cv.glmnet(X, y, alpha=0.25)
# 
# # Mnet
# cvfit[[5]] <- cv.ncvreg(X, y)
# cvfit[[6]] <- cv.ncvreg(X, y, alpha=0.75)
# cvfit[[7]] <- cv.ncvreg(X, y, alpha=0.5)
# cvfit[[8]] <- cv.ncvreg(X, y, alpha=0.25)

tab <- Tab4.1()

a <- stringr::str_replace(rownames(tab), '.*-', '')
rownames(tab) <- glue::glue('$\\alpha = {a}$')
colnames(tab) <- c('$R^2$', 'Size')
kbl(tab, 'latex', digits=c(2,0), escape=FALSE, booktabs=TRUE) %>%
  kable_styling(position = 'center') %>%
  pack_rows("Elastic Net", 1, 4) %>%
  pack_rows("Mnet", 5, 8)

Fig4.4()

x <- read_data(brca1)$X %>% ncvreg::std()
# Approximate average absolute correlation
aaac <- function(x, s = 10000) {
  n <- nrow(x)
  p <- ncol(x)
  rho <- numeric(s)
  for (i in 1:s) {
    idx <- sample(1:p, 2)
    rho[i] <- crossprod(x[, idx[1]], x[, idx[2]])/n
  }
  mean(abs(rho))
}
avg_cor <- aaac(x)

# Var explained
pca <- prcomp(x)
var_exp <- summary(pca)$importance[3,5]

# Effective rank
eigenvalues <- pca$sdev^2
p <- eigenvalues / sum(eigenvalues)
effective_rank <- exp(-sum(p * log(p)))

tab <- Tab4.2()

a <- stringr::str_replace(rownames(tab), '.*-', '')
rownames(tab) <- glue::glue('$\\alpha = {a}$')
colnames(tab) <- c('$R^2$', 'Size')
kbl(tab, 'latex', digits=c(2,0), escape=FALSE, booktabs=TRUE) %>%
  kable_styling(position = 'center') %>%
  pack_rows("Elastic Net", 1, 4) %>%
  pack_rows("Mnet", 5, 8)

x <- read_data(Scheetz2006)$X %>% ncvreg::std()
avg_cor <- aaac(x)

# Var explained
pca <- prcomp(x)
var_exp <- summary(pca)$importance[3,5]

# Effective rank
eigenvalues <- pca$sdev^2
p <- eigenvalues / sum(eigenvalues)
effective_rank <- exp(-sum(p * log(p)))