Introduction

Imagine you’re the performance analyst for a national Olympic weightlifting federation, and the head coach walks in asking a hard question: which athletes are actually worth the investment? That question is the whole reason for this project. It falls under Objective 4, recruitment and talent identification, and the goal is simple to state, harder to answer — figure out which athlete profiles tend to medal.

For data, I’m using the “120 Years of Olympic History: Athletes and Results” dataset. It originally came from Kaggle and sports-reference.com; my copy came through the TidyTuesday mirror. I filtered it down to Weightlifting only. Records run from 1896 through 2016, and cover the basics — age, height, weight, sex — along with country (NOC), event/weight class, and whether the athlete medaled.

My focus on this is because a federation with limited resources can’t scout everyone equally. If certain traits, or certain program-level factors, actually correlate with medaling, coaches can put scouting dollars where they’ll count before the next Olympic cycle starts.

Data Exploration and Preprocessing

Load and Filter Data

raw <- read_csv("data/athlete_events_raw.csv", show_col_types = FALSE) %>%
  clean_names()

wl <- raw %>% filter(sport == "Weightlifting")

cat("Total weightlifting records:", nrow(wl), "\n")
## Total weightlifting records: 3937

Handle Missing Values

Early Olympic Games didn’t track everything I’d want. Plenty of older records are missing age, height, or weight entirely. Since those three fields are core predictors here, I just drop any row missing one of them — a smaller clean dataset beats a bigger messy one.

wl <- wl %>% filter(!is.na(age), !is.na(height), !is.na(weight))
cat("Records after removing missing age/height/weight:", nrow(wl), "\n")
## Records after removing missing age/height/weight: 2974

Feature Engineering

A few new columns make the modeling step possible. BMI (weight divided by height in meters, squared) gives a body-composition measure that’s comparable across weight classes — something raw weight alone can’t do. I also build the actual target variable, medal_binary: did the athlete medal at all (gold, silver, or bronze), or not.

Two more additions round things out. An era flag splits records before and after 2000, since rules and equipment shifted enough over the decades to be worth capturing. And weight_class gets pulled straight out of the event name, so something like “Weightlifting Men’s Featherweight” just becomes “Featherweight”.

wl <- wl %>%
  mutate(
    bmi = weight / ((height / 100)^2),
    medal_binary = if_else(!is.na(medal) & medal %in% c("Gold", "Silver", "Bronze"), 1L, 0L),
    medal_binary = factor(medal_binary, levels = c(0, 1), labels = c("No_Medal", "Medal")),
    era = factor(if_else(year >= 2000, "Post_2000", "Pre_2000"), levels = c("Pre_2000", "Post_2000")),
    weight_class = str_remove(event, "^Weightlifting (Men's|Women's) "),
    sex = factor(sex)
  ) %>%
  select(id, name, sex, age, height, weight, bmi, noc, team, year, era,
         weight_class, event, medal, medal_binary)

wl_clean <- wl

Train / Test Split

Here I split the data 80/20 into training and test sets, stratified by medal outcome. That keeps the medal rate consistent across both sets, so I’m not accidentally testing on a skewed sample.

set.seed(42)
train_idx <- createDataPartition(wl_clean$medal_binary, p = 0.8, list = FALSE)
train_data <- wl_clean[train_idx, ]
test_data  <- wl_clean[-train_idx, ]

cat("Train rows:", nrow(train_data), " | Test rows:", nrow(test_data), "\n")
## Train rows: 2380  | Test rows: 594
cat("Train medal rate:", round(mean(train_data$medal_binary == "Medal"), 3),
    " | Test medal rate:", round(mean(test_data$medal_binary == "Medal"), 3), "\n")
## Train medal rate: 0.179  | Test medal rate: 0.178

Exploratory Summary Statistics

Before modeling anything, it’s worth just looking at the plain numbers first.

