To compute the test statistic $\chi^2$ the sum is typically performed over all cells in a table, see e.g. https://en.wikipedia.org/wiki/Pearson%27s_chi-squared_test
The $\chi^2$-distribution is the distribution of a sum of the squares of k independent standard normal random variables.
So, I wonder why the sum to compute the test statistics is not only over the cells that are independent: sum over as many cells that there are degrees of freedom. E.g. in the case of a homogeneity test: just sum $(row-1)(column-1)$ cells.
I stumbled over this question by simulating the $\chi^2$-distribution for some examples from the corresponding null hypothesis. The simulated sample distribution fits the $\chi^2$-distribution better if I only perform the sum over independent cells. However, in textbooks the sum is over all cells.
Thanks for your answer.
Here is the simulation (Is constraint correctly implemented?).
Contingency table for the expectation, if X and Y are independent:
Y
0 | 1 |
0 | 24 | 16 | 40
X 1 | 36 | 24 | 60
-------------------
60 40 100
t=10 #
data = np.array([[35,5],[25,35]]) *t # contingency table
# We want to test for independence, the expectation is
expected = np.array([[24,16],[36,24]]) *t
# simulation under the constraints
sim00 = np.random.binomial(n=100*t, p=24/100, size=10000000)
simulation = np.zeros((10000000, 2, 2))
simulation[:,0,0] = sim00
simulation[:,1,0] = 60*t - sim00
simulation[:,:,1] = np.array([40*t, 60*t]) - simulation[:,:,0]
# sum over all cells
sample_distribution = ((simulation[:,:,:]-expected[:,:])**2/expected[:,:]).sum(axis=1).sum(axis=1)
# the independet cell - you have to comment one out
sample_distribution = ((simulation[:,0,0]-expected[0,0])**2/expected[0,0])
# plot the simulation and the chi2 distribution
x = np.linspace(scipy.stats.chi2.ppf(0.3, df), scipy.stats.chi2.ppf(0.999, df), 100)
plt.plot(x, scipy.stats.chi2.pdf(x, df), 'r-', lw=2, alpha=0.6, label='$\chi^2$-pdf')
plt.hist(sample_distribution, bins=100, density=True, label='sampling')
plt.legend();
Update: The simulation (with the full sum) gives the expected result if the simulation is done with a hypergeometric_distribution.
# https://en.wikipedia.org/wiki/Hypergeometric_distribution#Working_example
nb_simulations=100000
s00 = np.random.hypergeometric(ngood=40*t, nbad=60*t, nsample=60*t, size=nb_simulations)
s = np.ndarray((nb_simulations, 2, 2))
s[:,0,0] = s00
s[:,0,1] = 40*t-s00
s[:,1,0] = 60*t-s00
s[:,1,1] = 60*t-s[:,1,0]
sample_distribution = ((s - expected)**2/expected).sum(axis=1).sum(axis=1)