How do I count thee? Let me count the ways?

Does state life expectancy correlate with political party voting?

Overview       Do people in red states live shorter lives?       This project examines whether state level life expectancy is st...

Thursday, July 23, 2026

Babe Ruth: How good was he? Mantle versus Ruth

      As a kid, my baseball hero was Mickey Mantle. In his 18-year Major League career, he amassed tremendous statistics, and surely he was one of the game's all-time greats. When Mantle retired, he was third on the all-time career home run list with 536 home runs, trailing only Babe Ruth (714) and Willie Mays (587). Unfortunately, he suffered numerous injuries, leaving fans to wonder what his numbers would look like if he had enjoyed a healthier career.

      For decades many people considered (and still consider) Babe Ruth to be the greatest all-around baseball player, because he was a star pitcher and a star hitter. The counter-argument to the Mantle "what-if" is: what if Ruth had not spent his first six years as a pitcher and instead had been a full-time hitter? (Similarly, what if Ted Williams had not spent some of his prime years in the military.)

      Nowadays Shohei Ohtani is certainly a once-in-a-generation two-way star. However, he spent his early prime years in Japan's professional baseball, and he may not have enough time left in US baseball to amass these historic lifetime stats. Aaron Judge is another exceptional modern hitter; hopefully, he stays healthy so we can watch his career numbers climb.

      To compare Mantle and Ruth’s lifetime stats objectively, I turned to data. The Sean Lahman baseball dataset contains Major League player stats back to 1871, and is available in the R library lahman. This was a good opportunity to practice data manipulation using dplyr.

      Ruth played 102 more games than Mantle over a 22-year career (2,503 versus 2,401). His lifetime statistics eclipse Mantle's in every offensive category except stolen bases — though Ruth did steal home 10 times and hit nearly twice as many triples, suggesting he was faster than most realize. His lifetime batting average was an impressive .342, sitting just behind Tris Speaker and Ted Williams.

      Lineup protection played a massive role for both men. Ruth generally batted directly before Lou Gehrig, and Mantle batted directly before Yogi Berra. While both Ruth and Mantle drew plenty of walks, presumably pitchers rarely chose to walk them intentionally just to face Gehrig or Berra.

      I was also curious about their defensive metrics. For a fair comparison, I filtered the data to isolate only their outfield appearances. There is no value in comparing an outfield throw to an assist Ruth made on a comeback grounder while pitching. Similarly, I excluded Mantle's infield appearances; he primarily played first base in his final two seasons and filled in briefly at other infield spots early in his career.

      Ruth played 222 more outfield games than Mantle (2,241 versus 2,019). Their total outfield putouts were nearly identical, likely because Mantle played centerfield and covered more ground. However, Ruth recorded nearly twice as many outfield assists (204 versus 117), which aligns with the arm strength expected of a former pitcher. While Ruth's 204 assists don't match stars such as Roberto Clemente's 266, it remains a highly respectable number.

      Ultimately, the data shows that even if you completely ignore his pitching career, Babe Ruth built an outstanding, standalone career as both a hitter and a fielder.


              
          MANTLE  RUTH
games       2401  2503
at_bats     8102  8398
runs        1677  2174
hits        2415  2873
doubles      344   506
triples       72   136
home_runs    536   714
BA         0.298 0.342
rbi         1509  2217
sb           153   123
bb          1733  2062
so          1710  1330
games_of    2019  2241
putouts     4438  4444
assists      117   204

      Here is my R code:


library(Lahman)
library(tidyverse)

data(People)
which(People$nameLast == "Mantle")   # 11773,  playerID = mantlmi01
which(People$nameLast == "Ruth")   # 16674 , playerID = ruthba01 
data(Batting)    
data(Fielding)   # for fielding, want games in outfield POS == 'OF'

df_fielding <- Fielding %>% 
  filter(playerID %in% c('mantlmi01', 'ruthba01'), POS == 'OF') %>% 
  group_by(playerID) %>% 
  summarize(games_of = sum(G, na.rm = TRUE), putouts = sum(PO, na.rm = TRUE),
    assists = sum(A, na.rm = TRUE)
  )

