Wrappers to generate the data, plots and animations for the post on polynomial regression.

polynomials.R

#### 2020-10-17-visualising-polynomial-regression -----

#### set up -------------------------------------------

# load packages
library(tidyverse)
library(gganimate)
library(data.table)
library(magick)
library(here)

my_theme <- theme(
  text = element_text(size = 12),
  axis.title = element_text(face = "bold", colour = "#e8b941"),
  axis.text = element_text(colour = "white", size = 12),
  panel.grid = element_line(colour = "grey"),
  panel.grid.major.x = element_blank(),
  panel.grid.minor.x = element_blank(),
  panel.grid.minor.y = element_blank(),
  strip.text = element_text(size = 15),
  legend.position = "top",
  legend.box = "vertical",
  legend.background = element_rect(fill = "#141414"),
  legend.title = element_text(colour = "#e8b941", face = "bold"),
  legend.key = element_rect(fill = "#3b3d42"),
  legend.text = element_text(colour = "white", size = 12),
  panel.background = element_rect(fill = "#3b3d42"),
  plot.background = element_rect(fill = "#141414", colour = "#141414")
)

# create functions
sim_poly <- function(x = seq(-1, 1, 0.1), degree = 1, coefs) {
  if (length(coefs)<degree) stop("You must provide as many coefficients as degrees, not including the intercept")
  cd <- coefs
  cd$x <- x
  c <- expand.grid(coefs)
  cd <- expand.grid(cd)
  p <- as.matrix(poly(x, degree = degree))
  d <- cbind(x, p)
  e <- full_join(as.data.frame(d), as.data.frame(cd), by = "x")
  e$y <- e[,2+degree] + rowSums(as.matrix(e[,2:(degree+1)] * e[,(degree+3):ncol(e)]))
  return(e)
}

gif_append <- function(a, b, duration = 6, rewind = TRUE){
  a_gif <- animate(a, duration = duration, rewind = rewind)
  b_gif <- animate(b, duration = duration, rewind = rewind)
  
  a_mgif <- image_read(a_gif)
  b_mgif <- image_read(b_gif)
  
  new_gif <- image_append(c(a_mgif[1], b_mgif[1]))
  for(i in 2:length(a_mgif)){
    combined <- image_append(c(a_mgif[i], b_mgif[i]))
    new_gif <- c(new_gif, combined)
  }
  return(new_gif)
}


#### example students ----------------------------------------------------------
attention <- sim_poly(0:59, 3, coefs = list(b0 = 0.5, b1 = 0.4, b2 = 2, b3 = 0.1)) %>% 
  mutate(y_rand = rnorm(n = nrow(.), mean = y, sd = 0.1))

ggplot(attention, aes(x = x, y = y_rand)) + 
  geom_point(size = 3, shape = 1, stroke = 1, colour = "orange") +
  #geom_line(aes(y = y), size = 1, colour = "orange") +
  labs(x = "Time (minutes)",
       y = "Proportion of student paying attention") +
  scale_y_continuous(limits = c(0, 1)) +
  my_theme
  

#### create animations ---------------------------------------------------------

# varying intercept
dat <- list(linear = sim_poly(degree = 1, coefs = list(b0 = seq(-1, 1, 0.1), b1 = 1)),
            quadratic = sim_poly(degree = 2, coefs = list(b0 = seq(-1, 1, 0.1), b1 = 1, b2 = 1)),
            cubic = sim_poly(degree = 3, coefs = list(b0 = seq(-1, 1, 0.1), b1 = 1, b2 = 1, b3 = 1))) %>% 
  bind_rows(.id = "degree") %>% 
  mutate(degree = str_to_sentence(degree)) %>% 
  mutate_all(function(x) ifelse(is.na(x), 0, x))

anim <- ggplot(dat, aes(x, y, colour = degree)) +
  geom_line(size = 1.5) +
  labs(x = "Input (x)", y = expression(y==beta[0]+beta[1]*x+beta[2]*x^2+beta[3]*x^3), colour = "Degree") +
  scale_colour_manual(values = c("#00FFFF", "#FFFF00", "#FF3030")) +
  theme_bw() +
  my_theme +
  theme(legend.position = c(0.05, 0.95),
        legend.justification = c(0, 1),
        axis.text = element_text(size = 20),
        text = element_text(size = 25),
        legend.direction = "horizontal") +
  transition_time(b0)

anim_bar <- dat %>% 
  select(degree, b0, b1, b2, b3) %>% 
  mutate(b0_time = b0) %>% 
  pivot_longer(-c(degree, b0_time), names_to = "coef", values_to = "values") %>% 
  filter(degree %in% "Cubic") %>% 
  ggplot(aes(coef, values, fill = coef, colour = coef)) +
  geom_col(position = position_identity(), width = 0.5) +
  geom_hline(yintercept = 0) +
  labs(x = "Coefficient", y = "Value", fill = "Coefficient") +
  scale_fill_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_colour_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_x_discrete(labels = c(expression(beta[0]), expression(beta[1]), expression(beta[2]), expression(beta[3]))) +
  theme_bw() +
  my_theme +
  theme(legend.position = "none",
        axis.text = element_text(size = 20),
        text = element_text(size = 25)) +
  transition_time(b0_time)

