Statistical Modeling with R (hoopR)

Using R and the hoopR package for advanced statistical analysis of post-up efficiency:

Setup and Data Loading

# Install required packages
install.packages(c("hoopR", "tidyverse", "lme4", "broom", "ggplot2", "corrplot"))

library(hoopR)
library(tidyverse)
library(lme4)
library(broom)
library(ggplot2)
library(corrplot)

# Load NBA play-by-play data
load_nba_data <- function(seasons = 2024) {
  # Get play-by-play data
  pbp_data <- load_nba_pbp(seasons = seasons)

  # Filter for post-up plays (using shot location and play description)
  post_plays <- pbp_data %>%
    filter(
      shooting_play == TRUE,
      shot_distance <= 10,
      str_detect(type_text, "(?i)(post|hook|layup|dunk)"),
      coordinate_x != "NA"
    ) %>%
    mutate(
      post_shot = case_when(
        shot_distance <= 3 ~ "At Rim",
        shot_distance <= 10 ~ "Paint",
        TRUE ~ "Other"
      ),
      shot_made = scoring_play == TRUE
    )

  return(post_plays)
}

# Load player tracking data
get_player_stats <- function(season = 2024) {
  # Get player box scores
  player_stats <- nba_leaguedashplayerstats(
    season = paste0(season-1, "-", substr(season, 3, 4))
  )

  return(player_stats)
}

# Execute data loading
post_data <- load_nba_data(seasons = 2024)
player_stats <- get_player_stats(season = 2024)

Regression Analysis: Post-Up Efficiency Predictors

# Prepare data for regression modeling
prepare_regression_data <- function(post_plays, player_stats) {
  # Aggregate post-up statistics by player
  player_post_stats <- post_plays %>%
    group_by(athlete_id_1, athlete_id_1_name) %>%
    summarize(
      post_attempts = n(),
      post_makes = sum(shot_made),
      fg_pct = mean(shot_made),
      avg_shot_distance = mean(shot_distance, na.rm = TRUE),
      at_rim_pct = sum(post_shot == "At Rim") / n(),
      .groups = 'drop'
    ) %>%
    filter(post_attempts >= 50) %>%  # Minimum volume threshold
    rename(player_id = athlete_id_1, player_name = athlete_id_1_name)

  # Merge with player physical/performance data
  merged_data <- player_post_stats %>%
    left_join(
      player_stats %>% select(player_id, height, weight, position,
                              usage_rate, assist_rate, turnover_rate),
      by = "player_id"
    ) %>%
    na.omit()

  # Convert height to inches
  merged_data <- merged_data %>%
    mutate(
      height_inches = as.numeric(str_extract(height, "\\d+")) * 12 +
                      as.numeric(str_extract(height, "-(\\d+)"))
    )

  return(merged_data)
}

# Build regression model
model_data <- prepare_regression_data(post_data, player_stats)

# Multiple linear regression: Predicting post-up FG%
post_efficiency_model <- lm(
  fg_pct ~ height_inches + weight + at_rim_pct + avg_shot_distance +
           usage_rate + assist_rate + position,
  data = model_data
)

# Model summary
summary(post_efficiency_model)

# Extract coefficients
model_coefs <- tidy(post_efficiency_model, conf.int = TRUE)
print(model_coefs)

# Visualize coefficients
ggplot(model_coefs %>% filter(term != "(Intercept)"),
       aes(x = reorder(term, estimate), y = estimate)) +
  geom_point(size = 4, color = "steelblue") +
  geom_errorbar(aes(ymin = conf.low, ymax = conf.high),
                width = 0.2, linewidth = 1) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
  coord_flip() +
  labs(
    title = "Post-Up Efficiency Predictors",
    subtitle = "Regression Coefficients with 95% Confidence Intervals",
    x = "Predictor Variables",
    y = "Coefficient Estimate (Effect on FG%)"
  ) +
  theme_minimal(base_size = 14)

ggsave("post_efficiency_predictors.png", width = 10, height = 8, dpi = 300)

Mixed Effects Model: Player and Team Effects

# Multi-level model accounting for player and team variation
library(lme4)

# Prepare hierarchical data
hierarchical_data <- post_plays %>%
  mutate(
    shot_value = case_when(
      score_value == 3 ~ 3,
      score_value == 2 ~ 2,
      TRUE ~ 0
    ),
    defender_distance = case_when(
      qtr <= 2 ~ "First Half",
      TRUE ~ "Second Half"
    )
  ) %>%
  filter(!is.na(athlete_id_1), !is.na(team_id))

# Mixed effects logistic regression
mixed_model <- glmer(
  shot_made ~ shot_distance + at_rim_pct + (1 | athlete_id_1) + (1 | team_id),
  data = hierarchical_data %>%
    filter(post_shot %in% c("At Rim", "Paint")),
  family = binomial(link = "logit"),
  control = glmerControl(optimizer = "bobyqa")
)

# Model summary
summary(mixed_model)

# Extract random effects (player-specific deviations)
player_effects <- ranef(mixed_model)$athlete_id_1 %>%
  rownames_to_column("player_id") %>%
  rename(player_effect = `(Intercept)`) %>%
  arrange(desc(player_effect))

print("Top 15 Players by Random Effect (Post-Up Skill Above Average):")
print(head(player_effects, 15))

# Visualize random effects distribution
ggplot(player_effects, aes(x = player_effect)) +
  geom_histogram(bins = 30, fill = "darkgreen", alpha = 0.7, color = "black") +
  geom_vline(xintercept = 0, linetype = "dashed", color = "red", linewidth = 1) +
  labs(
    title = "Distribution of Player-Specific Post-Up Skill",
    subtitle = "Random Effects from Mixed Model",
    x = "Player Effect (Log-Odds Scale)",
    y = "Frequency"
  ) +
  theme_minimal(base_size = 14)