df <- Batting %>% 
  filter(playerID %in% c('mantlmi01', 'ruthba01')) %>% 
   group_by(playerID) %>% 
   summarize(games = sum(G, na.rm = TRUE), at_bats = sum(AB, na.rm = TRUE),
                runs = sum(R, na.rm = TRUE), hits = sum(H, na.rm = TRUE),
                doubles = sum(X2B, na.rm = TRUE), 
                triples = sum(X3B, na.rm = TRUE), home_runs = sum(HR, na.rm = TRUE),
                BA = round(hits / at_bats,3), 
                rbi = sum(RBI, na.rm = TRUE),  sb = sum(SB, na.rm = TRUE),
                bb = sum(BB, na.rm = TRUE),  so = sum(SO, na.rm = TRUE),
                )  %>%
  merge(df_fielding, by = "playerID", all.x = TRUE) 

df_transposed <- as.data.frame(t(df))
colnames(df_transposed) <- c("MANTLE", "RUTH")
df_transposed <- df_transposed[-1, ]
df_transposed

End

Thursday, July 16, 2026

How much have prices increased?

      The Consumer Price Index (CPI) is a widely used measure of the prices of goods and services purchased by households. It’s the primary tool for tracking inflation and changes in the cost of living over time. The index is built from monthly price collections on a “basket” of goods and services from a sample of retail and service establishments. Historical CPI data is easy to download.

      A key feature of the CPI is that prices are adjusted for quality changes. If the price of a car’s side mirror rises by $200, but $120 of that increase reflects the mirror becoming “smart” rather than “dumb,” only the remaining $80 is counted as inflation. Similarly, if the price of a medical procedure rises because new equipment improves the quality of care, the portion attributable to improved quality is removed from the inflation calculation. Consumers still pay for these quality improvements, whether they want them or not, so in many cases the CPI understates pure cost increases.

      There are eight major categories of the CPI, and each category has its own index: Food & Beverages, Housing, Apparel, Transportation, Medical Care, Recreation, Education & Communication, and Other. These are weighted to form the overall CPI, with the largest weights as Housing at about 44% of the total, Transportation at 17%, Food & Beverages at 14%, and Medical Care at 8%. (Each of these is further sub-divided into its own index; for example, Other includes Personal Care, and Personal Care has seprate indices for Cosmetics, Perfume, Bath, and Nail Preparations.) Of course the weights will not reflect your percentages of what you buy.

      Downloading historical CPI data from FRED (Federal Reserve Bank of St. Louis) was easier than I expected. You need an API key which you can get from https://fredaccount.stlouisfed.org/login

      Here is a line graph showing cumulative CPI growth for the overall CPI and each of the eight major categories through June 2026, indexed to December 2016 = 1.00. The overall CPI has risen 37.1% since December 2016. Housing (rent, insurance, energy, etc.) has increased the most at 44.8%. Transportation (vehicle purchases, fuel, maintenance, insurance, public transit fares) is next at 43.5%. Food (groceries and restaurants) is up 39.9%. Medical Care is lower at 25.9%. Medical Care includes out‑of‑pocket spending on providers, hospitals, and insurance, but excludes employer‑paid and government‑paid health insurance premiums.

      The CPIs exclude a lot of things like the quality changes I mentioned above, and other items that you may pay but that the government does not classify as personal expenses. We all feel the cost of our groceries going up - I discussed this previously in groceries . Each category has its own reasons why it is increasing; I will leave that discussion to the economists, except to say that the cost of energy affects a lot of items in the cost of production and delivery.

      Like many broad measures, the CPI is an attempt to estimate the overall cost of goods and services. But what ultimately matters to you is the actual cost of the things you buy.

      Here is my R code:



