Review of Simple Linear Regression

require("knitr")

## Loading required package: knitr

knitr::opts_chunk$set(
    comment = "#",
    fig.width=8, 
    fig.height=6, 
    cache = TRUE
)

1 Input data

issu <- data.frame (
    driving = c(5, 2, 12, 9, 15, 6, 25, 16),
    premium = c(64, 87, 50, 71, 44, 56, 42, 60)
)

y <- issu$premium
x <- issu$driving
xbar <- mean(x)
ybar <- mean(y)
n <- length(y)
plot(x,y, xlab = "Driving", ylab = "Premium")

2 Estimating regression coefficients

fit.issu <- lm(y~x)
plot(x,y)
abline (fit.issu)

3 Fitted Values, Residuals and Sum Squares

beta0 <- fit.issu$coefficients[1]
beta1 <- fit.issu$coefficients[2]
fitted1 <- beta0+beta1*x
fitted0 <- rep(ybar, n)
residual1 <- y - fitted1
residual0 <- y - fitted0

data.frame (y, fitted0, residual0, fitted1, residual1, diff.residual=fitted1-fitted0)

#    y fitted0 residual0  fitted1  residual1 diff.residual
# 1 64   59.25      4.75 68.92243  -4.922425      9.672425
# 2 87   59.25     27.75 73.56519  13.434811     14.315189
# 3 50   59.25     -9.25 58.08931  -8.089309     -1.160691
# 4 71   59.25     11.75 62.73207   8.267927      3.482073
# 5 44   59.25    -15.25 53.44654  -9.446545     -5.803455
# 6 56   59.25     -3.25 67.37484 -11.374837      8.124837
# 7 42   59.25    -17.25 37.97066   4.029335    -21.279335
# 8 60   59.25      0.75 51.89896   8.101043     -7.351043

Visualize the fitted line, fitted values, and residuals

Definitions of SSR, SSE, SST:

SST <- sum((y-fitted0)^2); SST

# [1] 1557.5

SSE <- sum ((y-fitted1)^2);SSE

# [1] 639.0065

SSR <- SST-SSE;SSR

# [1] 918.4935

SSR can be computed directly with

sum((fitted1-fitted0)^2)

# [1] 918.4935

Plot of Residual Sum Squares

num.pars <- c(1,2,n)
rss <- c(SST, SSE, 0)
par(mfrow= c(1,1))
plot(rss~num.pars, type="b", 
     xlab = "Number of Parameters",
     ylab = "Residual Sum Squares")
abline(v=1:10, h=seq(0,1600,by=100), lty=3,col="grey")

Coefficient of Determination: \(R^2\)

R2 <- SSR/SST; R2

# [1] 0.5897229

Mean Sum Squares and F statististic

f <- (SSR/1)/(SSE/(n-2)); f

# [1] 8.624264

p-value to assess whether the relationship exists (statistical significance measure)

pvf <- pf(f, df1=1, df2=n-2, lower.tail = FALSE); pvf

# [1] 0.0260588

ANOVA table

Ftable <- data.frame(df = c(1,n-2), SS=c(SSR,SSE), 
                     MSS = c(SSR/1, SSE/(n-2)), Fstat=c(f,NA),
                     p.value = c(pvf, NA),
                     R2 = c(SSR, SSE)/SST)
Ftable

#   df       SS      MSS    Fstat   p.value        R2
# 1  1 918.4935 918.4935 8.624264 0.0260588 0.5897229
# 2  6 639.0065 106.5011       NA        NA 0.4102771

# A single function to compute F table:

anova(fit.issu)

# Analysis of Variance Table
# 
# Response: y
#           Df Sum Sq Mean Sq F value  Pr(>F)  
# x          1 918.49  918.49  8.6243 0.02606 *
# Residuals  6 639.01  106.50                  
# ---
# Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

4 Understanding the Sampling Distributions of SS and F

When the slope is 0

fit.issu <- lm (y~x)
fit.issu$coefficients

# (Intercept)           x 
#   76.660365   -1.547588

fit.issu$coefficients[2] <- 0
fit.issu

# 
# Call:
# lm(formula = y ~ x)
# 
# Coefficients:
# (Intercept)            x  
#       76.66         0.00

# Simulate a dataset with modified coeffients
sim.y <- predict(fit.issu, newdata = data.frame(x))+rnorm(n,0, sigma(fit.issu))

# Fit linear Model with Simulated Data
par(mfrow=c(1,2))
fit.sim.y <- lm(sim.y~x)
plot(sim.y~x, ylim=c(0,100)) 
abline(fit.sim.y,col=2)
abline(h=mean(sim.y),col=1)

anova.fit.sim.y <- anova(fit.sim.y); anova.fit.sim.y

# Analysis of Variance Table
# 
# Response: sim.y
#           Df Sum Sq Mean Sq F value Pr(>F)
# x          1   0.54   0.542  0.0073 0.9345
# Residuals  6 443.18  73.864

ss.sim.y <- anova.fit.sim.y$`Sum Sq`
rss.sim.y <- rev(cumsum(c(0, rev(ss.sim.y))))
num.par <- cumsum(c(1,anova.fit.sim.y$Df))
plot(rss.sim.y~num.par,xlab="Number of Parameters", ylab = "Residual Sum Squares", type ="b")
abline(v=1:25, h=(0:50)*100,lty=3,col="grey")

When the slope is non-zero

fit.issu <- lm (y~x)
fit.issu$coefficients

# (Intercept)           x 
#   76.660365   -1.547588

fit.issu$coefficients[2] <- -2
fit.issu

# 
# Call:
# lm(formula = y ~ x)
# 
# Coefficients:
# (Intercept)            x  
#       76.66        -2.00

# Simulate a dataset with modified coeffients
sim.y <- predict(fit.issu, newdata = data.frame(x))+rnorm(n,0, sigma(fit.issu))

# Fit linear Model with Simulated Data
par(mfrow=c(1,2))
fit.sim.y <- lm(sim.y~x)
plot(sim.y~x, ylim=c(0,80)) 
abline(fit.sim.y,col=2)
abline(h=mean(sim.y),col=1)

anova.fit.sim.y <- anova(fit.sim.y); anova.fit.sim.y

# Analysis of Variance Table
# 
# Response: sim.y
#           Df Sum Sq Mean Sq F value  Pr(>F)  
# x          1 874.22  874.22  7.4935 0.03385 *
# Residuals  6 699.98  116.66                  
# ---
# Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

ss.sim.y <- anova.fit.sim.y$`Sum Sq`
rss.sim.y <- rev(cumsum(c(0, rev(ss.sim.y))))
num.par <- cumsum(c(1,anova.fit.sim.y$Df))
plot(rss.sim.y~num.par,xlab="Number of Parameters", ylab = "Residual Sum Squares", type ="b")
abline(v=1:25, h=(0:50)*100,lty=3,col="grey")