summary(wl_clean %>% select(age, height, weight, bmi))
##       age            height          weight            bmi       
##  Min.   :15.00   Min.   :140.0   Min.   : 47.00   Min.   :14.68  
##  1st Qu.:22.00   1st Qu.:160.0   1st Qu.: 60.00   1st Qu.:23.88  
##  Median :25.00   Median :168.0   Median : 75.00   Median :26.53  
##  Mean   :25.29   Mean   :167.8   Mean   : 79.61   Mean   :27.75  
##  3rd Qu.:28.00   3rd Qu.:175.0   3rd Qu.: 90.00   3rd Qu.:30.35  
##  Max.   :43.00   Max.   :205.0   Max.   :176.50   Max.   :53.79
wl_clean %>%
  group_by(sex) %>%
  summarise(n = n(), medal_rate = round(mean(medal_binary == "Medal"), 3))
## # A tibble: 2 × 3
##   sex       n medal_rate
##   <fct> <int>      <dbl>
## 1 F       460      0.228
## 2 M      2514      0.17
wl_clean %>%
  filter(medal_binary == "Medal") %>%
  count(noc, sort = TRUE) %>%
  head(10)
## # A tibble: 10 × 2
##    noc       n
##    <chr> <int>
##  1 CHN      55
##  2 URS      51
##  3 BUL      36
##  4 POL      32
##  5 USA      30
##  6 RUS      26
##  7 HUN      20
##  8 KAZ      16
##  9 IRI      14
## 10 JPN      14

Exploratory Visualizations

A few charts make these patterns easier to spot than raw numbers ever could.

ggplot(wl_clean, aes(x = age, fill = medal_binary)) +
  geom_density(alpha = 0.5) +
  labs(title = "Age Distribution by Medal Outcome", x = "Age", y = "Density", fill = "Outcome") +
  theme_minimal()

ggplot(wl_clean, aes(x = medal_binary, y = bmi, fill = medal_binary)) +
  geom_boxplot(alpha = 0.7) +
  labs(title = "BMI by Medal Outcome", x = NULL, y = "BMI", fill = "Outcome") +
  theme_minimal()

wl_clean %>%
  group_by(year) %>%
  summarise(medal_rate = mean(medal_binary == "Medal")) %>%
  ggplot(aes(x = year, y = medal_rate)) +
  geom_line(color = "steelblue") +
  geom_point(color = "steelblue") +
  labs(title = "Medal Rate by Olympic Year", x = "Year", y = "Medal Rate") +
  theme_minimal()

ggplot(wl_clean, aes(x = height, y = weight, color = medal_binary)) +
  geom_point(alpha = 0.5) +
  labs(title = "Height vs Weight by Medal Outcome", x = "Height (cm)", y = "Weight (kg)", color = "Outcome") +
  theme_minimal()

top_classes <- wl_clean %>% count(weight_class, sort = TRUE) %>% slice_head(n = 15) %>% pull(weight_class)
wl_clean %>%
  filter(weight_class %in% top_classes) %>%
  count(weight_class, medal_binary) %>%
  ggplot(aes(x = reorder(weight_class, n), y = n, fill = medal_binary)) +
  geom_col() +
  coord_flip() +
  labs(title = "Records by Weight Class and Medal Outcome", x = "Weight Class", y = "Count", fill = "Outcome") +
  theme_minimal()

Analysis and Results

Model Feature Preparation

One problem before I can fit anything: noc, the country code, has 143 distinct levels. That’s far too many categories to throw into a factor variable directly — the models would either choke on it or overfit badly. So instead, each country gets collapsed into a single number: its historical weightlifting medal count, a rough stand-in for “program strength.” It’s not a perfect solution. This is a backward-looking, aggregate feature, not something you’d actually know ahead of time for a brand-new athlete with no international track record — worth flagging now, and I’ll come back to it in the discussion.

full_data <- bind_rows(train_data, test_data)
country_strength <- full_data %>%
  group_by(noc) %>%
  summarise(country_medal_count = sum(medal_binary == "Medal"), .groups = "drop")