library(tidyquant)
library(dplyr)
library(tidyr)
library(lubridate)
library(ggplot2)
library(ggrepel)   # repel overlapping text labels

# Set your API environment variable
Sys.setenv(FRED_API_KEY = "xxxx")


# Define the official FRED database tracking codes
cpi_series <- c(
  "CPIAUCSL", "CPIFABSL", "CPIHOSSL", "CPIAPPSL", 
  "CPITRNSL", "CPIMEDSL", "CPIRECNS", "CPIEDUNS", "CPIOGSNS"
)

# Download and clean data vectors
raw_data <- tq_get(cpi_series, get = "economic.data")   # get from FRED

# 2. Clean, Filter, and Perform Group-Indexing
cpi_processed <- raw_data %>%
  mutate(
    Year  = year(date),
    Month = month(date)
  ) %>%
  # Keep all Decembers from 2015 onward OR strictly isolate June 2026
  filter((Month == 12 & Year >= 2015) | (Year == 2026 & Month == 6)) %>%
  # Convert raw tracking codes into readable titles
  mutate(Category = case_when(
    symbol == "CPIAUCSL" ~ "Overall CPI",
    symbol == "CPIFABSL" ~ "1. Food & Bev",
    symbol == "CPIHOSSL" ~ "2. Housing",
    symbol == "CPIAPPSL" ~ "3. Apparel",
    symbol == "CPITRNSL" ~ "4. Transportation",
    symbol == "CPIMEDSL" ~ "5. Medical Care",
    symbol == "CPIRECNS" ~ "6. Recreation",
    symbol == "CPIEDUNS" ~ "7. Education & Comm",
    symbol == "CPIOGSNS" ~ "8. Other Goods"
  )) %>%
  # Chronologically sort each group, then anchor base-100 to the first row (Dec 2015)
  group_by(Category) %>%
  arrange(date, .by_group = TRUE) %>%
  mutate(Indexed_Value = (price / first(price)) * 100) %>% 
  ungroup() %>%
  # Convert Timeline to ordered categories for a clean discrete X-Axis
  mutate(Period = if_else(Month == 6, paste0(Year, " (June)"), as.character(Year))) %>%
  mutate(Period = factor(Period, levels = unique(Period[order(date)])))

# 3. Isolate final data point rows for the text tags
label_data <- cpi_processed %>%
  group_by(Category) %>%
  filter(date == max(date)) %>%
  ungroup()

# 4. Generate the Chart with the Categorical String Axis Baseline
ggplot(cpi_processed, aes(x = Period, y = Indexed_Value, color = Category, group = Category)) +
  geom_line(aes(linewidth = ifelse(Category == "Overall CPI", 1.5, 0.8))) +
  geom_point(size = 2) +
  
  # Non-overlapping direct text labels
  geom_text_repel(
    data = label_data,
    aes(label = paste0(Category, " (", round(Indexed_Value, 1), ")")),
    nudge_x = 0.5,             
    direction = "y",          
    hjust = 0,                
    segment.color = "grey50", 
    segment.size = 0.4,
    force = 2,
    fontface = "bold",
    size = 4   # millimters
  ) +
  
  # High-contrast visual color mapping matrix
  scale_color_manual(values = c(
    "Overall CPI"         = "#000000", 
    "1. Food & Bev"       = "#E64B35", 
    "2. Housing"          = "#56B4E9", 
    "3. Apparel"          = "#009E73", 
    "4. Transportation"   = "#4D8805", 
    "5. Medical Care"     = "#0072B2", 
    "6. Recreation"       = "#D55E00", 
    "7. Education & Comm" = "#CC79A7", 
    "8. Other Goods"      = "#999999"  
  )) +
  
  # Format plot margins to prevent label cropping
  scale_x_discrete(expand = expansion(mult = c(0.05, 0.35))) +
  scale_linewidth_identity() + 
  labs(
    title = "10.5-Year Cumulative CPI Growth Comparison",
    subtitle = "Base Index: December 2016 = 100",
    x = "Reporting Period",
    y = "Indexed Value (Relative to 100)"
  ) +
  theme_minimal(base_size = 12) +
  theme(
    legend.position = "none",          
    panel.grid.minor = element_blank(),
    text = element_text(face = "bold"),
    plot.title = element_text(face = "bold", size = 14),
    axis.text = element_text(face = "bold")
  )

