library(tidyverse)
library(broom)

data <- tribble(
  ~ age, ~ gender, ~ make, ~ claim_count, ~ severity,
  "18-25", "M", "Tesla", 0, NA,
  "18-25", "M", "Tesla", 0, NA,
  "18-25", "M", "Toyota", 0, NA,
  "18-25", "M", "Toyota", 0, NA,
  "18-25", "M", "Toyota", 0, NA,
  "18-25", "M", "Toyota", 0, NA,
  "18-25", "M", "Toyota", 0, NA,
  "18-25", "M", "Toyota", 1, 2000,
  "18-25", "M", "Toyota", 2, 1500,
  "18-25", "M", "BMW", 0, NA,
  "18-25", "M", "BMW", 0, NA,
  "18-25", "M", "BMW", 1, 8000,
  "18-25", "M", "BMW", 1, 1000000,
  "18-25", "M", "BMW", 2, 15000,
  "18-25", "F", "Tesla", 0, NA,
  "18-25", "F", "Tesla", 0, NA,
  "18-25", "F", "Tesla", 0, NA,
  "18-25", "F", "Toyota", 0, NA,
  "18-25", "F", "Toyota", 0, NA,
  "18-25", "F", "Toyota", 1, 500,
  "18-25", "F", "Toyota", 2, 1000,
  "18-25", "F", "Toyota", 2, 500,
  "18-25", "F", "BMW", 1, 20000,
  "26-30", "M", "Toyota", 0, NA,
  "26-30", "M", "Toyota", 0, NA,
  "26-30", "M", "Toyota", 0, NA,
  "26-30", "M", "Toyota", 0, NA,
  "26-30", "M", "Toyota", 0, NA,
  "26-30", "M", "Toyota", 0, NA,
  "26-30", "M", "Toyota", 0, NA,
  "26-30", "M", "Toyota", 2, 500,
  "26-30", "M", "Toyota", 2, 800,
  "26-30", "M", "BMW", 0, NA,
  "26-30", "M", "BMW", 0, NA,
  "26-30", "M", "BMW", 0, NA,
  "26-30", "M", "BMW", 0, NA,
  "26-30", "M", "BMW", 0, NA,
  "26-30", "F", "Tesla", 0, NA,
  "26-30", "F", "Tesla", 0, NA,
  "26-30", "F", "Toyota", 0, NA,
  "26-30", "F", "Toyota", 0, NA,
  "26-30", "F", "Toyota", 0, NA,
  "26-30", "F", "Toyota", 0, NA,
  "26-30", "F", "BMW", 0, NA,
  "26-30", "F", "BMW", 0, NA,
  "26-30", "F", "BMW", 0, NA,
  "26-30", "F", "BMW", 0, NA,
  "26-30", "F", "BMW", 0, NA
) %>%
  mutate(id = row_number())


glm_count <- glm(
  claim_count ~ age + gender + make,
  data = data,
  family = "poisson"
)

glm_severity <- glm(
  severity ~ age + gender + make,
  family = Gamma(link = "log"),
  data = data
)

glm_severity_random_effects <- lme4::glmer(
  severity ~ age + gender + (1|make),
  family = Gamma(link = "log"),
  data = data
)

tesla_new_data <- tibble(age = "18-25", gender = "M", make = "Tesla")

predict(
  glm_count,
  tesla_new_data
)

# errors
predict(
  glm_severity,
  tesla_new_data
)

predict(
  glm_severity_random_effects,
  tesla_new_data,
  type = "response",
  allow.new.levels = TRUE
)

# this gives pretty different results than the average relativity
# approach from before. maybe the average was missing a exp() or log()?