prep_features <- function(df) {
  df %>%
    left_join(country_strength, by = "noc") %>%
    mutate(
      sex = factor(sex),
      era = factor(era, levels = c("Pre_2000", "Post_2000")),
      weight_class = factor(weight_class),
      medal_binary = factor(medal_binary, levels = c("No_Medal", "Medal"))
    ) %>%
    select(medal_binary, sex, age, height, weight, bmi, era, weight_class, country_medal_count)
}

train_m <- prep_features(train_data)
test_m  <- prep_features(test_data)

all_levels <- union(levels(train_m$weight_class), levels(test_m$weight_class))
train_m$weight_class <- factor(train_m$weight_class, levels = all_levels)
test_m$weight_class  <- factor(test_m$weight_class, levels = all_levels)

Model 1: Logistic Regression

Let’s start simple, before trying anything fancier.

log_model <- glm(medal_binary ~ sex + age + height + weight + bmi + era +
                    weight_class + country_medal_count,
                  data = train_m, family = binomial)
summary(log_model)
## 
## Call:
## glm(formula = medal_binary ~ sex + age + height + weight + bmi + 
##     era + weight_class + country_medal_count, family = binomial, 
##     data = train_m)
## 
## Coefficients:
##                                    Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       29.470274 882.786132   0.033   0.9734    
## sexM                              -0.325407   0.269127  -1.209   0.2266    
## age                               -0.012250   0.015158  -0.808   0.4190    
## height                            -0.091966   0.049098  -1.873   0.0611 .  
## weight                             0.071614   0.043347   1.652   0.0985 .  
## bmi                               -0.188744   0.141176  -1.337   0.1812    
## eraPost_2000                      -0.044962   0.154364  -0.291   0.7708    
## weight_classBantamweight         -16.980872 882.743458  -0.019   0.9847    
## weight_classFeatherweight        -16.569562 882.743419  -0.019   0.9850    
## weight_classFlyweight            -16.886131 882.743483  -0.019   0.9847    
## weight_classHeavyweight          -16.113318 882.743429  -0.018   0.9854    
## weight_classHeavyweight I        -16.287872 882.743469  -0.018   0.9853    
## weight_classHeavyweight II       -16.545426 882.743496  -0.019   0.9850    
## weight_classLight-Heavyweight    -16.277454 882.743395  -0.018   0.9853    
## weight_classLightweight          -16.732942 882.743400  -0.019   0.9849    
## weight_classMiddle-Heavyweight   -16.386105 882.743408  -0.019   0.9852    
## weight_classMiddleweight         -16.378919 882.743392  -0.019   0.9852    
## weight_classSuper-Heavyweight    -16.377964 882.743520  -0.019   0.9852    
## weight_classUnlimited, One Hand    0.379180 984.471336   0.000   0.9997    
## weight_classUnlimited, Two Hands -13.439263 882.743790  -0.015   0.9879    
## country_medal_count                0.068376   0.004015  17.029   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 2236.5  on 2379  degrees of freedom
## Residual deviance: 1840.0  on 2359  degrees of freedom
## AIC: 1882
## 
## Number of Fisher Scoring iterations: 13
log_probs <- predict(log_model, newdata = test_m, type = "response")
log_pred  <- factor(if_else(log_probs > 0.5, "Medal", "No_Medal"), levels = c("No_Medal", "Medal"))

log_cm <- confusionMatrix(log_pred, test_m$medal_binary, positive = "Medal")
log_cm
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction No_Medal Medal
##   No_Medal      471    80
##   Medal          17    26
##                                           
##                Accuracy : 0.8367          
##                  95% CI : (0.8045, 0.8655)
##     No Information Rate : 0.8215          
##     P-Value [Acc > NIR] : 0.1817          
##                                           
##                   Kappa : 0.2742          
##                                           
##  Mcnemar's Test P-Value : 3.071e-10       
##                                           
##             Sensitivity : 0.24528         
##             Specificity : 0.96516         
##          Pos Pred Value : 0.60465         
##          Neg Pred Value : 0.85481         
##              Prevalence : 0.17845         
##          Detection Rate : 0.04377         
##    Detection Prevalence : 0.07239         
##       Balanced Accuracy : 0.60522         
##                                           
##        'Positive' Class : Medal           
## 
log_roc <- roc(response = test_m$medal_binary, predictor = log_probs,
                levels = c("No_Medal", "Medal"), direction = "<")