# 5. Save the final graphic output file
# ggsave("cpi_growth_comparison.png", width = 12, height = 7, dpi = 300, bg = "white")

# Extract and print the final 10.5-year cumulative values
final_column_summary <- cpi_processed %>%
  filter(date == max(date)) %>%
  select(Category, Indexed_Value) %>%
  mutate(Indexed_Value = round(Indexed_Value, 2)) %>%
  arrange(desc(Indexed_Value)) # Sorts from highest inflation to lowest

print(as.data.frame(final_column_summary))


End

Saturday, July 11, 2026

Testing racial predictions with BISG

      Regulators expect insurance companies, banks, and others not to discriminate by race, but businesses are not allowed to collect the race of their customers. One solution for testing discrimination in the aggregate is to use U.S. Census data to estimate the probable race of each customer based on the customer’s name and geographic location. (The Census race data is self identified and excludes a multi racial category; I will leave it to the sociologists to discuss issues with these things and whether Hispanic is a race.)

      The common algorithm for estimating race from names and geography is BISG (Bayesian Improved Surname Geocoding). As an example, suppose there is a customer Mary Johnson who lives in Essex County, NJ. If we want to predict her race only from her surname, we use nationwide Census probabilities such as P(Black ∣ Johnson) = 0.3441, P(White ∣ Johnson) = 0.5438, etc.

      To improve the prediction, we add her county, using geographic weights such as the relative density of each race in Essex compared with nationwide. The resulting probability P(Black | Johnson & Essex) = (Black national prob × Black geographic weight) / ∑race (race national prob × race geographic weight) = 0.7393.

      Finally, we add her first name. Under a Naive Bayes assumption that surname, location, and first name are conditionally independent given race, Bayes’ rule allows us to multiply the components: P(race ∣ surname & geo & first ) ≈ P(race) x P(surname ∣ race) x P(geo ∣ race) x P(first ∣ race). The probability P(Black ∣ Johnson & Mary & Essex) becomes 0.7158.

      A summary of these probabilities is:

Model Type White Black Hispanic Asian Other
Surname Only .5438 .3441 .0272 .0074 .0774
Surname+Geo .1766 .7393 .0244 .0046 .0551
Surname+First+Geo .2481 .7158 .0042 .0013 .0306

      These calculations can be done with the R package wru (Who are you?). You will need a Census API key to download the Census data, available from https://api.census.gov/data/key_ssignup.html .

      I was interested in measuring the accuracy of the BISG predictions. It is not easy to find real (not simulated) data linking name and race. However, several states make their voter registration records public. I requested and downloaded Florida voting registration data from Harvard Dataverse ( https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/UBIG3F ). I limited my study to a single Florida county, Lee, which had over 500,000 registered voters.

      The Florida voting data required considerable data cleaning. Among the amusing data issues is that some Lee County voters registered with a specialty post office box zip code, or a zip code with a typo, or an out-of-state zip code. The first two were mapped to county level, and the third category was dropped.

      There were 545,187 registered Lee County voters, distributed by race as follows:

Asian Black Hispanic Other White NA
6912 32407 66677 9435 422333 7423

      I performed the same BISG algorithm as above, against the Florida voter registration data for Lee County. The overall accuracy is 88%, but this is misleading because the Florida data is highly imbalanced toward the White race. A better accuracy measure is one that is calculated separately for each race, and I chose Balanced Accuracy = (Sensitivity + Specificity)/2. Sensitivity measures how well the algorithm finds people of that specific race, and specificity measures how well the algorithm avoids misclassifying people of other races into that specific race.

      My Balanced Accuracy results of comparing BISG race predictions versus Florida self-identified race data is as follows:

White Black Hispanic Asian Other
.8226 .7179 .9007 .7251 .5041

      I think these results are pretty good for Hispanic and White, not as good for Black and Asian, and of course pretty poor for Other. Other includes Native American, Native Hawaiian, and others, and are a small overall percentage of the county population.

The citation for the Florida voter data is: Sood, Gaurav, 2017, "Florida Voter Registration Data (2017 and 2022)", https://doi.org/10.7910/DVN/UBIG3F, Harvard Dataverse, V2.

      Here is my R code to predict Mary Johnson.


library(wru)
library(dplyr)

# 1. Create a 1-row dataset for Mary Johnson in Essex County, NJ (FIPS "013")
test_voter <- data.frame(
  surname = "JOHNSON",
  first = "MARY",
  state = "NJ",
  county = "013", # Essex County FIPS
  stringsAsFactors = FALSE
)

# 2. Grab the live Census geographic data for New Jersey
nj_census <- get_census_data(
  key = "xxxx",
  states = "NJ",
  age = FALSE,
  sex = FALSE
)

# 3. Generate the 3 prediction variations
pred_baseline  <- predict_race(voter.file = test_voter, surname.only = TRUE)
pred_geo       <- predict_race(voter.file = test_voter, census.data = nj_census, census.geo = "county")
pred_first_geo <- predict_race(voter.file = test_voter, census.data = nj_census, census.geo = "county", names.to.use = "surname, first")

# 4. Format and clean output to perfectly fit the RStudio layout
combined_predictions <- bind_rows(
  pred_baseline  %>% mutate(model_type = "Surname Only"),
  pred_geo       %>% mutate(model_type = "Surname+Geo"),
  pred_first_geo %>% mutate(model_type = "Surname+First+Geo")
) %>%
  rename(
    White = pred.whi,
    Black = pred.bla,
    Hispan = pred.his,
    Asian = pred.asi,
    Other = pred.oth
  ) %>%
  mutate(across(c(White, Black, Hispan, Asian, Other), ~ round(., 4))) %>%
  select(model_type, surname, first, White, Black, Hispan, Asian, Other)

# 5. Display the output summary
print(combined_predictions, row.names = FALSE)
# End of Mary Johnson prediction
##################################################################################

      Here is my R code to measure the accuracy of the BISG predictions against Lee County voters. A considerable portion of the code is for data cleaning of the Florida voter registration data.


library(data.table)
library(dplyr)
library(stringr)
library(caret)
library(wru)

voter_file_path <- "xxxx/LEE_20220621.txt"
Sys.setenv(CENSUS_API_KEY = "xxxx")

raw_lines <- fread(
  file = voter_file_path,
  sep = NULL, # Read the row as a single string to lock layout shifting
  header = FALSE,
  col.names = "raw_row"
)

wru_ready_data <- raw_lines %>%
  mutate(
    # Extract structural components securely out of text blocks
    surname = str_match(raw_row, "\\d{9}\\s+(\\w+)")[,2],
    first = str_match(raw_row, "\\d{9}\\s+\\w+\\s+(\\w+)")[,2],
    zipcode = str_extract(raw_row, "\\b3\\d{4}\\b"),
    raw_race = str_match(raw_row, "\\b([MF])\\s+(\\d)\\b")[,3]
  ) %>%
  # Filter out empty records or text extraction errors
  filter(!is.na(surname) & !is.na(raw_race)) %>%
  # Map Florida's state administrative codes to wru categories
  mutate(
    true_race = case_when(
      raw_race == "5" ~ "white",
      raw_race == "3" ~ "black",
      raw_race == "4" ~ "hispanic",
      raw_race == "2" ~ "asian",
      raw_race %in% c("1", "6", "7") ~ "other",
      TRUE ~ NA_character_ # Automatically standardizes missing responses to NA
    )
  )

