Use the dataset assigned to you in the LAB Instructions.
Use the depvar indicated for your dataset.
Select 5-8 independent variables from your dataset as predictors, and estimate a linear regression model.
Note that not all variables in your dataset make sense as predictors. For example, categorical variables with many levels or all components of a derived variable are unreasonable.
Work through the EDA to identified good candidates. Your model should be linear in the parameters. Consider alternative functional forms for your depvar or continuous variables.
Include interpretations and comments after each chunk.
Conclude with a final assessment of the model and a reflection of the value of data visualization to aid model development.
Example Code for FIFA with depvar=Wage
DATA
options(scipen=999)library(socviz) # use for the round_df() functionload("NBA_viz.RData")library(psych)describe(NBA_viz)
Look at the shape of the relationship between the dependent variable and all of the continuous potential independent variables.
plot_scatterplot(test, by="GP")
Correlation Matrix
plot_correlation(train, type="continuous")
Variables with stronger correlation are FGM AND W which can make us infer that will be more significant in the model. The plot help us see that some variables are correlated between them like FTM and FTA and also FGM and PTS
MODEL
Linear Regression Model
From the above EDA, we chose the following variables to start. Note that, given the shape of the relationship between Wage and Age, we entered Age as a quadratic.
From experience and prior research, it is common to specify a dependent variable that is a currency variable (e.g., sales, revenue, wages) in log form.
Estimate the following model: \(log(gp) = Age + Age^2 + Potential + International.Reputation + Value.n + Special + Height + RightFoot + Skill.Moves\)\(log(gp)= PTS + FGM + FGA + FTM + FTA\)
Estimate Coefficients and show coefficients table
model <-lm(log(GP) ~ PTS + W + MIN + FGM + AGE, data=train)model %>% as_flextable
F-statistic: 116.1 on 313 and 5 DF, p-value: 0.0000
Among the factors considered, both “W” (wins) and “MIN” (minutes played) significantly influence the number of games played. For instance, for every additional win, the expected increase in games played is about 0.030, and for each extra minute played, there’s a rise in games played. On the other hand, indicators like “PTS” (points), “FGM” (field goals made), and “AGE” (age) show that were not statistically significant. The model shows a good ability to explain the variation in games played, with an adjusted R-squared value of 0.644, meaning that about 64.4% of the variability in games played can be understood by the predictors. Moreover, the low p-value of 0.0000 for the F-statistic (116.1) indicates that the overall model is highly statistically significant.
Coefficient Magnitude Plot
# note this simple code using coefplot package produces chart like Figure 6.6 in one short line. library(coefplot)coefplot(model, soft="magnitude", intercept=FALSE)
The chart helps visualize the 3 predictors that were not statistically significant because they are the ones that has the line that pass on zero, it means that the predictor variable does not have a significant impact on the outcome variable in the model.
Check for predictor independence
Using Variance Inflation Factors (VIF)
library(car)vif(model)
PTS W MIN FGM AGE
13.825569 2.570092 2.875727 13.286659 1.042462
Since a VIF value greater than 10 indicates a high degree of multicollinearity. “PTS” and “FGM” have VIF values exceeding 10, indicating potential multicollinearity issues which mean that these variables are highly correlated with other predictor variables in the model, which can affect the reliability of coefficient estimates. However, the model was made with those predictors because they were the ones that helped the model to explain the variablity more.
Residual Analysis
Residual Range
quantile(round(residuals(model, type ="deviance"),2))
We are looking for: - Random distribution of residuals vs fitted values - Normally distributed residuals : Normal Q-Q plot with values along line - Homoskedasticity with a Scale-Location line that is horizontal and no residual pattern - Minimal influential obs - that is, those outside the borders of Cook’s distance
plot(model)
Since the points closely follow the line, it can indicate that the model is performing well, and there are no apparent issues with heteroscedasticity or non-linearity. If the points in the Normal Q-Q plot deviate from the diagonal line, it suggests that the distribution of residuals does not perfectly follow a normal distribution.
Plot Fitted Value by Actual Value
library(broom)train2 <-augment(model, data=train) # this appends predicted to original dataset to build plots on your ownp <-ggplot(train2, mapping =aes(y=.fitted, x=GP))p +geom_point()
The model’s predictions can mean they are accurate, and there is a strong linear relationship between the predicted and actual values.
Plot Residuals by Fitted Values
p <-ggplot(train2, mapping =aes(y=.resid, x=.fitted))p +geom_point()
Since the points are not randomly distributed it can indicate heteroscedasticity ### Performance Evaluation
Use Model to Score test dataset (Display First 10 values - depvar and fitted values only)
library(dplyr)pred <-predict(model, newdata = test, interval="predict") #score test dataset wtih modeltest_w_predict <-cbind(test, pred) # append score to original test datasettest_w_predict %>%head(10) %>%select(PLAYER, GP, fit, lwr,upr) %>%as_flextable(show_coltype =FALSE)
# more in depth plots of performance p <-ggplot(data =subset(test_w_predict, TEAM %in%c("LAC", "GSW")),aes(x = GP,y = fit, ymin = lwr, ymax = upr,color = TEAM,fill = TEAM,group = TEAM))p +geom_point(alpha =0.5) +geom_line() +geom_ribbon(alpha =0.2, color =FALSE) +labs(title="Actual vs Fitted with Upper and Lower CI",subtitle="Golden State Warriors and Los Angeles Clippers",caption="NBA{datasetsICR}") +xlab(NULL) +ylab("fitted") +theme(legend.position ="bottom")
Since both lines in the plot (representing the teams Golden State Warriors and Los Angeles Clippers) are ascending or going up, it suggests a positive relationship between the number of games played (“GP”) and the fitted values generated by the model.
ggplot2 explore (Intro Section)
This is just replication of the first section within Chapter 6 (as FYI)
pip <-lm(log(GP) ~ W, data=train)summary(pip)
Call:
lm(formula = log(GP) ~ W, data = train)
Residuals:
Min 1Q Median 3Q Max
-2.44684 -0.30764 0.07194 0.41148 1.23955
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.399525 0.066373 36.15 <0.0000000000000002 ***
W 0.047317 0.002205 21.46 <0.0000000000000002 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.6323 on 317 degrees of freedom
Multiple R-squared: 0.5923, Adjusted R-squared: 0.591
F-statistic: 460.6 on 1 and 317 DF, p-value: < 0.00000000000000022
# knitr does not like autoplot#autoplot(lm(log(Wage.n) ~ Potential, data=train), data=train, colour = 'Nationality_Continent')library(ggplot2)p <-ggplot(data=train, mapping=aes(y=log(GP), x=W))p +geom_point(alpha =0.2) +geom_smooth(method ="lm", aes(color ="OLS", fill ="OLS"))
