Project 2: Data Mining, Classification, Prediction
SDS322E
Mining, Classification, Prediction
Name: Viren Velacheri, EID: vv6898
Introduction
Paragraph or two introducing your datasets and variables, why they are interesting to you, etc. See instructions for more information
I am a big sports lover and one of my favorite sports is basketball. I played competitively through middle school and was a part of the school team and still play it today as a side hobby. In basketball, there are various metrics/stats used to evaluate players and teams.
The dataset I chose contains some common variables such as the players and number of games played, but instead of some of the basic statistics that are commonly brought up, this dataset contains more advanced metrics such as win shares, box plus/minus, etc. I found this dataset through the usual google search and the reason this specfic dataset stood out to me, as implied above, is because of the different metrics of data they contain. I am a basketball lover and look forward to exploring this dataset! I expect there to be some obvious associations between some of the metrics, but my goal is to more look at by position and see potentially if these advanced metrics tend to favor certain positions over others. There are over 703 observations/rows and for the binary variable regarding guard or not guard, 289 of them are guards, while the rest or 414 are non guards.
library(tidyverse)
library(stringr)
# read your datasets in here, e.g., with read_csv()
# if your dataset needs tidying, do so here
dataset <- read_csv("project2_dataset.csv")
# removed some of the unnecessary columns
dataset <- dataset %>% select(-"X20") %>% select(-"X25")
# function that I use in creating binary variable for
# determining whether player is a guard or not (1 means
# player, 0 means it is not)
guardOrNot <- function(posVector) {
guardBinary <- vector()
guardCheck <- c("SG", "PG")
for (pos in posVector) {
if (pos %in% guardCheck) {
guardBinary <- c(guardBinary, 1)
} else {
guardBinary <- c(guardBinary, 0)
}
}
return(guardBinary)
}
# add guard or not column to dataset
dataset <- dataset %>% mutate(guardOrNot = guardOrNot(Pos))
# removed Rk column since it is irrelevant to what I am doing
dataset <- dataset %>% select(-"Rk")
# removed any rows with NA values
dataset <- dataset %>% na.omit
# gets number of players that are guards
dataset %>% filter(guardOrNot == 1) %>% count()## # A tibble: 1 x 1
## n
## <int>
## 1 289Cluster Analysis
# your code here
library(cluster)
# clustering code here
cluster_selection <- dataset %>% select(PER, WS, VORP)
sil_width1 <- vector()
for (i in 2:10) {
kms <- kmeans(cluster_selection, centers = i) #compute k-means solution for each k
sil <- silhouette(kms$cluster, dist(cluster_selection)) #get sil widths
sil_width1[i] <- mean(sil[, 3]) #take averages (higher is better)
}
ggplot() + geom_line(aes(x = 1:10, y = sil_width1)) + scale_x_continuous(name = "k",
breaks = 1:10)
ts_pam <- cluster_selection %>% pam(k = 2)
ts_pam$silinfo$avg.width## [1] 0.4551869library(GGally)
cluster_selection <- cluster_selection %>% ungroup %>% mutate(cluster = as.factor(ts_pam$clustering))
ggpairs(cluster_selection, aes(color = cluster))
Contrary to the first graph, the k that gives the best silhouette width is the value 2. Even though the graph shows the best value as k=3, the average silhouette width for this k value is like 0.39. The average silhouette width here though is still only 0.4551869, so the structure present is weak and could be artificial. Since this produced this k value produced the largest average silhouette width though, I picked 2 clusters. The reason I chose the 3 numeric variables that I did was because they seemed the most intriguing of the lot and I didn’t want to pick more as then the plots would be less interpretable.
Based on the clusters, the red cluster has higher WS (win shares) and VORP than the blue cluster. This appears to suggest that the red cluster represents players with higher win shares and VORP, while the blue cluster represents players with lower win shares and VORP. There is not too much overlap between the clusters for any of the variables. There is also some very strong correlation between the variables VORP and WS (0.883). There is also some moderate correlation between the variables WS and PER (0.575) and VORP and PER (0.560). Overall, there is some good correlation between these respective variables.
Dimensionality Reduction with PCA
# PCA code here
cluster_selection <- dataset %>% select(PER, WS, VORP)
pca1 <- princomp(cluster_selection, cor = T)
summary(pca1)## Importance of components:
## Comp.1 Comp.2 Comp.3
## Standard deviation 1.5352403 0.7251372 0.34236424
## Proportion of Variance 0.7856543 0.1752746 0.03907109
## Cumulative Proportion 0.7856543 0.9609289 1.00000000scaledability <- data.frame(scale(cluster_selection))
summary(pca1, loadings = "T")## Importance of components:
## Comp.1 Comp.2 Comp.3
## Standard deviation 1.5352403 0.7251372 0.34236424
## Proportion of Variance 0.7856543 0.1752746 0.03907109
## Cumulative Proportion 0.7856543 0.9609289 1.00000000
##
## Loadings:
## Comp.1 Comp.2 Comp.3
## PER 0.509 0.861
## WS 0.610 -0.344 -0.714
## VORP 0.607 -0.376 0.700corr_matrix <- pca1$scores %>% cor %>% round(10)
matrix_scores <- pca1$scores %>% as.data.frame
matrix_scores %>% ggplot(aes(pca1$scores[, 1], pca1$scores[,
2])) + geom_point() + xlab("PCA1") + ylab("PCA2") + ggtitle("Principal Component Analysis")
library(factoextra)
fviz_pca_biplot(pca1)
I keep the first 2 PCS as they account for over 96% of the variance. PC1 accounts for the highest variability amongst the three variables at 0.78565. It has positive values for all three of the variables, so this means if someone scores high on PC1, they score high on all three of the variables (PER, WS, VORP) and vice versa if someone scores low on PC1, they score low on all three of the variables. PC2 has not nearly as high a variability at 0.175 and has only a positive value for the variable PER. This means if someone scores high on PC2, they score high only on PER and low on the other two variables (WS and VORP). Likewise, if someone scores low on PC2, they score low on PER and high on the other two variables (WS and VORP).
Linear Classifier
# linear classifier code here
linear_fit <- lm(guardOrNot ~ PER + WS + BPM + VORP + FTr + OWS +
DWS + OBPM + DBPM + MP, data = dataset)
prob_reg <- predict(linear_fit)
class_diag(prob_reg, dataset$guardOrNot, positive = 1)## acc sens spec ppv f1 ba auc
## Metrics 0.6373 0.4152 0.7923 0.5825 0.4848 0.6037 0.7119# cross-validation of linear classifier here
set.seed(322)
k = 10
data <- sample_frac(dataset) #randomly order rows
folds <- rep(1:k, length.out = nrow(data)) #create folds
diags <- NULL
i = 1
for (i in 1:k) {
# create training and test sets
train <- data[folds != i, ]
test <- data[folds == i, ]
truth <- test$guardOrNot
# train model
fit <- lm(guardOrNot ~ PER + WS + BPM + VORP + FTr + OWS +
DWS + OBPM + DBPM + MP, data = train, family = "binomial")
# test model
probs <- predict(fit, test) ### GET PREDICTIONS FROM THE TRAINED MODEL ON THE TEST SET HERE
# get performance metrics for each fold
diags <- rbind(diags, class_diag(probs, truth, positive = 1))
}
# average performance metrics across all folds
summarize_all(diags, mean)## acc sens spec ppv f1 ba auc
## 1 0.619 0.37624 0.78265 0.54658 0.44191 0.57943 0.68567The AUC for the linear regression is 0.7119, while for the 10-fold CV it is 0.68567. This appears to imply that the linear regression model is doing a slightly better job of predicting new observations and that the 10-fold CV might be possibly showing slight signs of overfitting since the AUC is just a tad bit smaller.
Non-Parametric Classifier
library(caret)
# non-parametric classifier code here
knn_fit <- knn3(guardOrNot ~ PER + WS + BPM + VORP + FTr + OWS +
DWS + OBPM + DBPM + MP, data = dataset)
prob_knn <- predict(knn_fit, dataset)
class_diag(prob_knn[, 2], dataset$guardOrNot, positive = 1)## acc sens spec ppv f1 ba auc
## Metrics 0.7112 0.5744 0.8068 0.6748 0.6206 0.6906 0.7706# cross-validation of np classifier here
set.seed(322)
k = 10
data <- sample_frac(dataset) #randomly order rows
folds <- rep(1:k, length.out = nrow(data)) #create folds
diags <- NULL
i = 1
for (i in 1:k) {
# create training and test sets
train <- data[folds != i, ]
test <- data[folds == i, ]
truth <- test$guardOrNot
# train model
fit <- knn3(guardOrNot ~ PER + WS + BPM + VORP + FTr + OWS +
DWS + OBPM + DBPM + MP, data = train) ### SPECIFY THE LOGISTIC REGRESSION MODEL FIT TO THE TRAINING SET HERE
# test model
probs <- predict(fit, test) ### GET PREDICTIONS FROM THE TRAINED MODEL ON THE TEST SET HERE
# get performance metrics for each fold
diags <- rbind(diags, class_diag(probs[, 2], truth, positive = 1))
}
# average performance metrics across all folds
summarize_all(diags, mean)## acc sens spec ppv f1 ba auc
## 1 0.542 0.3429 0.67901 0.42121 0.37292 0.51095 0.51547The AUC for KNN was 0.7706, while for the 10-fold CV it had a significantly lower AUC at a value of 0.51547. This means that KNN not only did a better job, but also that there are some major signs of overfitting since the AUC for the 10-fold CV is way smalller.
Regression/Numeric Prediction
# regression model code here
fit <- lm(PER ~ WS + BPM + VORP + FTr + OWS + DWS + OBPM + DBPM +
MP, data = dataset)
yhat <- predict(fit)
cbind(yhat, y = dataset$PER)## yhat y
## 1 11.1794581 14.2
## 2 -8.8394545 -6.5
## 3 15.5824850 15.1
## 4 22.6694666 22.7
## 5 14.4576419 15.7
## 6 14.2260214 15.1
## 7 16.2123126 18.2
## 8 7.1210914 4.2
## 9 12.0843329 12.5
## 10 14.5286465 12.8
## 11 19.3871981 20.1
## 12 19.8465028 21.5
## 13 18.5962947 19.8
## 14 8.4369557 8.9
## 15 9.3676394 9.9
## 16 4.9302979 3.0
## 17 18.6321051 17.2
## 18 27.9926355 29.2
## 19 -0.5829472 0.0
## 20 10.6445790 10.3
## 21 13.4044221 14.6
## 22 10.0475266 12.1
## 23 15.8228849 15.4
## 24 11.6334541 9.8
## 25 13.3874279 11.9
## 26 13.6174885 12.1
## 27 13.0657231 10.2
## 28 14.8956141 15.5
## 29 9.3169306 7.6
## 30 19.4824113 20.3
## 31 10.1466425 10.3
## 32 7.3063471 9.9
## 33 12.5503741 16.3
## 34 16.5142247 17.5
## 35 15.7374194 15.2
## 36 15.8462308 19.3
## 37 13.3671621 12.2
## 38 15.1641909 15.4
## 39 13.4370681 13.4
## 40 12.4512020 11.8
## 41 12.0663467 12.0
## 42 16.6235596 12.9
## 43 9.2645701 10.5
## 44 11.6074078 10.6
## 45 7.5234899 9.7
## 46 20.6709392 22.7
## 47 14.7543688 14.6
## 48 7.5306683 8.0
## 49 7.8611962 8.5
## 50 7.8960327 5.6
## [ reached getOption("max.print") -- omitted 653 rows ]mean((dataset$PER - yhat)^2)## [1] 4.152874# cross-validation of regression model here
set.seed(322)
k = 5
data <- sample_frac(dataset) #randomly order rows
folds <- rep(1:k, length.out = nrow(data)) #create folds
diags <- NULL
i = 1
for (i in 1:k) {
# create training and test sets
train <- data[folds != i, ]
test <- data[folds == i, ]
truth <- test$guardOrNot
# train model
fit <- lm(PER ~ WS + BPM + VORP + FTr + OWS + DWS + OBPM +
DBPM + MP, data = train) ### SPECIFY THE LOGISTIC REGRESSION MODEL FIT TO THE TRAINING SET HERE
# test model
yhat <- predict(fit, newdata = test)
# get performance metrics for each fold
diags <- mean((dataset$PER - yhat)^2)
}
# average performance metrics across all folds
mean(diags)## [1] 69.85123For linear regression the MSE value of 4.152874 was calculated which is very good and implies that the prediction error is not too bad. The 5-fold CV had a bigger MSE at a value of 69.85123 which is definitely bigger, but still not too bad. The prediction error is definitely bigger, but still alright.
Python
library(reticulate)
use_python("/usr/bin/python3", required = F)
ball <- "Ball is"
cat(c(ball, py$life))## Ball is# python code here
life = 'Life!'
print(r.ball, life)## Ball is Life!First, in R, ball is defined as “Ball is”. From there, variable life in python is defined as “Life!”. They are each defined in different environments, but still can be used in each others respective environment thanks to the reticulate package. In R environment, the python variable life is accessed with the use of the py$ syntax and with the use of the cat command, the “Ball is Life!” line is printed. Meanwhile, in the python environment, the ball variable from r is accessed using the r. syntax and so with the use of the python print statement, the same saying is printed out. This shows how we are able to grab variables from the different environments and use them accordingly.
Concluding Remarks
I had a great semester and really learned a lot from this class. I look forward to utilizing all these skills in the future! Thank you!