1-grad.R

# Sun Nov 14 19:58:50 2021 ------------------------------

# Regression Oranges
library(RColorBrewer)

data("Orange")

# Determine number of iterations
niter <- 25
# Determine learning rate / step size
alpha <- 1e-4
set.seed(1)
b0 <- rnorm(1) # intercept
b1 <- rnorm(1) # slope

# Set palette
cols <- colorRampPalette(rev(brewer.pal(n = 7, name = "RdBu")))(niter)

# Plot
layout(matrix(c(1,1,2,2,3), nrow = 1, ncol = 5))
plot(age ~ circumference, data = Orange, pch = 16,
     xlab = "Circumference (mm)",
     ylab = "Age (days)",
     col = rgb(0,0,1,.5))

# Perform gradient descent
slopes <- rep(NA, niter)
intercepts <- rep(NA, niter)
for(i in 1:niter){
  # prediction
  y <- b0 + b1 * Orange$circumference
  # b0 = b0 - dJ/da * alpha
  b0 <- b0 - sum(y - Orange$age) / nrow(Orange) * alpha
  # b1 = b1 - dJ/db * alpha
  b1 <- b1 -  sum((y - Orange$age) * Orange$circumference) / nrow(Orange) * alpha
  abline(a = b0, b = b1, col = cols[i], lty = 2)
  # Save estimates over all iterations
  intercepts[i] <- b0
  slopes[i] <- b1
}
title("Regression Fit")

# Cost function contour
allCombs <- expand.grid(b0 = seq(-50, 50, length.out = 100),
                        b1 = seq(7.5, 8.5, length.out = 100))
res <- matrix(NA, 100, 100)
# a by rows, b by cols
for(i in 1:nrow(allCombs)){
  y <- allCombs$b0[i] + allCombs$b1[i] * Orange$circumference
  res[i] <- sum((y - Orange$age)**2)/(2*nrow(Orange))
}

# Plot contour
contour(t(res), xlab = expression(beta[1]),
        ylab = expression(beta[0]), axes = F, nlevels = 25)
axis(1, at = c(0, .5, 1), labels = c(7.5, 8, 8.5))
axis(2, at = c(0, .5, 1), labels = c(-50, 0, 50))
points((slopes - min(allCombs$b1)) / diff(range(allCombs$b1)),
       (intercepts - min(allCombs$b0)) / diff(range(allCombs$b0)),
       pch = 4, col = cols)
title("MSE contour")

# Add colorbar
z = matrix(1:niter, nrow = 1)
image(1, 1:niter, z,
      col = cols, axes = FALSE, xlab = "", ylab = "")
title("Iteration No.")
axis(2, at = c(1, seq(5, niter, by=5)))