2by3by3 Chi-square test thanks for reading my question.
I planned to do a 2by2by3 Chi-square test , and ask how

*

*variable A (condition, 2 levels: Cued and Uncued) and

*variable B (group, 3 levels: English, Spa-Eng, and Chi-Eng) affected

*learning of variable C--the number of labels (label, 3 levels: 2 labels, 1 label, and 0 labels).

See the photo below. The sum of variable C should be 16 per row.
I intended to use R function chisq.test(). But I cannot find a way to deal with 3-way data. Can anyone help with this issue? Thanks a lot.


The Pearson $X^2$ test is the score test for association vs independence in a two-way table.  A corresponding test for multiway tables would be the score test for a saturated model vs independence in a loglinear model. In R, the loglm function in the MASS package is one way to do it.
For example, from the help page for loglm
> library(MASS)
> ## the data -- a four-way table
> minn38a <- xtabs(f ~ ., minn38)
> str(minn38a)
 'xtabs' int [1:3, 1:4, 1:7, 1:2] 53 163 309 13 28 38 7 30 17 76 ...
 - attr(*, "dimnames")=List of 4
  ..$ hs : chr [1:3] "L" "M" "U"
  ..$ phs: chr [1:4] "C" "E" "N" "O"
  ..$ fol: chr [1:7] "F1" "F2" "F3" "F4" ...
  ..$ sex: chr [1:2] "F" "M"
 - attr(*, "call")= language xtabs(formula = f ~ ., data = minn38)
> ## the model
> loglm(~ phs + hs+fol+sex, minn38a)
Call:
loglm(formula = ~phs + hs + fol + sex, data = minn38a)

Statistics:
                      X^2  df P(> X^2)
Likelihood Ratio 3711.914 155        0
Pearson          4161.597 155        0
```

