Whenever working with new data it’s good practice to run descriptive statistics and visualizations
# run in your R to get clearer output
describe(GPA) vars n mean sd median trimmed mad min max range skew kurtosis
year1gpa 1 20 2.62 0.64 2.72 2.65 0.49 1.1 3.86 2.76 -0.41 0.04
act 2 20 24.25 2.59 24.00 24.25 2.97 20.0 30.00 10.00 0.15 -0.66
se
year1gpa 0.14
act 0.58The two variables seem fairly normally distributed.
GPA %>%
# When creating multiple plots at once, long data works best
pivot_longer(cols = colnames(GPA),
names_to = "variable",
values_to = "value") %>%
ggplot(aes(x = value)) +
geom_density() +
# the facet_wrap() fucntion generates plots for each value in "variable"
facet_wrap(~variable, scales = "free") +
theme(axis.title.y=element_blank(),
axis.text.y=element_blank(),
axis.ticks.y=element_blank(),
axis.line.y = element_blank())
A scatterplot is ideal for visualizing bivariate relations
We have seen this code before
ggplot(GPA,
aes(x = act, y = year1gpa)) +
geom_point()
There seems to be a positive trend (as ACT scores increase, first year GPA also increases). Let’s draw a regression line to confirm that.
A scatterplot is ideal for visualizing bivariate relations
ggplot(GPA,
aes(x = act, y = year1gpa)) +
geom_point() +
geom_smooth(method = "lm",
se = FALSE)
The regression line confirms our observation.
A scatterplot is ideal for visualizing bivariate relations
ggplot(GPA,
aes(x = act, y = year1gpa)) +
geom_point() +
geom_smooth(method = "lm",
se = FALSE) +
geom_smooth(method = "loess",
color = "red",
se = FALSE)
For many reasons, once we have a regression model, we may want to predict values of \(Y\) based on some values of \(X\).
ggplot(GPA,
aes(x = act, y = year1gpa)) +
geom_point(alpha = .5, shape = 1) +
geom_smooth(method = "lm",
se = FALSE) +
geom_segment(x = 23,
xend = 23,
y = 0,
yend = coef(reg)[1] + coef(reg)[2]*23,
linetype = 2) +
geom_segment(x = 0,
xend = 23,
y = coef(reg)[1] + coef(reg)[2]*23,
yend = coef(reg)[1] + coef(reg)[2]*23,
linetype = 2) +
geom_segment(x = 27,
xend = 27,
y = 0,
yend = coef(reg)[1] + coef(reg)[2]*27,
linetype = 2) +
geom_segment(x = 0,
xend = 27,
y = coef(reg)[1] + coef(reg)[2]*27,
yend = coef(reg)[1] + coef(reg)[2]*27,
linetype = 2) 
ggplot(GPA,
aes(x = act, y = year1gpa)) +
geom_point(shape = 1) 
ggplot(GPA,
aes(x = act, y = year1gpa)) +
geom_point(shape = 1) +
geom_point(aes(y = predictions),
colour = "red")
ggplot(GPA,
aes(x = act, y = year1gpa)) +
geom_point(shape = 1) +
geom_point(aes(y = predictions),
color = "red") +
geom_segment(y = predictions,
x = GPA$act,
xend = GPA$act,
yend = GPA$year1gpa,
linetype = 2,
alpha = .5)
The fact that residuals must be normally distributed is an assumption built into regression (why? The appendix explains it). QQplots can help us infer how much residuals deviate from normality.
The qqPlot() function from the car package is helpful here:
car::qqPlot(reg, id = FALSE)As long as roughly 95% of the residuals (dots) are within the band, there should be no cause for concern.
Another way of checking that residuals are evenly distributed around 0, and do not show any type of relationship with the predictor.
ggplot(GPA,
aes(x = act, y = residuals)) +
geom_point() +
geom_hline(yintercept = 0) +
geom_smooth(method = "loess",
se = FALSE)Once again, the loess line does not do too well because there are few points. The points look reasonably distributed around the 0 line to me.
On the previous slide you saw we had \(Y_i\), \(\hat{Y_i}\), \(\hat{Y}\), and crucially, just \(Y\) (). One important thing to keep in mind is that, whenever we run statistics, the question that we are asking is really…
Likewise, regression tries to answer how your observed sample of the \(Y\) random variable comes into existence.
The \(Y \sim N(b_0 + b_1X, \epsilon)\) (for simplicity, I dropped the \(\mu\) and \(\sigma\)) is how linear regression thinks \(Y\) crosses over the dotted line and ends up in our data.
We can visualize what our \(Y\) variable should look like according to our regression model in a perfect world.
ggplot(GPA,
aes(x = act, y = year1gpa)) +
# assume we have no data yet
geom_point(alpha = 0) +
geom_smooth(method = "lm", se = FALSE) +
ylim(0, 5)The line is what the mean of year1gpa should be at every value of act. However, there is no data yet.
We can visualize what our \(Y\) variable should look like according to our regression model in a perfect world.
# Generate idealized data
set.seed(7757)
# generate X values
act <- rep(seq(min(GPA$act), max(GPA$act), by =.1), 150)
# the code below generates data according to the regression model
intercept <- coef(reg)[1]
slope <- coef(reg)[2]
residual_var <- sigma(reg)
ideal_data <- rnorm(length(act), mean = intercept + act*slope, sd = sigma(reg))
# plot the idealized data
ggplot(mapping = aes(x = act, y = ideal_data)) +
geom_point(alpha = .2)+
geom_smooth(method = "lm", se = FALSE) +
ylim(0, 5)Each “column of dots” is normally distributed around the line, with standard deviation of \(.61\), the residual variance (\(\epsilon\)). The dots are essentially our residuals, \(e\).
The \(Y\) variable looking like this is the only “assumption” of regression.
The \(Y\) variable never looks like this.
This is what the \(Y\) variable looks like usually (our data).
Whether you think your observed data it “too different” from the expected data is usually a very subjective process.