cat("Logistic Regression AUC:", round(auc(log_roc), 3), "\n")
## Logistic Regression AUC: 0.806

Model 2: Random Forest

Now for something with a bit more flexibility.

set.seed(42)
rf_model <- randomForest(medal_binary ~ sex + age + height + weight + bmi + era +
                            weight_class + country_medal_count,
                          data = train_m, ntree = 500, importance = TRUE)
rf_model
## 
## Call:
##  randomForest(formula = medal_binary ~ sex + age + height + weight +      bmi + era + weight_class + country_medal_count, data = train_m,      ntree = 500, importance = TRUE) 
##                Type of random forest: classification
##                      Number of trees: 500
## No. of variables tried at each split: 2
## 
##         OOB estimate of  error rate: 16.68%
## Confusion matrix:
##          No_Medal Medal class.error
## No_Medal     1874    80  0.04094166
## Medal         317   109  0.74413146
rf_probs <- predict(rf_model, newdata = test_m, type = "prob")[, "Medal"]
rf_pred  <- predict(rf_model, newdata = test_m, type = "response")

rf_cm <- confusionMatrix(rf_pred, test_m$medal_binary, positive = "Medal")
rf_cm
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction No_Medal Medal
##   No_Medal      468    76
##   Medal          20    30
##                                           
##                Accuracy : 0.8384          
##                  95% CI : (0.8063, 0.8671)
##     No Information Rate : 0.8215          
##     P-Value [Acc > NIR] : 0.1542          
##                                           
##                   Kappa : 0.3051          
##                                           
##  Mcnemar's Test P-Value : 1.984e-08       
##                                           
##             Sensitivity : 0.28302         
##             Specificity : 0.95902         
##          Pos Pred Value : 0.60000         
##          Neg Pred Value : 0.86029         
##              Prevalence : 0.17845         
##          Detection Rate : 0.05051         
##    Detection Prevalence : 0.08418         
##       Balanced Accuracy : 0.62102         
##                                           
##        'Positive' Class : Medal           
## 
rf_roc <- roc(response = test_m$medal_binary, predictor = rf_probs,
              levels = c("No_Medal", "Medal"), direction = "<")
cat("Random Forest AUC:", round(auc(rf_roc), 3), "\n")
## Random Forest AUC: 0.795

Feature Importance

So which variables is the random forest actually leaning on?

importance_df <- as.data.frame(importance(rf_model)) %>%
  rownames_to_column("feature") %>%
  arrange(desc(MeanDecreaseGini))

importance_df
##               feature  No_Medal      Medal MeanDecreaseAccuracy
## 1 country_medal_count 43.405925  79.456958            75.156968
## 2                 bmi 28.958916  -3.527162            30.096129
## 3                 age  6.290835  -1.749097             4.605189
## 4              height 25.879523  -1.160072            27.450039
## 5              weight 27.047223 -11.636999            27.949419
## 6        weight_class 35.902754 -18.892882            34.395820
## 7                 era  4.369466   5.509961             6.629381
## 8                 sex  9.080627   3.203887            11.163930
##   MeanDecreaseGini
## 1       187.456757
## 2       103.232336
## 3        89.598863
## 4        81.090303
## 5        74.839551
## 6        54.942566
## 7        11.938028
## 8         6.920498
ggplot(importance_df, aes(x = reorder(feature, MeanDecreaseGini), y = MeanDecreaseGini)) +
  geom_col(fill = "steelblue") +
  coord_flip() +
  labs(title = "Random Forest Feature Importance", x = "Feature", y = "Mean Decrease in Gini") +
  theme_minimal()