ggsave("player_effects_distribution.png", width = 10, height = 6, dpi = 300)

Shot Location Heatmap Analysis

# Create shot location heatmap for post-up attempts
library(ggplot2)
library(hexbin)

# Filter for post area shots
post_shots_spatial <- post_plays %>%
  filter(
    !is.na(coordinate_x),
    !is.na(coordinate_y),
    shot_distance <= 12
  ) %>%
  mutate(
    shot_outcome = ifelse(shot_made, "Made", "Missed")
  )

# Create heatmap
ggplot(post_shots_spatial, aes(x = coordinate_x, y = coordinate_y)) +
  stat_density_2d(aes(fill = ..level..), geom = "polygon", alpha = 0.6) +
  scale_fill_gradient(low = "yellow", high = "red", name = "Shot Density") +
  facet_wrap(~shot_outcome) +
  labs(
    title = "Post-Up Shot Location Heatmap",
    subtitle = "Spatial Distribution of Made vs Missed Shots",
    x = "Court X Coordinate",
    y = "Court Y Coordinate"
  ) +
  theme_minimal(base_size = 14) +
  theme(legend.position = "bottom")

ggsave("post_shot_heatmap.png", width = 12, height = 6, dpi = 300)

# Calculate efficiency by court zone
zone_efficiency <- post_shots_spatial %>%
  mutate(
    zone = case_when(
      abs(coordinate_x) < 4 & coordinate_y < 5 ~ "Restricted Area",
      abs(coordinate_x) < 8 & coordinate_y < 10 ~ "Low Post",
      coordinate_y < 10 ~ "Mid Post",
      TRUE ~ "High Post"
    )
  ) %>%
  group_by(zone) %>%
  summarize(
    attempts = n(),
    fg_pct = mean(shot_made),
    avg_distance = mean(shot_distance),
    .groups = 'drop'
  ) %>%
  arrange(desc(fg_pct))

print("Post-Up Efficiency by Court Zone:")
print(zone_efficiency)

# Visualize zone efficiency
ggplot(zone_efficiency, aes(x = reorder(zone, fg_pct), y = fg_pct, fill = zone)) +
  geom_col(color = "black", linewidth = 0.5) +
  geom_text(aes(label = sprintf("%.1f%%", fg_pct * 100)),
            vjust = -0.5, size = 5, fontface = "bold") +
  scale_y_continuous(labels = scales::percent, limits = c(0, 0.7)) +
  labs(
    title = "Post-Up Shooting Efficiency by Court Zone",
    x = "Court Zone",
    y = "Field Goal Percentage"
  ) +
  theme_minimal(base_size = 14) +
  theme(legend.position = "none")

ggsave("zone_efficiency.png", width = 10, height = 6, dpi = 300)

Time Series Analysis: Within-Game Trends

# Analyze post-up efficiency across game situations
game_situation_analysis <- post_plays %>%
  mutate(
    quarter = qtr,
    time_remaining = period_seconds_remaining,
    score_margin = score_value - lag(score_value, default = 0)
  ) %>%
  group_by(quarter) %>%
  summarize(
    attempts = n(),
    fg_pct = mean(shot_made),
    avg_shot_distance = mean(shot_distance, na.rm = TRUE),
    .groups = 'drop'
  )

print("Post-Up Efficiency by Quarter:")
print(game_situation_analysis)

# Visualize quarter-by-quarter trends
ggplot(game_situation_analysis, aes(x = factor(quarter), y = fg_pct, group = 1)) +
  geom_line(linewidth = 1.5, color = "darkblue") +
  geom_point(size = 4, color = "darkred") +
  geom_text(aes(label = sprintf("%.1f%%", fg_pct * 100)),
            vjust = -1, size = 4, fontface = "bold") +
  scale_y_continuous(labels = scales::percent, limits = c(0.3, 0.6)) +
  labs(
    title = "Post-Up Efficiency Across Game Quarters",
    subtitle = "How does post-up effectiveness change throughout the game?",
    x = "Quarter",
    y = "Field Goal Percentage"
  ) +
  theme_minimal(base_size = 14)

ggsave("quarter_efficiency.png", width = 10, height = 6, dpi = 300)

Correlation Analysis

# Examine correlations between post-up metrics and team success
library(corrplot)

# Aggregate team-level post statistics
team_post_stats <- post_plays %>%
  group_by(team_id) %>%
  summarize(
    post_freq = n(),
    post_fg_pct = mean(shot_made),
    at_rim_rate = sum(post_shot == "At Rim") / n(),
    avg_shot_distance = mean(shot_distance, na.rm = TRUE),
    .groups = 'drop'
  )

# Merge with team performance metrics
team_performance <- player_stats %>%
  group_by(team_id) %>%
  summarize(
    team_off_rating = mean(off_rating, na.rm = TRUE),
    team_def_rating = mean(def_rating, na.rm = TRUE),
    team_net_rating = team_off_rating - team_def_rating,
    .groups = 'drop'
  )

team_analysis <- team_post_stats %>%
  left_join(team_performance, by = "team_id") %>%
  na.omit()

# Correlation matrix
cor_matrix <- cor(team_analysis %>% select(-team_id))

# Visualize correlations
corrplot(cor_matrix, method = "color", type = "upper",
         addCoef.col = "black", number.cex = 0.8,
         tl.col = "black", tl.srt = 45,
         title = "Correlation: Post-Up Metrics vs Team Performance",
         mar = c(0, 0, 2, 0))