gif_append(anim, anim_bar) 


# linear
dat <- sim_poly(degree = 1, coefs = list(b0 = 0, b1 = seq(-1, 1, 0.05)))

anim <- ggplot(dat, aes(x, y, group = b1)) +
  geom_line(colour = "#00FFFF", size = 1.5) +
  labs(x = "Input (x)", y = expression(y==beta[0]+beta[1]*x+beta[2]*x^2+beta[3]*x^3)) +
  theme_bw() +
  my_theme +
  theme(legend.position = "none",
        axis.text = element_text(size = 20),
        text = element_text(size = 25)) +
  transition_time(b1)

anim_bar <- dat %>% 
  mutate(b2 = 0, b3 = 0) %>% 
  select(b0, b1, b2, b3) %>% 
  mutate(b1_time = b1) %>% 
  pivot_longer(-b1_time, names_to = "coef", values_to = "values") %>% 
  ggplot(aes(coef, values, fill = coef, colour = coef)) +
  geom_col(position = position_identity(), width = 0.5) +
  geom_hline(yintercept = 0) +
  labs(x = "Coefficient", y = "Value", fill = "Coefficient") +
  scale_fill_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_colour_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_x_discrete(labels = c(expression(beta[0]), expression(beta[1]), expression(beta[2]), expression(beta[3]))) +
  theme_bw() +
  my_theme +
  theme(legend.position = "none",
        axis.text = element_text(size = 20),
        text = element_text(size = 25)) +
  transition_time(b1_time)

gif_append(anim, anim_bar)

# quadratic
dat <- sim_poly(degree = 2, coefs = list(b0 = 0, b1 = seq(-1, 1, 0.01), b2 = seq(-1, 1, 0.05)))
anim <- ggplot(dat, aes(x, y, colour = b1, group = b1)) +
  geom_line(size = 1) +
  labs(x = "Input (x)", y = expression(y==beta[0]+beta[1]*x+beta[2]*x^2+beta[3]*x^3), colour = expression(beta[1])) +
  scale_colour_gradient2(low = "#00FFFF", mid = "#FFFF00", high = "#FF3030") +
  theme_bw() +
  my_theme +
  theme(legend.position = c(0.05, 0.95),
        legend.justification = c(0, 1),
        axis.text = element_text(size = 20),
        text = element_text(size = 25),
        legend.direction = "horizontal") +
  transition_time(b2)

anim_bar <- dat %>% 
  mutate(b3 = 0) %>% 
  select(b0, b1, b2, b3) %>% 
  mutate(b2_time = b2) %>% 
  pivot_longer(-b2_time, names_to = "coef", values_to = "values") %>% 
  ggplot(aes(coef, values, fill = coef, colour = coef)) +
  geom_col(position = position_identity(), width = 0.5) +
  geom_hline(yintercept = 0) +
  labs(x = "Coefficient", y = "Value", fill = "Coefficient") +
  scale_fill_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_colour_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_x_discrete(labels = c(expression(beta[0]), expression(beta[1]), expression(beta[2]), expression(beta[3]))) +
  theme_bw() +
  my_theme +
  theme(legend.position = "none",
        axis.text = element_text(size = 20),
        text = element_text(size = 25)) +
  transition_time(b2_time)

gif_append(anim, anim_bar)


# cubic
dat <- sim_poly(degree = 3, coefs = list(b0 = 0, b1 = seq(-1, 1, 0.05), b2 = seq(-1, 1, 0.05), b3 = seq(-1, 1, 0.5)))
anim <- dat %>% 
  mutate(b3 = paste0("beta[3]==", b3)) %>% 
  ggplot(aes(x, y, colour = b1, group = b1)) +
  facet_wrap(~b3, labeller = label_parsed) +
  geom_line(size = 1) +
  labs(x = "Input (x)", y = expression(y==beta[0]+beta[1]*x+beta[2]*x^2+beta[3]*x^3), colour = expression(beta[1])) +
  scale_colour_gradient2(low = "#00FFFF", mid = "#FFFF00", high = "#FF3030") +
  theme_bw() +
  my_theme +
  theme(legend.position = c(0.85, 0.25)) +
  transition_time(b2)

anim_bar <- dat %>% 
  select(b0, b1, b2, b3) %>% 
  mutate(b2_time = b2) %>% 
  pivot_longer(-b2_time, names_to = "coef", values_to = "values") %>% 
  ggplot(aes(coef, values, fill = coef, colour = coef)) +
  geom_col(position = position_identity(), width = 0.5) +
  geom_hline(yintercept = 0) +
  labs(x = "Coefficient", y = "Value", fill = "Coefficient") +
  scale_fill_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_colour_manual(values = c("white", "#00FFFF", "#FFFF00", "#FF3030")) +
  scale_x_discrete(labels = c(expression(beta[0]), expression(beta[1]), expression(beta[2]), expression(beta[3]))) +
  theme_bw() +
  my_theme +
  theme(legend.position = "none",
        axis.text = element_text(size = 20),
        text = element_text(size = 25)) +
  transition_time(b2_time)

gif_append(anim, anim_bar)