ROC Curve Comparison

Side by side, here’s how the two models stack up.

plot(log_roc, col = "darkorange", main = "ROC Curve: Logistic Regression vs Random Forest")
plot(rf_roc, col = "steelblue", add = TRUE)
legend("bottomright",
       legend = c(paste0("Logistic (AUC=", round(auc(log_roc), 3), ")"),
                  paste0("Random Forest (AUC=", round(auc(rf_roc), 3), ")")),
       col = c("darkorange", "steelblue"), lwd = 2)

Metrics Summary

And the numbers, all in one place:

metrics_summary <- tibble(
  model = c("Logistic Regression", "Random Forest"),
  accuracy = c(log_cm$overall["Accuracy"], rf_cm$overall["Accuracy"]),
  sensitivity = c(log_cm$byClass["Sensitivity"], rf_cm$byClass["Sensitivity"]),
  specificity = c(log_cm$byClass["Specificity"], rf_cm$byClass["Specificity"]),
  auc = c(auc(log_roc), auc(rf_roc))
)
metrics_summary
## # A tibble: 2 × 5
##   model               accuracy sensitivity specificity   auc
##   <chr>                  <dbl>       <dbl>       <dbl> <dbl>
## 1 Logistic Regression    0.837       0.245       0.965 0.806
## 2 Random Forest          0.838       0.283       0.959 0.795

Discussion and Recommendations

So what did I actually learn? Both models land around an AUC of 0.80, solidly above chance. That tells me athlete profile data, combined with country-level program information, really does carry a signal about who’s likely to medal. This isn’t noise.

The single biggest predictor, by a wide margin, is country program strength — essentially a country’s historical medal count in weightlifting. Athletes coming out of established programs (China, the former Soviet Union, Bulgaria, Poland, USA) are systematically more likely to medal. Makes sense once you think about it: those countries bring decades of coaching infrastructure and a deep competitive pool behind every athlete they send.

After that, BMI, weight class, and weight round out the next tier. Body composition relative to weight class clearly matters — an athlete with a favorable strength-to-mass ratio for their class tends to be more competitive.

One thing worth being upfront about: sensitivity sits around 25–30% for both models. In plain terms, they’re much better at correctly ruling out non-medalists (specificity near 96%) than at confidently picking out who will medal. Given the medal rate is only about 18% of the field, and elite competition is inherently unpredictable, that’s not a shock.

A few limitations deserve a mention too. country_medal_count looks backward — it summarizes history, not something you’d know for a brand-new athlete with no international track record yet. The dataset also only includes athletes who already qualified for the Olympics; there’s no data on the ones who didn’t make it. So really, the model tells us “who medals among Olympians,” not “who becomes an Olympian” in the first place. And records stretching back to 1896 don’t map cleanly onto today’s training, nutrition, or competitive standards.

Given all that, here’s what I’d tell stakeholders:

  1. Put money into program infrastructure, not just individual scouting. Country-level effects dominate the prediction, which says coaching and training systems are the bigger lever to pull.
  2. Use BMI and weight-class fit as a scouting screen for recruits — paired with real technique assessment, not as a replacement for it.
  3. Bring in data on athletes who never qualified (junior or national championship results, say) so future versions of this model can predict who becomes a medal contender, not just who medals after already getting there.

Conclusion

Stepping back, this project took a talent-identification angle on Olympic weightlifting. I built and compared two classifiers in RStudio, logistic regression and random forest, to predict medal outcomes from athlete and country-level data. Both landed in the same ballpark, around an AUC of 0.80, with national program strength and body composition doing most of the predictive work. Those results back two concrete recommendations: invest in program infrastructure, and use weight-class-specific scouting criteria.

None of this would have been nearly as rigorous in Tableau alone. RStudio’s modeling tools — glm, randomForest, pROC — gave me a reproducible, interpretable way to actually test predictions rather than just visualize patterns, which fits a prediction-driven objective like this one much better.