# Breusch-Pagan Test
# Module 8 - The Breusch-Pagan Test
# Dataset: MODULE8_INCOME_SAVINGS_SYNTHETIC

# Breusch-Pagan Test
#
# Module 8 notebook lab. This notebook uses original Ceteris Lab teaching data and does not report real empirical findings.

# Learning goal
# Run and manually audit the Breusch-Pagan workflow.
#
# Dataset: MODULE8_INCOME_SAVINGS_SYNTHETIC. Variables: household_id, income, savings, age, education, family_size.

# %% Cell 3
import pandas as pd
import statsmodels.api as sm
from statsmodels.stats.diagnostic import het_breuschpagan

df = pd.read_csv("/data/module-8/module8_income_savings_heteroskedastic.csv")
X = sm.add_constant(df[["income", "age", "education"]])
model = sm.OLS(df["savings"], X).fit()
lm, lm_pvalue, fvalue, f_pvalue = het_breuschpagan(model.resid, X)
aux = sm.OLS(model.resid ** 2, X).fit()
print({"LM": round(lm, 4), "LM p-value": round(lm_pvalue, 4), "F": round(fvalue, 4), "F p-value": round(f_pvalue, 4)})
print("Manual nR^2:", round(len(df) * aux.rsquared, 4))

# Reflection
# Write two sentences: one sentence explaining what the diagnostic or robust result says, and one sentence explaining a limitation or next step.