Merge pull request #220 from Ziyu-Mu/enet

enet

slides/regularization/rsrc/data-augmentation-images-cat.py → slides/regularization/attic/data-augmentation-images-cat.py

@@ -1,53 +1,53 @@
-# -*- coding: utf-8 -*-
-"""
-Created on Tue Dec 12 17:55:57 2023
-
-@author: chris
-"""
-
-from keras.datasets import cifar10
-import matplotlib.pyplot as plt
-import numpy as np
-from tensorflow.keras.preprocessing.image import ImageDataGenerator
-
-# Load CIFAR-10 dataset
-(x_train, y_train), (_, _) = cifar10.load_data()
-
-# Selecting a random dog image
-# In CIFAR-10, the label for dogs is 5
-dog_indices = np.where(y_train == 3)[0] # cat is 3
-random_index = np.random.choice(dog_indices)
-dog_image = x_train[random_index]
-
-# Data augmentation techniques
-datagen = ImageDataGenerator(
-    rotation_range=20,
-    width_shift_range=0.2,
-    height_shift_range=0.2,
-    shear_range=0.2,
-    zoom_range=0.2,
-    horizontal_flip=True,
-    fill_mode='nearest'
-)
-
-# Preparing the image for augmentation
-dog_image = dog_image.reshape((1,) + dog_image.shape)
-
-# Applying the augmentation and plotting
-fig, axs = plt.subplots(1, 5, figsize=(15, 3))
-axs[0].imshow(dog_image[0])
-axs[0].axis('off')
-axs[0].set_title("Original", fontsize=20)
-
-# Generate 4 augmented images
-i = 1
-for batch in datagen.flow(dog_image, batch_size=1):
-    axs[i].imshow(batch[0].astype('uint8'))
-    axs[i].axis('off')
-    axs[i].set_title(f"Augmented {i}", fontsize=20)
-    i += 1
-    if i > 4:
-        break
-
-plt.tight_layout()
-plt.show()
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Dec 12 17:55:57 2023
+
+@author: chris
+"""
+
+from keras.datasets import cifar10
+import matplotlib.pyplot as plt
+import numpy as np
+from tensorflow.keras.preprocessing.image import ImageDataGenerator
+
+# Load CIFAR-10 dataset
+(x_train, y_train), (_, _) = cifar10.load_data()
+
+# Selecting a random dog image
+# In CIFAR-10, the label for dogs is 5
+dog_indices = np.where(y_train == 3)[0] # cat is 3
+random_index = np.random.choice(dog_indices)
+dog_image = x_train[random_index]
+
+# Data augmentation techniques
+datagen = ImageDataGenerator(
+    rotation_range=20,
+    width_shift_range=0.2,
+    height_shift_range=0.2,
+    shear_range=0.2,
+    zoom_range=0.2,
+    horizontal_flip=True,
+    fill_mode='nearest'
+)
+
+# Preparing the image for augmentation
+dog_image = dog_image.reshape((1,) + dog_image.shape)
+
+# Applying the augmentation and plotting
+fig, axs = plt.subplots(1, 5, figsize=(15, 3))
+axs[0].imshow(dog_image[0])
+axs[0].axis('off')
+axs[0].set_title("Original", fontsize=20)
+
+# Generate 4 augmented images
+i = 1
+for batch in datagen.flow(dog_image, batch_size=1):
+    axs[i].imshow(batch[0].astype('uint8'))
+    axs[i].axis('off')
+    axs[i].set_title(f"Augmented {i}", fontsize=20)
+    i += 1
+    if i > 4:
+        break
+
+plt.tight_layout()
+plt.show()

slides/regularization/rsrc/enet_exp.R

@@ -1,10 +1,24 @@
+# ------------------------------------------------------------------------------
+# enetlogreg
+
+# DATA: generate regression data y = X(n*q ~Normal)*beta(q=500) + eps(n ~Normal)
+#   (1) beta is sparse, only 5 is non-zero out of 500
+#   (2) beta is non-sparse
+#       then calculate R-squared with enet, lasso and ridge regression
+# ------------------------------------------------------------------------------
+
 library(mlr3)
 library(glmnet)
 library(mlr3learners)
 library(mlr3tuning)
 library(mlr3misc)
 library(pracma)
 library(mvtnorm)
+library(future)
+
+set.seed(123)
+
+# DATA -------------------------------------------------------------------------
 
 n_train = 100
 n_test = 10000
@@ -17,30 +31,41 @@ p = 500
 q_seq = c(5, 500)
 x_corr = 0.8
 
-tuner1 = tnr("grid_search", resolution = gs1_grid)
-tuner2 = tnr("grid_search", 
-  param_resolutions = c(alpha = gs2_grid[1], lambda = gs2_grid[2]))
-inner = rsmp ("cv", folds = n_folds)
+# Initialize grid search tuners
+tuner1 = tnr("grid_search", resolution = gs1_grid) # Tuner for lambda only
+tuner2 = tnr("grid_search",                       # Tuner for both alpha and lambda
+             param_resolutions = c(alpha = gs2_grid[1], lambda = gs2_grid[2])
+)
+
+inner = rsmp("cv", folds = n_folds)
 mm = msr("regr.mse")
 
 l1 = lrn("regr.glmnet", alpha = 0, id = "ridge")
 l2 = lrn("regr.glmnet", alpha = 1, id = "lasso")
 l3 = lrn("regr.glmnet", id = "enet")
+
 ss1 = ps(
-  lambda = p_dbl(1e-3, 1e2, logscale = TRUE)        
-)        
+  lambda = p_dbl(1e-3, 1e2, logscale = TRUE) # Log-scaled lambda search space
+)
+
 l1 = auto_tuner(tuner1, l1, inner, mm, search_space = ss1)
 l2 = auto_tuner(tuner1, l2, inner, mm, search_space = ss1)
+
 ss2 = ps(
-  alpha = p_dbl(0, 1), 
-  lambda = p_dbl(1e-3, 1e2, logscale = TRUE)        
-)        
+  alpha = p_dbl(0, 1),
+  lambda = p_dbl(1e-3, 1e2, logscale = TRUE)
+) # Search space
+
 l3 = auto_tuner(tuner2, l3, inner, mm, search_space = ss2)
 
 mylearners = list(l1, l2, l3)
+
 myrsmp = rsmp("holdout", ratio = n_train / n)
-# lrn_order = c("LM", "ridge", "lasso") 
+# lrn_order = c("LM", "ridge", "lasso")
+
+# FUNC -------------------------------------------------------------------------
 
+# Simulate data based on the given parameters and return regression task
 make_simul_data = function(rep_i, q) {
   sigma = x_corr^(0:(p-1))
   sigma = Toeplitz(sigma)
@@ -50,34 +75,39 @@ make_simul_data = function(rep_i, q) {
   y = X %*% theta + eps
   d = as.data.frame(X)
   colnames(d) = sprintf("x%03i", 1:p)
-  d$y = y
+  d$y = y               
   tt = as_task_regr(d, target = "y", id = sprintf("q:%i", q))
   return(tt)
 }
 
+# Function to run benchmarking
 run_bm = function(n_reps) {
   simul_grid = expand.grid(q = q_seq, rep_i = 1:n_reps)
   mytasks = lapply(1:nrow(simul_grid), function(i) {
-    row = simul_grid[i,]                 
+    row = simul_grid[i,]           
     make_simul_data(rep_i = row$rep_i, q = row$q)
   })
   bg = benchmark_grid(mytasks, mylearners, myrsmp)
   bmr = benchmark(bg, store_models = TRUE)
   ba = bmr$aggregate(msr("regr.rsq"))
-  list(bmr = bmr, ba = ba)
+  list(bmr = bmr, ba = ba)# detailed and aggregated benchmark result
 }
 
-future::plan("multicore")
+# DATA -------------------------------------------------------------------------
+
+# Execute benchmarking in parallel using multiple cores
+plan("multicore")
 z = run_bm(n_reps)
 ba = z$ba
 bmr = z$bmr
 
+# Extract and save model coefficients (betas)
 nn = length(bmr$uhashes)
-betas = lapply(1:nn, function(i) {
+betas = lapply(1:nn, function(i){
   at = bmr$resample_results$resample_result[[i]]$learners[[1]]
-  gmod = at$learner$model
-  as.numeric(gmod$beta)
+  gmod = at$learner$model                                 
+  as.numeric(gmod$beta)               
 })
-ba$betas = betas
+ba$betas = betas      
 ba$resample_result = NULL
 save(file = "enet_exp.RData", bmr_aggr = ba)

slides/regularization/rsrc/enet_lasso_ridge_r2.R

@@ -0,0 +1,53 @@
+# ------------------------------------------------------------------------------
+# enetlogreg
+
+# FIG: boxplot of R-squared for elasticnet, lasso and ridge
+# LEFT: linear model with 5 non-Zero coefficients (sparse)
+# RIGHT: linear model with 500 non-Zero coefficients
+# ------------------------------------------------------------------------------
+
+library(ggplot2)
+library(gridExtra)
+load("enet_exp.RData")
+
+# PLOT -------------------------------------------------------------------------
+
+q_values <- sapply(ba$task_id, function(task) {
+  as.numeric(sub("q:(\\d+)", "\\1", task))
+})
+
+performance_df <- as.data.frame(ba)
+performance_df$q <- q_values
+performance_df$learner_id <- as.factor(gsub("\\.tuned", "", performance_df$learner_id))
+
+# linear model with sparse features
+df_5 <- performance_df[performance_df['q']==5,]
+
+p1 <- ggplot(data = df_5, aes(x = regr.rsq, y = learner_id)) +
+  geom_boxplot() +
+  coord_flip() +
+  ylab("") +
+  labs(title="sparse") +
+  xlab("R-squared")+
+  xlim(0.5,1)+
+  theme_minimal(base_size = 10) +
+  theme(legend.position="none",
+        axis.title.x=element_blank())
+
+# linear model with non-sparse features
+df_500 <- performance_df[performance_df['q']==500,]
+
+p2 <- ggplot(data = df_500, aes(x = regr.rsq, y = learner_id)) +
+  geom_boxplot() +
+  coord_flip() +
+  ylab("") +
+  xlab("R-squared")+
+  labs(title="non-sparse") +
+  xlim(0.5,1)+
+  theme_minimal(base_size = 10) +
+  theme(legend.position="none",
+        axis.title.x=element_blank())
+
+p <- grid.arrange(p1, p2, nrow= 1)
+
+ggsave("../figure/enet_lasso_ridge_r2.png", plot = p, width = 6, height = 2)

slides/regularization/rsrc/make_l1_l2_regupaths_mse.R → slides/regularization/rsrc/l1_l2_regupaths_mse.R

@@ -1,42 +1,42 @@
-library(glmnet)
-library(datasets)
-
-set.seed(42)
-
-# Load mtcars
-data(mtcars)
-x <- as.matrix(mtcars[, -1])  
-y <- mtcars$mpg 
-
-# Normalize data
-x <- scale(x, center = TRUE, scale = TRUE)
-y <- scale(y, center = TRUE, scale = FALSE)
-
-# Lasso with cross-validation
-cvfit_lasso <- cv.glmnet(x, y, alpha = 1)
-
-# Ridge with cross-validation
-cvfit_ridge <- cv.glmnet(x, y, alpha = 0)
-
-# Plot layout
-par(mfrow = c(2, 2), oma = c(0, 0, 0, 0), mar = c(5, 4, 4, 2) + 0.1)
-par(cex.main = 1.7, cex.axis = 1.5, cex.lab = 1.5, lwd=1.2)
-
-# Plot regularization path for Lasso
-plot(cvfit_lasso$glmnet.fit, xvar = "lambda", label = TRUE, main = "")
-title(main = "lasso coefficients path", line = 2.5)
-
-# Plot regularization path for Ridge
-plot(cvfit_ridge$glmnet.fit, xvar = "lambda", label = TRUE, main = "")
-title(main = "ridge coefficients path", line = 2.5)
-
-# Plot RMSE vs. Lambda for Lasso
-plot(cvfit_lasso, main = "")
-title(main = "MSE vs. lambda for lasso", line = 2.5)
-
-# Plot RMSE vs. Lambda for Ridge
-plot(cvfit_ridge, main = "")
-title(main = "MSE vs. lambda for ridge", line = 2.5)
-
-# Reset layout to default
-par(mfrow = c(1, 1), oma = c(0, 0, 0, 0), mar = c(5, 4, 4, 2) + 0.1)
+library(glmnet)
+library(datasets)
+
+set.seed(42)
+
+# Load mtcars
+data(mtcars)
+x <- as.matrix(mtcars[, -1])  
+y <- mtcars$mpg 
+
+# Normalize data
+x <- scale(x, center = TRUE, scale = TRUE)
+y <- scale(y, center = TRUE, scale = FALSE)
+
+# Lasso with cross-validation
+cvfit_lasso <- cv.glmnet(x, y, alpha = 1)
+
+# Ridge with cross-validation
+cvfit_ridge <- cv.glmnet(x, y, alpha = 0)
+
+# Plot layout
+par(mfrow = c(2, 2), oma = c(0, 0, 0, 0), mar = c(5, 4, 4, 2) + 0.1)
+par(cex.main = 1.7, cex.axis = 1.5, cex.lab = 1.5, lwd=1.2)
+
+# Plot regularization path for Lasso
+plot(cvfit_lasso$glmnet.fit, xvar = "lambda", label = TRUE, main = "")
+title(main = "lasso coefficients path", line = 2.5)
+
+# Plot regularization path for Ridge
+plot(cvfit_ridge$glmnet.fit, xvar = "lambda", label = TRUE, main = "")
+title(main = "ridge coefficients path", line = 2.5)
+
+# Plot RMSE vs. Lambda for Lasso
+plot(cvfit_lasso, main = "")
+title(main = "MSE vs. lambda for lasso", line = 2.5)
+
+# Plot RMSE vs. Lambda for Ridge
+plot(cvfit_ridge, main = "")
+title(main = "MSE vs. lambda for ridge", line = 2.5)
+
+# Reset layout to default
+par(mfrow = c(1, 1), oma = c(0, 0, 0, 0), mar = c(5, 4, 4, 2) + 0.1)