# Clean data of invalid zip codes and reroute PO Boxes to residential ZCTAs
wru_ready_data <- wru_ready_data %>%
  # Exclude out-of-state ZIP codes (Georgia 30xxx/31xxx, Alabama 35xxx/36xxx)
  filter(!substr(zipcode, 1, 2) %in% c("30", "31", "35", "36")) %>%
  # Map missing/typos and PO boxes directly to matching valid residential ZCTAs
  mutate(
    zipcode = case_when(
      # Defunct / Typos / Missing data map to Central Lee County Baseline
      zipcode %in% c("fl-NA", "33929", "33932", "33945", "33970") ~ "33901",
     
      # PO Box explicit routing to residential census counterparts
      zipcode == "33902" ~ "33901",  # Fort Myers PO Box -> Fort Myers Residential
      zipcode == "33906" ~ "33907",  # Fort Myers PO Box -> South Fort Myers Residential
      zipcode == "33910" ~ "33904",  # Cape Coral PO Box -> Cape Coral Residential
      zipcode == "33915" ~ "33919",  # Fort Myers PO Box -> Cypress Lake Residential
      zipcode == "33918" ~ "33903",  # N. Fort Myers PO Box -> N. Fort Myers Residential
      zipcode == "33994" ~ "33928",  # Bonita Springs PO Box -> Estero/Bonita Residential
      zipcode == "34106" ~ "34102",  # Naples PO Box -> Naples Residential
      zipcode == "34133" ~ "34135",  # Bonita Springs PO Box -> Bonita Residential
      zipcode == "34136" ~ "34135",  # Bonita Springs PO Box -> Bonita Residential
     
      TRUE ~ zipcode # Keep all other valid residential ZCTAs as they are
    ),
    # Ensure wru recognizes the geographic county boundary for the fallback imputation
    county = "12071" # FIPS code for Lee County, FL
  )

nrow(wru_ready_data)
table(wru_ready_data$true_race, useNA = "always")

lee_county_test <- as.data.table(wru_ready_data) %>%
  filter(!is.na(true_race)) %>%
  mutate(
    state = "fl",
    surname = as.character(surname),
    first = as.character(first),
    zcta = as.character(zipcode) # Required column label mapping for ZCTA use
  )

florida_zcta_layers <- wru::get_census_data(
  key = Sys.getenv("CENSUS_API_KEY"),
  state = "FL",
  age = FALSE,
  sex = FALSE,
  census.geo = "zcta",
  county.list = NULL # Required parameter boundary choice for independent ZCTAs
)

predicted_lee_county <- wru::predict_race(
  voter.file = lee_county_test,
  census.surname = TRUE,
  surname.only = FALSE,
  census.geo = "zcta",
  census.data = florida_zcta_layers,
  impute.missing = TRUE, # Forces wru to fall back to County priors if needed
  skip_bad_geos = TRUE
)

# "53233 (10.1%) individuals' last names were not matched, but ZCTA's baseline geographic racial demographics to calculate the prediction.

evaluation_ready_data <- predicted_lee_county %>%
  rowwise() %>%
  mutate(
    # Isolate which of the five columns held the highest posterior probability string
    max_prob_col = c("white", "black", "hispanic", "asian", "other")[
      which.max(c(pred.whi, pred.bla, pred.his, pred.asi, pred.oth))
    ]
  ) %>%
  ungroup() %>%
  # Convert fields into factor arrays to prevent tracking sequence errors
  mutate(
    true_race = factor(true_race, levels = c("white", "black", "hispanic", "asian", "other")),
    predicted_map = factor(max_prob_col, levels = c("white", "black", "hispanic", "asian", "other"))
  )

# Compute the final metric confusion matrix tracking loop
accuracy_report <- confusionMatrix(
  data = evaluation_ready_data$predicted_map,
  reference = evaluation_ready_data$true_race
)

print(accuracy_report$overall["Accuracy"]) # Overall model precision
print(accuracy_report$byClass[, "Balanced Accuracy"]) # Success broken down per group

End