Final Project Sports Analytics
getwd()## [1] "/cloud/project"Exploratory Data Analysis
setwd("/cloud/project")
mlb_stats <- read.csv("MLB Players-hittingstats-ss.csv", header = TRUE)
#Data structure
str(mlb_stats)## 'data.frame': 47 obs. of 17 variables:
## $ Player : chr "Trea Turner" "Bo Bichette" "Amed Rosario" "Xander Bogaerts" ...
## $ Pos : chr "SS" "SS" "SS" "SS" ...
## $ Team : chr "LAD" "TOR" "CLE" "BOS" ...
## $ GS : int 160 158 151 148 161 133 148 151 129 138 ...
## $ AB : int 652 652 637 557 630 522 591 593 481 563 ...
## $ H : int 194 189 180 171 170 152 150 145 135 134 ...
## $ X2B : int 39 43 26 38 25 24 31 24 22 31 ...
## $ X3B : int 4 1 9 0 5 1 6 1 5 0 ...
## $ HR : int 21 24 11 15 26 22 20 33 10 31 ...
## $ RBI : int 100 93 71 73 107 64 80 83 55 98 ...
## $ AVG : num 0.298 0.29 0.283 0.307 0.27 0.291 0.254 0.245 0.281 0.238 ...
## $ OBP : num 0.343 0.333 0.312 0.377 0.339 0.366 0.294 0.317 0.327 0.298 ...
## $ SLG : num 0.466 0.469 0.403 0.456 0.449 0.467 0.428 0.455 0.41 0.458 ...
## $ OPS : num 0.809 0.802 0.715 0.833 0.788 0.834 0.722 0.772 0.736 0.756 ...
## $ WAR : num 4.84 3.44 3.95 5.42 5.4 5.55 1.05 4.04 4.5 4.42 ...
## $ Cash2023: chr "$27,272,727 " "$6,100,000 " "$7,800,000 " "$30,000,000 " ...
## $ Age : int 29 24 26 29 28 27 22 28 25 26 ...names(mlb_stats)## [1] "Player" "Pos" "Team" "GS" "AB" "H"
## [7] "X2B" "X3B" "HR" "RBI" "AVG" "OBP"
## [13] "SLG" "OPS" "WAR" "Cash2023" "Age"#Turn the variable Cash2023 numeric
mlb_stats$Cash2023 <- as.numeric(gsub("[$, ]", "", mlb_stats$Cash2023))
str(mlb_stats)## 'data.frame': 47 obs. of 17 variables:
## $ Player : chr "Trea Turner" "Bo Bichette" "Amed Rosario" "Xander Bogaerts" ...
## $ Pos : chr "SS" "SS" "SS" "SS" ...
## $ Team : chr "LAD" "TOR" "CLE" "BOS" ...
## $ GS : int 160 158 151 148 161 133 148 151 129 138 ...
## $ AB : int 652 652 637 557 630 522 591 593 481 563 ...
## $ H : int 194 189 180 171 170 152 150 145 135 134 ...
## $ X2B : int 39 43 26 38 25 24 31 24 22 31 ...
## $ X3B : int 4 1 9 0 5 1 6 1 5 0 ...
## $ HR : int 21 24 11 15 26 22 20 33 10 31 ...
## $ RBI : int 100 93 71 73 107 64 80 83 55 98 ...
## $ AVG : num 0.298 0.29 0.283 0.307 0.27 0.291 0.254 0.245 0.281 0.238 ...
## $ OBP : num 0.343 0.333 0.312 0.377 0.339 0.366 0.294 0.317 0.327 0.298 ...
## $ SLG : num 0.466 0.469 0.403 0.456 0.449 0.467 0.428 0.455 0.41 0.458 ...
## $ OPS : num 0.809 0.802 0.715 0.833 0.788 0.834 0.722 0.772 0.736 0.756 ...
## $ WAR : num 4.84 3.44 3.95 5.42 5.4 5.55 1.05 4.04 4.5 4.42 ...
## $ Cash2023: num 27272727 6100000 7800000 30000000 27000000 ...
## $ Age : int 29 24 26 29 28 27 22 28 25 26 ...#Eliminating non-numeric columns
mlb_stats_num <- mlb_stats[sapply(mlb_stats, is.numeric)]
str(mlb_stats_num)## 'data.frame': 47 obs. of 14 variables:
## $ GS : int 160 158 151 148 161 133 148 151 129 138 ...
## $ AB : int 652 652 637 557 630 522 591 593 481 563 ...
## $ H : int 194 189 180 171 170 152 150 145 135 134 ...
## $ X2B : int 39 43 26 38 25 24 31 24 22 31 ...
## $ X3B : int 4 1 9 0 5 1 6 1 5 0 ...
## $ HR : int 21 24 11 15 26 22 20 33 10 31 ...
## $ RBI : int 100 93 71 73 107 64 80 83 55 98 ...
## $ AVG : num 0.298 0.29 0.283 0.307 0.27 0.291 0.254 0.245 0.281 0.238 ...
## $ OBP : num 0.343 0.333 0.312 0.377 0.339 0.366 0.294 0.317 0.327 0.298 ...
## $ SLG : num 0.466 0.469 0.403 0.456 0.449 0.467 0.428 0.455 0.41 0.458 ...
## $ OPS : num 0.809 0.802 0.715 0.833 0.788 0.834 0.722 0.772 0.736 0.756 ...
## $ WAR : num 4.84 3.44 3.95 5.42 5.4 5.55 1.05 4.04 4.5 4.42 ...
## $ Cash2023: num 27272727 6100000 7800000 30000000 27000000 ...
## $ Age : int 29 24 26 29 28 27 22 28 25 26 ...#mean
average <- mean(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("average:", format(round(average, 0), scientific = FALSE, big.mark = ","))## [1] "average: 6,855,709"#variance
variance <- var(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("variance:", format(round(variance, 0), scientific = FALSE, big.mark = ","))## [1] "variance: 93,571,304,657,024"#standard deviation
stdv <- sd(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("standard deviation:", format(round(stdv, 0), scientific = FALSE, big.mark = ","))## [1] "standard deviation: 9,673,226"#median
median1 <- median(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("median:", format(round(median1, 0), scientific = FALSE, big.mark = ","))## [1] "median: 2,000,000"#min number
minimum <- min(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("minimum:", format(round(minimum, 0), scientific = FALSE, big.mark = ","))## [1] "minimum: 410,326"#max number
maximum <- max(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("maximum:", format(round(maximum, 0), scientific = FALSE, big.mark = ","))## [1] "maximum: 36,000,000"#range of the numbers
range1 <- range(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("range:", format(round(range1, 0), scientific = FALSE, big.mark = ","))## [1] "range: 410,326" "range: 36,000,000"#difference between max and min
difference <- diff(range(mlb_stats_num$Cash2023, na.rm = TRUE))
paste("difference:", format(round(difference, 0), scientific = FALSE, big.mark = ","))## [1] "difference: 35,589,674"#IQR
IQR1 <- IQR(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("IQR:", format(round(IQR1, 0), scientific = FALSE, big.mark = ","))## [1] "IQR: 7,525,400"#quantile
quant <- quantile(mlb_stats_num$Cash2023, na.rm = TRUE)
paste("quantile:", format(round(quant, 0), scientific = FALSE, big.mark = ","))## [1] "quantile: 410,326" "quantile: 724,600" "quantile: 2,000,000"
## [4] "quantile: 8,250,000" "quantile: 36,000,000"names(mlb_stats_num)## [1] "GS" "AB" "H" "X2B" "X3B" "HR"
## [7] "RBI" "AVG" "OBP" "SLG" "OPS" "WAR"
## [13] "Cash2023" "Age"#correlation of the variables
cor(mlb_stats_num)## GS AB H X2B X3B HR
## GS 1.0000000 0.9878832 0.9387030 0.8804607 0.366752978 0.7339877
## AB 0.9878832 1.0000000 0.9741364 0.9156261 0.389106791 0.7840134
## H 0.9387030 0.9741364 1.0000000 0.9356879 0.375379829 0.7724508
## X2B 0.8804607 0.9156261 0.9356879 1.0000000 0.293089007 0.7399494
## X3B 0.3667530 0.3891068 0.3753798 0.2930890 1.000000000 0.1576920
## HR 0.7339877 0.7840134 0.7724508 0.7399494 0.157692045 1.0000000
## RBI 0.8833421 0.9243508 0.9298298 0.8896238 0.315943566 0.8973054
## AVG 0.2527911 0.3350714 0.4999607 0.4540903 0.115910827 0.2808246
## OBP 0.1430023 0.1935521 0.3283211 0.2978383 0.001986825 0.2055263
## SLG 0.1628489 0.2577958 0.3830350 0.3966754 0.102071214 0.5288822
## OPS 0.1672546 0.2520635 0.3908039 0.3880660 0.069776549 0.4423801
## WAR 0.7585507 0.7838553 0.8122353 0.7383154 0.308012465 0.7154940
## Cash2023 0.4710708 0.5099290 0.5628422 0.4634702 0.030049924 0.6053461
## Age 0.2663626 0.2643532 0.2483782 0.2137952 -0.098017716 0.1988217
## RBI AVG OBP SLG OPS WAR
## GS 0.8833421 0.2527911 0.143002251 0.1628489 0.16725465 0.7585507
## AB 0.9243508 0.3350714 0.193552094 0.2577958 0.25206353 0.7838553
## H 0.9298298 0.4999607 0.328321103 0.3830350 0.39080388 0.8122353
## X2B 0.8896238 0.4540903 0.297838268 0.3966754 0.38806597 0.7383154
## X3B 0.3159436 0.1159108 0.001986825 0.1020712 0.06977655 0.3080125
## HR 0.8973054 0.2808246 0.205526332 0.5288822 0.44238007 0.7154940
## RBI 1.0000000 0.3961871 0.285046646 0.4542913 0.42255309 0.7653890
## AVG 0.3961871 1.0000000 0.807340495 0.7975364 0.86254468 0.4335819
## OBP 0.2850466 0.8073405 1.000000000 0.7032172 0.87380311 0.3843565
## SLG 0.4542913 0.7975364 0.703217214 1.0000000 0.96018521 0.4288262
## OPS 0.4225531 0.8625447 0.873803113 0.9601852 1.00000000 0.4439624
## WAR 0.7653890 0.4335819 0.384356543 0.4288262 0.44396237 1.0000000
## Cash2023 0.5789837 0.3434524 0.373311756 0.3875104 0.41241558 0.6341681
## Age 0.2358260 0.1043517 0.054119374 0.1016015 0.09065536 0.1841738
## Cash2023 Age
## GS 0.47107079 0.26636259
## AB 0.50992895 0.26435320
## H 0.56284220 0.24837822
## X2B 0.46347021 0.21379517
## X3B 0.03004992 -0.09801772
## HR 0.60534606 0.19882166
## RBI 0.57898374 0.23582598
## AVG 0.34345236 0.10435169
## OBP 0.37331176 0.05411937
## SLG 0.38751040 0.10160148
## OPS 0.41241558 0.09065536
## WAR 0.63416813 0.18417379
## Cash2023 1.00000000 0.44225191
## Age 0.44225191 1.00000000If pairs are correlating with each other at 0.7+, that’s multicollinearity, and including both in the same regression could make the coefficients unstable or misleading. In that case, it’s often better to pick just one or two representative predictors rather than throwing all four in together. WAR is a strong single choice since it already tries to summarize overall value in one number.
#box plot chart
options(scipen = 999)
boxplot(mlb_stats_num$Cash2023, main="Boxplot of Salaries", ylab="Price ($)")
#histogram chart
options(scipen = 999)
hist(mlb_stats_num$Cash2023, main = "Histogram of Player Prices", xlab = "Price ($)")
#table
table(mlb_stats_num$Cash2023)##
## 410326 520429 536130 541940 632766 654193 661941 720000
## 1 1 1 1 1 1 1 1
## 720100 722000 723200 724200 725000 727600 730000 734500
## 1 1 1 1 1 1 1 1
## 738600 745750 754900 850000 950000 1800000 2000000 2525000
## 1 1 1 1 1 1 2 1
## 2662000 3000000 5000000 5585000 6000000 6100000 6500000 7000000
## 1 1 1 1 2 1 1 1
## 7800000 8700000 9000000 10000000 10250000 12500000 16000000 22000000
## 1 1 1 1 1 1 1 1
## 27000000 27272727 30000000 35000000 36000000
## 1 1 1 1 1# scatterplot
plot(x = mlb_stats_num$RBI, y = mlb_stats_num$Cash2023,
main = "Scatterplot of RBI vs. Salary",
xlab = "RBI",
ylab = "Price ($)")
#scatter plot against Cash2023
par(mfrow = c(3, 5))
for (col in names(mlb_stats_num)) {
if (col != "Cash2023") {
plot(mlb_stats_num[[col]], mlb_stats_num$Cash2023,
main = paste(col, "vs. Salary"),
xlab = col,
ylab = "Price ($)")
}
}
par(mfrow = c(1, 1)) 
Prediction Model
library(rpart)
library(randomForest)## randomForest 4.7-1.2## Type rfNews() to see new features/changes/bug fixes.- Check for multicollinearity among predictors
pairs(mlb_stats_num[, c("WAR", "OPS", "HR", "RBI", "AVG", "OBP", "Age")])
If two predictors are highly correlated with each other (like AVG and
OBP are), keeping both adds little and can destabilize the model.
- Split into training and test sets
set.seed(123)
n <- nrow(mlb_stats_num)
train_idx <- sample(1:n, size = 0.8 * n)
train <- mlb_stats[train_idx, ]
test <- mlb_stats[-train_idx, ]- Choose the best model
Examine model results using summary function
Best Subset Selection
#install.packages("leaps", repos = "https://cran.r-project.org")
library(leaps)
#Best Subset Selection on the train dataset
best_subset <- regsubsets(Cash2023 ~ H + HR + RBI + WAR + AVG + OBP + SLG + OPS + Age + GS + AB + X2B + X3B,
data = train, nvmax = 10) # max number of predictors to consider
summary(best_subset)## Subset selection object
## Call: regsubsets.formula(Cash2023 ~ H + HR + RBI + WAR + AVG + OBP +
## SLG + OPS + Age + GS + AB + X2B + X3B, data = train, nvmax = 10)
## 13 Variables (and intercept)
## Forced in Forced out
## H FALSE FALSE
## HR FALSE FALSE
## RBI FALSE FALSE
## WAR FALSE FALSE
## AVG FALSE FALSE
## OBP FALSE FALSE
## SLG FALSE FALSE
## OPS FALSE FALSE
## Age FALSE FALSE
## GS FALSE FALSE
## AB FALSE FALSE
## X2B FALSE FALSE
## X3B FALSE FALSE
## 1 subsets of each size up to 10
## Selection Algorithm: exhaustive
## H HR RBI WAR AVG OBP SLG OPS Age GS AB X2B X3B
## 1 ( 1 ) " " " " " " "*" " " " " " " " " " " " " " " " " " "
## 2 ( 1 ) " " " " " " "*" " " " " " " " " "*" " " " " " " " "
## 3 ( 1 ) " " "*" " " "*" " " " " " " " " "*" " " " " " " " "
## 4 ( 1 ) " " "*" " " "*" " " "*" " " " " "*" " " " " " " " "
## 5 ( 1 ) "*" "*" " " "*" " " " " " " " " "*" " " "*" " " " "
## 6 ( 1 ) "*" "*" " " " " " " " " "*" "*" "*" " " "*" " " " "
## 7 ( 1 ) "*" "*" " " "*" " " " " "*" "*" "*" " " "*" " " " "
## 8 ( 1 ) "*" "*" " " "*" " " "*" "*" "*" "*" " " "*" " " " "
## 9 ( 1 ) "*" "*" " " "*" " " "*" "*" "*" "*" " " "*" " " "*"
## 10 ( 1 ) "*" "*" " " "*" " " "*" "*" "*" "*" "*" "*" " " "*"#best values for the linear regression model
results <- summary(best_subset)
results$adjr2 # adjusted R-squared for each model size## [1] 0.3187883 0.4179625 0.4608279 0.4641609 0.4661644 0.4773790 0.4948131
## [8] 0.4884056 0.4858591 0.4688147which.max(results$adjr2) # which size is best## [1] 7coef(best_subset, which.max(results$adjr2)) # variables in that best model## (Intercept) H HR WAR SLG
## -17178396.63 215714.47 862406.66 1382414.00 -132522548.11
## OPS Age AB
## 67932216.23 941931.66 -76017.15Liner regression equation: Salary = -1,717,8396.63 + 215,714.47 * H + 862,406.66 * HR + 1,382,414 * WAR - 132,522,548.11 * SLG + 67,932,216.23 * OPS + 941,931.66 * AGE -76,017.15 * AB
#training the liner regression on the train dataset
lm_best <- lm(Cash2023 ~ H + HR + WAR + SLG + OPS + Age + AB, data = train)
summary(lm_best)##
## Call:
## lm(formula = Cash2023 ~ H + HR + WAR + SLG + OPS + Age + AB,
## data = train)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11046717 -3573723 -107538 3681824 13752118
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -17178397 15141144 -1.135 0.26585
## H 215714 131153 1.645 0.11082
## HR 862407 289454 2.979 0.00579 **
## WAR 1382414 968998 1.427 0.16436
## SLG -132522548 72390068 -1.831 0.07744 .
## OPS 67932216 46099503 1.474 0.15136
## Age 941932 318256 2.960 0.00608 **
## AB -76017 39010 -1.949 0.06107 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6198000 on 29 degrees of freedom
## Multiple R-squared: 0.593, Adjusted R-squared: 0.4948
## F-statistic: 6.037 on 7 and 29 DF, p-value: 0.0002109#R squared result on the test dataset prediction
lm_pred <- predict(lm_best, newdata = test)
R_squared <- 1 - sum((test$Cash2023 - lm_pred)^2) / sum((test$Cash2023 - mean(train$Cash2023))^2)
paste("R-squared:", round(R_squared, 4))## [1] "R-squared: 0.7858"prediction_examples <- data.frame(Player = mlb_stats$Player[as.numeric(rownames(test))],
Actual = test$Cash2023,
Predicted = round(lm_pred, 0),
Difference = round(test$Cash2023 - lm_pred, 0)
)
prediction_examples# RMSE formula
rmse <- function(actual, predicted) sqrt(mean((actual - predicted)^2))
# Generate Predictions
lm_pred <- predict(lm_best, newdata = test)
cat("Linear Regression RMSE: $", rmse(test$Cash2023, lm_pred), "\n")## Linear Regression RMSE: $ 5828064Linear Regression: $5,828,064 average prediction error
Conclusion: The linear regression model’s predictions are, on average, about $5.8 million off from the actual salary.
mape <- function(actual, predicted) mean(abs((actual - predicted) / actual)) * 100
lm_mape <- mape(test$Cash2023, lm_pred)
cat("Linear Regression MAPE:", lm_mape,"%\n")## Linear Regression MAPE: 192.3953 %# Decision Tree Regression Model
dt_model <- rpart(Cash2023 ~ H + HR + RBI + WAR + AVG + OBP + SLG + OPS + Age + GS + AB + X2B + X3B,
data = train)
summary(dt_model)## Call:
## rpart(formula = Cash2023 ~ H + HR + RBI + WAR + AVG + OBP + SLG +
## OPS + Age + GS + AB + X2B + X3B, data = train)
## n= 37
##
## CP nsplit rel error xerror xstd
## 1 0.49471047 0 1.0000000 1.0458792 0.3919677
## 2 0.05948364 1 0.5052895 0.8603924 0.2418243
## 3 0.01000000 2 0.4458059 0.7827374 0.2255773
##
## Variable importance
## AB RBI H GS SLG HR Age OPS X3B WAR
## 23 19 17 16 11 10 3 1 1 1
##
## Node number 1: 37 observations, complexity param=0.4947105
## mean=6225058, MSE=7.398353e+13
## left son=2 (30 obs) right son=3 (7 obs)
## Primary splits:
## AB < 545 to the left, improve=0.4947105, (0 missing)
## RBI < 66.5 to the left, improve=0.4412836, (0 missing)
## Age < 27.5 to the left, improve=0.3793108, (0 missing)
## GS < 140 to the left, improve=0.3692215, (0 missing)
## WAR < 2.64 to the left, improve=0.3685210, (0 missing)
## Surrogate splits:
## RBI < 66.5 to the left, agree=0.973, adj=0.857, (0 split)
## H < 140 to the left, agree=0.946, adj=0.714, (0 split)
## GS < 145.5 to the left, agree=0.946, adj=0.714, (0 split)
## HR < 14.5 to the left, agree=0.892, adj=0.429, (0 split)
## SLG < 0.4455 to the left, agree=0.892, adj=0.429, (0 split)
##
## Node number 2: 30 observations, complexity param=0.05948364
## mean=3302713, MSE=1.229556e+13
## left son=4 (17 obs) right son=5 (13 obs)
## Primary splits:
## Age < 26.5 to the left, improve=0.4414330, (0 missing)
## H < 99.5 to the left, improve=0.1603620, (0 missing)
## WAR < 1.14 to the left, improve=0.1396855, (0 missing)
## X2B < 8.5 to the left, improve=0.1156420, (0 missing)
## X3B < 0.5 to the right, improve=0.1125486, (0 missing)
## Surrogate splits:
## SLG < 0.407 to the right, agree=0.700, adj=0.308, (0 split)
## OPS < 0.703 to the right, agree=0.700, adj=0.308, (0 split)
## X3B < 0.5 to the right, agree=0.700, adj=0.308, (0 split)
## WAR < 0.94 to the left, agree=0.667, adj=0.231, (0 split)
## H < 99.5 to the left, agree=0.633, adj=0.154, (0 split)
##
## Node number 3: 7 observations
## mean=1.874939e+07, MSE=1.449011e+14
##
## Node number 4: 17 observations
## mean=1265419, MSE=2.319694e+12
##
## Node number 5: 13 observations
## mean=5966867, MSE=1.281554e+13printcp(dt_model)##
## Regression tree:
## rpart(formula = Cash2023 ~ H + HR + RBI + WAR + AVG + OBP + SLG +
## OPS + Age + GS + AB + X2B + X3B, data = train)
##
## Variables actually used in tree construction:
## [1] AB Age
##
## Root node error: 2737390566136008/37 = 73983528814487
##
## n= 37
##
## CP nsplit rel error xerror xstd
## 1 0.494710 0 1.00000 1.04588 0.39197
## 2 0.059484 1 0.50529 0.86039 0.24182
## 3 0.010000 2 0.44581 0.78274 0.22558Variables actually used: AB and Age only. Even though you gave the tree 7 predictors (H, HR, WAR, OBP, SLG, Age, AB), it only found AB (at-bats) and Age useful enough to split on. This is an important finding that suggests that in your training data, playing time (AB) and player age/experience are doing most of the work in explaining salary, while things like OBP, SLG, and even WAR weren’t chosen as strong enough splitters. This can happen when a few variables are highly correlated with each other (remember, AB correlates with H, HR, etc.) — the tree just picks whichever one creates the cleanest split first, and once it does, the others become redundant.
plotcp(dt_model)
the cross-validated error (X-val Relative Error) is still declining as
tree size increases, all the way to size 3 (the largest tree). It hasn’t
started curving back up, which is the sign of overfitting to normally
watch for.
# Random Forest Regression Model
rf_model <- randomForest(Cash2023 ~ H + HR + RBI + WAR + AVG + OBP + SLG + OPS + Age + GS + AB + X2B + X3B,
data = train)
rf_model##
## Call:
## randomForest(formula = Cash2023 ~ H + HR + RBI + WAR + AVG + OBP + SLG + OPS + Age + GS + AB + X2B + X3B, data = train)
## Type of random forest: regression
## Number of trees: 500
## No. of variables tried at each split: 4
##
## Mean of squared residuals: 51195680143859
## % Var explained: 30.8# We can evaluate the performance of each model using metrics such as RMSE, R-squared, etc.
# Function to calculate RMSE
rmse <- function(actual, predicted) {sqrt(mean((actual - predicted)^2))}
# Generate Predictions
dt_pred <- predict(dt_model, newdata = test)
rf_pred <- predict(rf_model, newdata = test)
cat("Linear Regression RMSE: $", rmse(test$Cash2023, lm_pred), "\n")## Linear Regression RMSE: $ 5828064cat("Decision Tree RMSE: $", rmse(test$Cash2023, dt_pred), "\n")## Decision Tree RMSE: $ 11312926cat("Random Forest RMSE: $", rmse(test$Cash2023, rf_pred), "\n")## Random Forest RMSE: $ 9334539Linear Regression RMSE: $5,828,064 prediction error Decision Tree: $11,312,926 average prediction error Random Forest: $9,335,072 average prediction error
The predictive model with the lowest RMSE is the Linear Regression making it th ebest and most accurate
# Analysis of Variance (ANOVA) for linear regression model
anova(lm_best)# Alternative Linear Regression model for comparison using ANOVA
lm_final <- lm(Cash2023 ~ H + HR + WAR + Age + AB, data = train)
anova(lm_final, lm_best)p = 0.171, not significant so adding OBP and SLG to the simpler model does not significantly improve the fit
# Save the trained model
saveRDS(lm_best, "linear_regression_model.rds")