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).


    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)
        X_fit = sm.add_constant(X)
        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):
        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))
    plt.figure(figsize=(10, 8))
    plt.hist(statistics, bins=50);
    plt.axvline(threshold_95, color="red", label="$\\chi^2$ threshold");
    plt.title("Distribution of computed $\\chi^2$ statistic");

Any help appreciated!


enter image description here

  • $\begingroup$ 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. $\endgroup$ – whuber Apr 8 at 14:33
  • $\begingroup$ 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. $\endgroup$ – Arnold Davidson Apr 8 at 14:46
  • $\begingroup$ Why do you use np.var(predictions) for $\sigma^2$? $\endgroup$ – stefgehrig Apr 8 at 15:45
  • $\begingroup$ Your $\sigma$ seems to be 1 $\endgroup$ – stefgehrig Apr 8 at 15:53
  • $\begingroup$ Perhaps my question should be: what is the correct way to carry out this test? (I'm aware I'm not doing something right.) $\endgroup$ – Arnold Davidson Apr 8 at 16:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.