ION606
266 lines
10 KiB
R
266 lines
10 KiB
R
# lab 3: (A) bologna
|
|
|
|
# I NEEDED TO INSTALL FORTRAN BS FOR THIS LMAO
|
|
# sudo pacman -S --needed base-devel gcc-fortran lapack openblas libxml2 curl openssl
|
|
|
|
# packages (install if necessary)
|
|
required_pkgs <- c("dplyr", "ggplot2", "caret", "class", "cluster", "factoextra", "gridExtra")
|
|
for (p in required_pkgs) {
|
|
if (!requireNamespace(p, quietly = TRUE)) {
|
|
install.packages(p, repos = "https://cloud.r-project.org", dependencies=TRUE)
|
|
}
|
|
library(p, character.only = TRUE)
|
|
}
|
|
|
|
# path handling: prefer the uploaded path if present
|
|
uploaded_path <- "/home/ion606/Desktop/Data Analytics/Lab 3/abalone_dataset.csv"
|
|
fallback_path <- "/home/ion606/Desktop/Data Analytics/Lab 3/abalone_dataset.csv"
|
|
data_path <- if (file.exists(uploaded_path)) uploaded_path else fallback_path
|
|
|
|
# read dataset
|
|
abalone.data <- read.csv(data_path, stringsAsFactors = FALSE)
|
|
|
|
# canonicalize column names to predictable lower-case tokens
|
|
names(abalone.data) <- tolower(gsub("[[:space:]]+", ".", names(abalone.data)))
|
|
|
|
# if rings column was named differently, try to find it
|
|
if (!"rings" %in% names(abalone.data)) {
|
|
stop("could not find 'rings' column in dataset. column names found: ", paste(names(abalone.data), collapse = ", "))
|
|
}
|
|
|
|
print(names(abalone.data))
|
|
|
|
# old code but I left it here anyways
|
|
# abalone.data$age.group[abalone.data$rings <= 8] <- "young"
|
|
# abalone.data$age.group[abalone.data$rings > 8 & abalone.data$rings <= 11] <- "adult"
|
|
# abalone.data$age.group[abalone.data$rings > 11 & abalone.data$rings <= 35] <- "old"
|
|
# abalone.data$age.group <- factor(abalone.data$age.group, levels = c("young", "adult", "old"))
|
|
|
|
|
|
# new code
|
|
abalone.data$age.group <- cut(abalone.data$rings, breaks = c(0, 8, 11, 35),
|
|
labels = c("young", "adult", "old"),
|
|
right = TRUE, include.lowest = TRUE)
|
|
|
|
if ("sex" %in% names(abalone.data)) {
|
|
abalone.data$sex <- as.factor(abalone.data$sex)
|
|
}
|
|
|
|
# preview
|
|
cat("dataset dims:", dim(abalone.data), "\n")
|
|
cat("columns:", paste(names(abalone.data), collapse = ", "), "\n")
|
|
|
|
expected_num_cols <- c("length", "diameter", "height",
|
|
"whole.weight", "shucked.weight",
|
|
"viscera.weight", "shell.weight")
|
|
|
|
|
|
num_cols_present <- intersect(expected_num_cols, names(abalone.data))
|
|
if (length(num_cols_present) < 3) {
|
|
stop("expected at least three numeric measurement columns; found ", paste(num_cols_present, collapse = ", "))
|
|
}
|
|
|
|
# feature subsets
|
|
features_full <- num_cols_present # numeric
|
|
features_small <- intersect(c("length","diameter","height"), names(abalone.data)) # subset lmao
|
|
|
|
cat("using features (full):", paste(features_full, collapse = ", "), "\n")
|
|
cat("using features (small):", paste(features_small, collapse = ", "), "\n")
|
|
|
|
# data split
|
|
set.seed(123)
|
|
train_index <- createDataPartition(abalone.data$age.group, p = 0.7, list = FALSE)
|
|
train_df <- abalone.data[train_index, , drop = FALSE]
|
|
test_df <- abalone.data[-train_index, , drop = FALSE]
|
|
|
|
# helper to scale numeric features / return matrix + labels
|
|
scale_features <- function(df, feature_names, center = NULL, scale = NULL) {
|
|
mat <- as.data.frame(df[, feature_names, drop = FALSE])
|
|
|
|
# compute center/scale from provided if present (for train/test separation)
|
|
if (is.null(center)) {
|
|
center <- sapply(mat, mean, na.rm = TRUE)
|
|
}
|
|
|
|
if (is.null(scale)) {
|
|
scale <- sapply(mat, sd, na.rm = TRUE)
|
|
# avoid zero sd
|
|
scale[scale == 0] <- 1
|
|
}
|
|
|
|
scaled <- as.data.frame(scale(mat, center = center, scale = scale))
|
|
list(scaled = scaled, center = center, scale = scale)
|
|
}
|
|
|
|
# scale train/test for both feature sets
|
|
train_full_scaled <- scale_features(train_df, features_full)
|
|
test_full_scaled <- scale_features(test_df, features_full, center = train_full_scaled$center, scale = train_full_scaled$scale)
|
|
|
|
train_small_scaled <- scale_features(train_df, features_small)
|
|
test_small_scaled <- scale_features(test_df, features_small, center = train_small_scaled$center, scale = train_small_scaled$scale)
|
|
|
|
# labels for knn
|
|
train_labels <- train_df$age.group
|
|
test_labels <- test_df$age.group
|
|
|
|
# 2 kNN models (initial comparison)
|
|
library(class) # knn()
|
|
|
|
# pick an initial k (odd)
|
|
k_init <- 5
|
|
|
|
knn_predict_and_confmat <- function(train_mat, test_mat, train_labels, test_labels, k) {
|
|
pred <- knn(train = as.matrix(train_mat), test = as.matrix(test_mat), cl = train_labels, k = k)
|
|
cm <- confusionMatrix(pred, test_labels)
|
|
list(pred = pred, confmat = cm)
|
|
}
|
|
|
|
res_full_init <- knn_predict_and_confmat(train_full_scaled$scaled, test_full_scaled$scaled, train_labels, test_labels, k_init)
|
|
res_small_init <- knn_predict_and_confmat(train_small_scaled$scaled, test_small_scaled$scaled, train_labels, test_labels, k_init)
|
|
|
|
cat("\ninitial results (k =", k_init, ")\n")
|
|
cat("full-features accuracy:", res_full_init$confmat$overall["Accuracy"], "\n")
|
|
print(res_full_init$confmat$table)
|
|
cat("small-features accuracy:", res_small_init$confmat$overall["Accuracy"], "\n")
|
|
print(res_small_init$confmat$table)
|
|
|
|
# choose better performing feature subset (by accuracy)
|
|
acc_full <- as.numeric(res_full_init$confmat$overall["Accuracy"])
|
|
acc_small <- as.numeric(res_small_init$confmat$overall["Accuracy"])
|
|
if (acc_full >= acc_small) {
|
|
best_features <- features_full
|
|
best_train_scaled <- train_full_scaled
|
|
best_test_scaled <- test_full_scaled
|
|
chosen_tag <- "full"
|
|
} else {
|
|
best_features <- features_small
|
|
best_train_scaled <- train_small_scaled
|
|
best_test_scaled <- test_small_scaled
|
|
chosen_tag <- "small"
|
|
}
|
|
cat("\nchosen feature subset for tuning:", chosen_tag, "(", paste(best_features, collapse = ", "), ")\n")
|
|
|
|
# optimal k for best performing subset
|
|
k_values <- seq(1, 25, by = 2) # odd ks
|
|
accuracy_by_k <- numeric(length(k_values))
|
|
names(accuracy_by_k) <- k_values
|
|
|
|
for (i in seq_along(k_values)) {
|
|
k <- k_values[i]
|
|
tmp <- knn(train = as.matrix(best_train_scaled$scaled),
|
|
test = as.matrix(best_test_scaled$scaled),
|
|
cl = train_labels,
|
|
k = k)
|
|
cm <- confusionMatrix(tmp, test_labels)
|
|
accuracy_by_k[i] <- as.numeric(cm$overall["Accuracy"])
|
|
}
|
|
|
|
best_k_idx <- which.max(accuracy_by_k)
|
|
best_k <- k_values[best_k_idx]
|
|
cat("\naccuracy_by_k:\n")
|
|
print(round(accuracy_by_k, 4))
|
|
cat("\nbest k:", best_k, "with accuracy", round(accuracy_by_k[best_k_idx], 4), "\n")
|
|
|
|
# final model with best_k
|
|
final_knn <- knn(train = as.matrix(best_train_scaled$scaled),
|
|
test = as.matrix(best_test_scaled$scaled),
|
|
cl = train_labels,
|
|
k = best_k)
|
|
|
|
final_cm <- confusionMatrix(final_knn, test_labels)
|
|
cat("\nfinal confusion matrix (best k):\n")
|
|
print(final_cm)
|
|
|
|
# per-class
|
|
print(final_cm$byClass)
|
|
|
|
# summary
|
|
output_pdf <- "lab3_output.pdf"
|
|
pdf(output_pdf, width = 10, height = 7)
|
|
|
|
# accuracy vs k
|
|
plot(k_values, accuracy_by_k, type = "b", pch = 19, xlab = "k (odd)", ylab = "accuracy", main = paste("k-NN accuracy (chosen subset:", chosen_tag, ")"))
|
|
grid()
|
|
|
|
# Exercise 2: clustering (k-means and pam) using best feature subset
|
|
|
|
# use scaled dataset (all observations) for clustering
|
|
# scale using full population mean/sd
|
|
all_scaled_res <- scale_features(abalone.data, best_features)
|
|
all_scaled <- all_scaled_res$scaled
|
|
|
|
# use fviz_nbclust for optimal K using silhouette
|
|
# < (k = 10 or sqrt(n))
|
|
k_max <- min(10, floor(sqrt(nrow(all_scaled)) * 2))
|
|
|
|
# silhouette
|
|
factoextra::fviz_nbclust(all_scaled, kmeans, method = "silhouette") + ggtitle("fviz_nbclust: silhouette (kmeans)")
|
|
|
|
# show elbow
|
|
factoextra::fviz_nbclust(all_scaled, kmeans, method = "wss") + ggtitle("fviz_nbclust: wss (kmeans)")
|
|
|
|
# pick K by the maximum average silhouette
|
|
avg_sil <- numeric(k_max - 1)
|
|
for (k in 2:k_max) {
|
|
km_tmp <- kmeans(all_scaled, centers = k, nstart = 25)
|
|
sil <- cluster::silhouette(km_tmp$cluster, dist(all_scaled))
|
|
avg_sil[k - 1] <- mean(sil[, 3])
|
|
}
|
|
k_values_clust <- 2:k_max
|
|
best_k_clust <- k_values_clust[which.max(avg_sil)]
|
|
cat("\navg silhouette by k (kmeans):\n")
|
|
print(data.frame(k = k_values_clust, avg_silhouette = round(avg_sil, 4)))
|
|
cat("\nchosen best k for kmeans (max avg silhouette):", best_k_clust, "\n")
|
|
|
|
# run kmeans with best_k_clust and plot silhouette
|
|
km_final <- kmeans(all_scaled, centers = best_k_clust, nstart = 25)
|
|
sil_km <- cluster::silhouette(km_final$cluster, dist(all_scaled))
|
|
factoextra::fviz_silhouette(sil_km) + ggtitle(paste("kmeans silhouette (k=", best_k_clust, ")", sep = ""))
|
|
|
|
# run pam with same range and pick best k for pam by avg silhouette
|
|
avg_sil_pam <- numeric(k_max - 1)
|
|
for (k in 2:k_max) {
|
|
pam_tmp <- cluster::pam(all_scaled, k = k)
|
|
avg_sil_pam[k - 1] <- mean(pam_tmp$silinfo$avg.width)
|
|
}
|
|
|
|
best_k_pam <- k_values_clust[which.max(avg_sil_pam)]
|
|
cat("\navg silhouette by k (pam):\n")
|
|
print(data.frame(k = k_values_clust, avg_silhouette = round(avg_sil_pam, 4)))
|
|
cat("\nchosen best k for pam (max avg silhouette):", best_k_pam, "\n")
|
|
|
|
# run pam for best_k_pam/show silhouette plot
|
|
pam_final <- cluster::pam(all_scaled, k = best_k_pam)
|
|
factoextra::fviz_silhouette(pam_final) + ggtitle(paste("pam silhouette (k=", best_k_pam, ")", sep = ""))
|
|
|
|
# also plot cluster centers (2-d PCA scatter with cluster colors) for kmeans and pam
|
|
pca_res <- prcomp(all_scaled, center = TRUE, scale. = FALSE)
|
|
pcs <- data.frame(pca_res$x[, 1:2])
|
|
pcs$kmeans_cluster <- factor(km_final$cluster)
|
|
pcs$pam_cluster <- factor(pam_final$clustering)
|
|
|
|
# kmeans PCA plot
|
|
ggplot(pcs, aes(x = PC1, y = PC2, color = kmeans_cluster)) + geom_point(alpha = 0.6) + ggtitle(paste("kmeans clusters (k=", best_k_clust, ")", sep = ""))
|
|
|
|
# pam PCA plot
|
|
ggplot(pcs, aes(x = PC1, y = PC2, color = pam_cluster)) + geom_point(alpha = 0.6) + ggtitle(paste("pam clusters (k=", best_k_pam, ")", sep = ""))
|
|
|
|
# kill the pdf device
|
|
dev.off()
|
|
|
|
cat("plots and clustering/kNN visuals saved to", output_pdf, "\n")
|
|
cat("final chosen k for kNN:", best_k, "\n")
|
|
cat("final chosen k for kmeans:", best_k_clust, "\n")
|
|
cat("final chosen k for pam:", best_k_pam, "\n")
|
|
|
|
# yes I am lazy thx
|
|
summary_txt <- paste0(
|
|
"lab 3 results\n\n",
|
|
"chosen feature subset (by initial k=", k_init, "): ", chosen_tag, " features: ", paste(best_features, collapse = ", "), "\n",
|
|
"best k (k-NN tuning): ", best_k, " (accuracy = ", round(accuracy_by_k[best_k_idx], 4), ")\n",
|
|
"kmeans best k (silhouette): ", best_k_clust, "\n",
|
|
"pam best k (silhouette): ", best_k_pam, "\n"
|
|
)
|
|
|
|
writeLines(summary_txt, con = "lab3_summary.txt")
|