# Chi-squared test for regression residuals

Per Cosma Shalizi, here, equation 11.2, I understood the chi-squared test for regression residuals to involve calculating:

$$\sum_{i=1}^{n} \frac{(y_i - f(x_i, \theta))^2}{\sigma_i^2}$$

which will then be $$\chi^2$$ distributed with $$n - k$$ degrees of freedom.

I wanted to test this by simulation but what I get doesn't seem to make sense (I was expecting to reject the null hypothesis in line with the confidence level).

Note on the simulation, performed n_sims times:

1. Generate a dataset of size n using the true parameter values.
2. Fit a linear regression and calculate the residuals.
3. From the residuals attempt to calculate the $$\chi^2$$ statistic (note I've tried several variations on this to no avail so clearly have some misunderstanding).

Code:

    import matplotlib.pyplot as plt
import numpy as np
import statsmodels.api as sm
from scipy import stats

sigma = 1  # true sigma
beta = 0.5  # true beta
n = 25  # number of data points

def generate_fit_single_regression(beta, sigma, n):
X = np.random.uniform(50, 80, size=n)
y = X * beta + np.random.normal(0, sigma, size=n)
return X, y, sm.OLS(y, X_fit).fit()

df = n - 2
threshold_95 = stats.chi2(df).ppf(0.975)

n_sims = 10_000
statistics = []
for i in range(n_sims):
np.random.seed(i)
X, y, results = generate_fit_single_regression(beta, sigma, n)
predictions = results.predict()
residuals = y - predictions
std_residuals = np.std(residuals)
chi_squared_statistic = np.sum(residuals**2 / np.var(predictions))
statistics.append(chi_squared_statistic)

plt.figure(figsize=(10, 8))
plt.hist(statistics, bins=50);
plt.axvline(threshold_95, color="red", label="$$\\chi^2$$ threshold");
plt.legend();
plt.title("Distribution of computed $$\\chi^2$$ statistic");


Any help appreciated!

Arnie • You seem to confuse $\sigma_i$ with some kind of estimate of it. That's not what Shalizi means. Second, Shalizi explicitly says the degrees of freedom is $1,$ but you have used n - 2 instead: that is, your "$\chi^2$ threshold" is wrong. – whuber Apr 8 at 14:33
• Thanks for the reply, @whuber. Re. $\sigma_i$, what should I use instead? Re. degrees of freedom I thought, given we are summing over $n$ $\chi^2(1)$ random variables we would have a $\chi^2(n-k)$ random variable? Clearly I'm probably just highlighting my confusion again at this point so any instruction welcome. – Arnold Davidson Apr 8 at 14:46
• Why do you use np.var(predictions) for $\sigma^2$? – stefgehrig Apr 8 at 15:45
• Your $\sigma$ seems to be 1 – stefgehrig Apr 8 at 15:53
• Perhaps my question should be: what is the correct way to carry out this test? (I'm aware I'm not doing something right.) – Arnold Davidson Apr 8 at 16:05