Post hoc $\chi^2$ test with R

Question

This is a quite simple question but I don't find any good, clear, precise answers: I'm looking for a way to perform post hoc test on a chi$^2$ test.

I have 2 variables : var1 : good/fair/poor and var2: a/b/c

This is the contingency table :

      a   b  c
good 120  70 13
fair 230 130 26
poor  84  83 18

with R :

mat <- matrix(c(120,230,84,70,130,83,13,26,18),3)
dimnames(mat) <- list(c("good","fair","poor"),c("a","b","c"))
mat

I perform $\chi^2$

chisq.test(mat)

I obtain p =0.022, so there is a relation between var1 and var2... but which one? How can I handle this? (contribution analysis?, pairwise prop test?)

Answer 1

I like this question because too often, people do omnibus tests and then don't ask more specific questions about what is happening.

If the goal is to compare "treatments" a, b, and c, I would suggest summarizing the data showing the percentages within each column, so you can see more clearly how they differ. Then to test these comparisons, one simple idea is to do the $\chi^2$ test on each pair of columns:

> for (j in 1:3) print(chisq.test(mat[, -j]))

    Pearson's Chi-squared test

data:  mat[, -j]
X-squared = 0.1542, df = 2, p-value = 0.9258


    Pearson's Chi-squared test

data:  mat[, -j]
X-squared = 4.5868, df = 2, p-value = 0.1009


    Pearson's Chi-squared test

data:  mat[, -j]
X-squared = 9.5653, df = 2, p-value = 0.008374

Since 3 tests are done, a Bonferroni correction is advised (multiply each $P$ value by 3). The last test, where column 3 is omitted, has a very low $P$ value, so you can conclude that the distributions of (good, fair, poor) are different for conditions a and b. Note, however, that condition c does not have much data, and that's largely why the other two results are nonsignificant.

You could use a similar strategy to do pairwise comparisons of the rows.

Answer 2

In case anyone still comes across this ancient thread - the procedure suggested by Guilherme is now implemented in the chisq.posthoc.test package, which also offers specific p-values based on the residuals.

Answer 3

Another way of doing this is by means of pearson standardized residuals, as suggested by Agresti, A. (2007) in his book Categorical Data Analysis section 3.3 Following-up chi-squared tests.

The Pearson standardized residuals ($e_{ij}$) measures how large is the deviation from each cell to the null hypothesis (in this case, independence between row and column's).

Positive residuals indicate positive association between row and column variables. Negative residuals indicate negative association between row and column variables.

Once $e_{ij} \sim N(0,1)$, $e_{ij} > |2|$ are indicative of association. They are obtained this way:

tab2 <- chisq.test(mat)
tab2$stdres
           a        b        c
good  1,0164 -0,71661 -0,60995
fair  1,9643 -1,66201 -0,66782
poor -3,3512  2,68760  1,41203

The category poor is negatively associated with "a", therefore, "a" is associated with few poor people. The category poor is associated positively with "b", therefore, "b" is associated with large number of poor people.