# ==============================================================================
# BRFSS 2024 DATA CLEANING & DESCRIPTIVE REPORT
# Author: Herve Fossou / Samar Ikram
# Focus: Menopausal Status and Osteoarthritis
# ==============================================================================
# 1. Clear environment memory
rm(list = ls())
# 2. Load required packages
library(tidyverse)
library(haven)
# 3. Optimized Data Ingestion
brfss_subset <- read_xpt(
"LLCP2024.XPT",
col_select = c("SEXVAR", "_LLCPWT", "_STATE", "_PSU", "HAVARTH4", "_AGEG5YR", "_BMI5")
) %>%
# Filter immediately for female respondents (SEXVAR == 2)
filter(SEXVAR == 2)
# 4. Handling Missing Values and Exclusions
brfss_missing_handled <- brfss_subset %>%
filter(
HAVARTH4 %in% c(1, 2), # Keep only explicit 'Yes' or 'No'
`_AGEG5YR` %in% 1:14, # Exclude 'Refused/Missing' age categories
`_BMI5` > 0 & `_BMI5` != 99999 # Exclude missing BMI values
)
# 5. Variable Recoding and Collapsing Categories
brfss_clean <- brfss_missing_handled %>%
mutate(
# Recode Osteoarthritis binary outcome
osteoarthritis = case_when(
HAVARTH4 == 1 ~ "Yes",
HAVARTH4 == 2 ~ "No"
),
osteoarthritis = factor(osteoarthritis, levels = c("No", "Yes")),
# Collapse 14 raw age brackets into 3 biological lifecycle periods
menopause_status = case_when(
`_AGEG5YR` %in% 1:5 ~ "Pre-Menopausal Age (18-44)",
`_AGEG5YR` %in% 6:7 ~ "Peri/Early Post-Menopausal Age (45-54)",
`_AGEG5YR` %in% 8:14 ~ "Post-Menopausal Age (55+)"
),
menopause_status = factor(menopause_status,
levels = c("Pre-Menopausal Age (18-44)",
"Peri/Early Post-Menopausal Age (45-54)",
"Post-Menopausal Age (55+)")),
# Correct the BRFSS 2-decimal implicit multiplier for BMI
bmi = `_BMI5` / 100
) %>%
# Retain only final structural variables for the analytic file
select(`_LLCPWT`, `_STATE`, `_PSU`, osteoarthritis, menopause_status, bmi)
# 6. Applying Clear Variable Metadata Labels
attr(brfss_clean$osteoarthritis, "label") <- "Doctor-Diagnosed Osteoarthritis Status"
attr(brfss_clean$menopause_status, "label") <- "Biological Menopausal Age Categories"
attr(brfss_clean$bmi, "label") <- "Calculated Body Mass Index (kg/m²)"
# ==============================================================================
# SECTION 7: STUDY POPULATION CHARACTERISTICS (DEMOGRAPHICS & COVARIATES)
# ==============================================================================
# 1. Summary of Categorical Variables (Frequencies & Percentages)
cat("\n--- Osteoarthritis Frequencies ---\n")
##
## --- Osteoarthritis Frequencies ---
oa_table <- table(brfss_clean$osteoarthritis)
print(oa_table)
##
## No Yes
## 125226 84300
print(prop.table(oa_table) * 100)
##
## No Yes
## 59.76633 40.23367
cat("\n--- Menopausal Age Categories Frequencies ---\n")
##
## --- Menopausal Age Categories Frequencies ---
meno_table <- table(brfss_clean$menopause_status)
print(meno_table)
##
## Pre-Menopausal Age (18-44) Peri/Early Post-Menopausal Age (45-54)
## 57543 28387
## Post-Menopausal Age (55+)
## 123596
print(prop.table(meno_table) * 100)
##
## Pre-Menopausal Age (18-44) Peri/Early Post-Menopausal Age (45-54)
## 27.46342 13.54820
## Post-Menopausal Age (55+)
## 58.98838
# 2. Summary of Continuous Variable (BMI Mean & Standard Deviation)
cat("\n--- Continuous BMI Summary ---\n")
##
## --- Continuous BMI Summary ---
bmi_summary <- brfss_clean %>%
summarise(
Mean_BMI = mean(bmi, na.rm = TRUE),
SD_BMI = sd(bmi, na.rm = TRUE),
Median_BMI = median(bmi, na.rm = TRUE),
IQR_BMI = IQR(bmi, na.rm = TRUE)
)
print(bmi_summary)
## # A tibble: 1 × 4
## Mean_BMI SD_BMI Median_BMI IQR_BMI
## <dbl> <dbl> <dbl> <dbl>
## 1 28.5 7.13 27.4 8.63
# 3. Visualizations: Bar Chart for Osteoarthritis by Menopausal Status
ggplot(brfss_clean, aes(x = menopause_status, fill = osteoarthritis)) +
geom_bar(position = "fill") +
scale_y_continuous(labels = scales::percent) +
theme_minimal() +
labs(
title = "Osteoarthritis Prevalence Across Menopausal Age Groups",
x = "Biological Menopausal Age Categories",
y = "Percentage (%)",
fill = "Osteoarthritis"
) +
theme(axis.text.x = element_text(angle = 15, hjust = 1))

# 4. Visualizations: Histogram of Body Mass Index (BMI)
ggplot(brfss_clean, aes(x = bmi)) +
geom_histogram(binwidth = 1, fill = "steelblue", color = "white", alpha = 0.8) +
geom_vline(aes(xintercept = mean(bmi, na.rm = TRUE)), color = "red", linetype = "dashed", size = 1) +
xlim(15, 50) + # Limit x-axis to standard ranges for clean viewing
theme_minimal() +
labs(
title = "Distribution of Body Mass Index (BMI) in Study Population",
subtitle = "Red dashed line indicates the population mean",
x = "Calculated BMI (kg/m²)",
y = "Count"
)
