# --- PSM_PROTOCOL.R ---
# R version 4.1.2
# Required packages: MatchIt, ggplot2, tableone

# STEP 1: Propensity Score Estimation
ps_model <- glm(female ~ age + bmi + hypertension + diabetes + prior_fracture,
                data = op_data,
                family = binomial())

# STEP 2: 1:1 Nearest Neighbor Matching
matched_data <- matchit(ps_model,
                        method = "nearest",
                        ratio = 1,
                        caliper = 0.05,
                        replace = FALSE)

# STEP 3: Balance Assessment
# Generate Love plot
love_plot <- love.plot(matched_data, 
                       thresholds = c(0.1, 0.2),
                       colors = c("red", "darkgreen"))

# STEP 4: Outcome Analysis
matched_df <- match.data(matched_data)
psm_model <- clogit(op_diagnosis ~ female + strata(subclass),
                    data = matched_df)

# SENSITIVITY ANALYSIS
# Varying caliper distances
sensitivity <- matchit(ps_model, caliper = c(0.01, 0.1))

# OUTPUT REPORT
sink("PSM_Results.txt")
print(summary(psm_model))
print(bal.tab(matched_data))
sink()