2.4 K-Way Tables
These notes rely on PSU STATS 504 course notes.
A k-way table has k independent variables. Each cell in the k-way table has joint probability \(\pi_{ij..k}\). You can collapse a k-way table along a dimension to create a marginal table and the resulting joint probabilities in the marginal table are referred to as marginal distributions. Alternatively, you can consider a single level of a dimension, creating a conditional distribution.
Here are two case studies to illustrate the concepts.
Study 1: Death Penalty. Data is collected for n = 326 murder cases along three dimensions: defendant’s color, victim’s color, and whether or not the defendant received the death penalty.
deathp <- array(c(19, 132, 11, 52, 0, 9, 6, 97), dim=c(2,2,2))
dimnames(deathp) <- list(
DeathPen = c("yes","no"),
Defendant=c("white","black"),
Victim=c("white","black")
)
ftable(deathp, row.vars = c("Defendant", "Victim"), col.vars = "DeathPen")## DeathPen yes no
## Defendant Victim
## white white 19 132
## black 0 9
## black white 11 52
## black 6 97Study 2: Boy Scouts. Data is collected for n = 800 Boy Scouts along three dimensions: socioeconomic status, scout status, and delinquency status.
scout <- expand.grid(
delinquent = c("no","yes"),
scout = c("no", "yes"),
SES = c("low", "med","high")
)
scout <- cbind(scout, count = c(169,42,43,11,132,20,104,14,59,2,196,8))
scout_ct <- xtabs(count ~ ., scout)
ftable(scout_ct, row.vars = c("SES", "scout"), col.vars = "delinquent")## delinquent no yes
## SES scout
## low no 169 42
## yes 43 11
## med no 132 20
## yes 104 14
## high no 59 2
## yes 196 82.4.1 Odds Ratio
In the DeathP study, the marginal odds ratio is 1.18, meaning the odds of death penalty for a white defendant are 1.18 times as high as they are for a black defendant.
deathp_m <- margin.table(deathp, margin = c(1, 2))
deathp_p <- prop.table(deathp_m, margin = 2)
deathp_odds <- deathp_p[1, ] / deathp_p[2, ]
deathp_or <- deathp_odds[1] / deathp_odds[2]
print(deathp_or)## white
## 1.2The conditional odds ratio given the victim is white is 0.68, and 0.79 given the victim is black, meaning the odds of death penalty for a white defendant are 0.68 times as high as they are for a black defendent if the victim is white, and 0.79 times as high as they are for a black defendent if the victim is black.
deathp_m <- margin.table(deathp[, , 1], margin = c(2, 1))
deathp_p <- prop.table(deathp_m, margin = 2)
deathp_odds <- deathp_p[1, ] / deathp_p[2, ]
deathp_or <- deathp_odds[1] / deathp_odds[2]
print(deathp_or)## yes
## 0.68deathp_m <- margin.table(deathp[, , 2], margin = c(2, 1)) + 0.5
deathp_p <- prop.table(deathp_m, margin = 2)
deathp_odds <- deathp_p[1, ] / deathp_p[2, ]
deathp_or <- deathp_odds[1] / deathp_odds[2]
print(deathp_or)## yes
## 0.79Notice above that the second margin table had a zero in one cell. To get an odds ratio in this case, the convention is to add 0.5 to all cells.
Interesting that the marginal odds of a white defendent receiving the death penalty are >1, but the conditional odds of a white defendent receiving the death penalty given the race of the victim are <1. Simpson’s paradox is the phenomenon that a pair of variables can have marginal association and partial (conditional) associations in opposite direction. Another way to think about this is that the nature and direction of association changes due to presence or absence of a third (possibly confounding) variable.
2.4.2 Chi-Square Independence Test
There are several ways to think about “independence” when dealing with a k-way table.
- Mutual independence. All variables are independent from each other, (A, B, C).
- Joint independence. Two variables are jointly independent of the third, (AB, C).
- Marginal independence. Two variables are independent if you ignore the third, (A, B).
- Conditional independence Two variables are independent given the third, (AC, BC).
- Homogeneous associations Conditional (partial) odds-ratios are not related on the value of the third, (AB, AC, BC).
Under the assumption that the model of independence is true, once you know the marginal probability values, we have enough information to estimate all unknown cell probabilities. Because each of the marginal probability vectors must add up to one, the number of free parameters in the model is (I − 1) + (J − 1) + (K −I ). This is exactly like the two-way table
2.4.2.1 Mutual Independence
The simplest model is that ALL variables are independent of one another (A, B, C). The joint probabilities are equal to the product of the marginal probabilities, \(\pi_{ijk} = \pi_i \pi_j \pi_k\). In terms of odds ratios, the model (A, B, C) implies that the odds ratios in the marginal tables A × B, B × C, and A × C are equal to 1.
The chi-square independence test tests whether observed joint frequency counts \(O_{ijk}\) differ from expected frequency counts \(E_{ijk}\) under the independence model \(\pi_{ijk} = \pi_{i+} \pi_{+j} \pi_{k+}\). \(H_0\) is \(O_{ijk} = E_{ijk}\). This is essentially a one-way table chi-squre goodness-of-fit test.
For the Boy Scouts study, \(X^2\) and its associated p-value is
scout_e <- array(NA, dim(scout_ct))
for (i in 1:dim(scout_ct)[1]) {
for (j in 1:dim(scout_ct)[2]) {
for (k in 1:dim(scout_ct)[3]) {
scout_e[i,j,k] <- (
margin.table(scout_ct, 3)[k] *
margin.table(scout_ct, 2)[j] *
margin.table(scout_ct, 1)[i]) /
(sum(scout_ct))^2
}
}
}
scout_df <- (prod(dim(scout_ct)) - 1) - sum(dim(scout_ct) - 1)
scout_x2 <- sum((scout_ct - scout_e)^2 / scout_e)
print(scout_x2)## [1] 215## [1] 0.00000000000000000000000000000000000000000079You can also get X^2 this way, although longer code.
x <- scout %>% group_by(delinquent) %>%
summarize(n = sum(count)) %>% ungroup() %>% mutate(p = n / sum(n))## `summarise()` ungrouping output (override with `.groups` argument)p_d <- x$p
names(p_d) <- x$delinquent
x <- scout %>% group_by(scout) %>%
summarize(n = sum(count)) %>% ungroup() %>% mutate(p = n / sum(n))## `summarise()` ungrouping output (override with `.groups` argument)p_b <- x$p
names(p_b) <- x$scout
x <- scout %>% group_by(SES) %>%
summarize(n = sum(count)) %>% ungroup() %>% mutate(p = n / sum(n))## `summarise()` ungrouping output (override with `.groups` argument)p_s <- x$p
names(p_s) <- x$SES
scout_e <- scout %>%
mutate(e = as.numeric(p_d[delinquent] * p_b[scout] * p_s[SES] * sum(count)))
x2 <- sum((scout_e$count - scout_e$e)^2 / scout_e$e)
print(x2)## [1] 215The deviance statistic is
## [1] 219## [1] 0.00000000000000000000000000000000000000000013Safe to say the mutual independence model does not fit.
chisq.test() does not test the complete independence model for k-way tables. Instead, conduct tests on all \({3 \choose 2} = 3\) embedded 2-way tables. For the Boy Scouts study, you can conduct a series of 2-way chi-sq tests.
SES*scout:
## scout
## SES no yes
## low 211 54
## med 152 118
## high 61 204##
## Pearson's Chi-squared test
##
## data: scout_ses_scout
## X-squared = 172, df = 2, p-value <0.0000000000000002SES*delinquent
## delinquent
## SES no yes
## low 212 53
## med 236 34
## high 255 10##
## Pearson's Chi-squared test
##
## data: scout_ses_delinquent
## X-squared = 33, df = 2, p-value = 0.00000007delinquent*scout
## scout
## delinquent no yes
## no 360 343
## yes 64 33##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: scout_delinquent_scout
## X-squared = 7, df = 1, p-value = 0.009The odds ratio for scout delinquency vs non-scout delinquency is 0.54 with 95% CI (0.347, 0.845), so reject the null hypothesis that boy scout status and delinquent status are independent of one another (OR = 1), and thus that boy scout status, delinquent status, and socioeconomic status are not mutually independent.
## odds ratios for delinquent and scout
##
## [1] 0.54## 2.5 % 97.5 %
## no:yes/no:yes 0.35 0.842.4.2.2 Joint Independence
The joint independence model is that two variables are jointly independent of a third (AB, C). The joint probabilities are equal to the product of the AB marginal probability and the C marginal probability, \(\pi_{ijk} = \pi_{ij} \pi_k\).
The degrees of freedom in the associated chi-sq test are the number of free parameters in the saturated model minus the number of free parameters in the joint independence model. \(df = (IJK-1) - ((IJ-1) + (K-1)) = (IJ-1)(K-1)\).
From the Boy Scouts study, suppose juvenile delinquency is the response variable. Test the null hypothesis that juvenile delinquency is independent of boy scout status and socioeconomic status.
scout_d_bs <- ftable(
scout_ct,
row.vars = c("SES", "scout"),
col.vars = "delinquent"
)
print(scout_d_bs)## delinquent no yes
## SES scout
## low no 169 42
## yes 43 11
## med no 132 20
## yes 104 14
## high no 59 2
## yes 196 8##
## Pearson's Chi-squared test
##
## data: scout_d_bs
## X-squared = 33, df = 5, p-value = 0.000004Notice how the contingency table is constructed to not indicate whether SES and scout are related.
2.4.2.3 Marginal Independence
The marginal independence model is that two variables are independent while ignoring the third (A, B). The joint probabilities are equal to the product of the two marginal probabilities, \(\pi_{ij} = \pi_{i} \pi_j\), just like a 2-way table.
The degrees of freedom in the associated chi-sq test are the number of free parameters in the saturated model minus the number of free parameters in the marginal independence model. \(df = (I - 1) + (J - 1) = (I-1)(J-1)\).
From the Boy Scouts study, suppose juvenile delinquency is the response variable. Test the null hypothesis that juvenile delinquency is independent of boy scout status and socioeconomic status.
## SES
## scout low med high
## no 211 152 61
## yes 54 118 204##
## Pearson's Chi-squared test
##
## data: ct
## X-squared = 172, df = 2, p-value <0.00000000000000022.4.2.4 Conditional Independence
The conditional independence model is that two variables are independent given a third (AB, AC). The joint probabilities are equal to the product of the conditional probabilities (B|A and C|A) and the marginal probability of A, \(\pi_{ijk} = \pi_{j|i}\pi_{k|i}\pi_{i++}\).
You can test for conditional independence (AB, AC) by individually testing the independence (B, C) for each level of A. \(X^2\) and \(G^2\) are the sum of the individual values and the degrees of freedom are \(I(J-1)(K-1)\). A better way to do it is with the Cochran-Mantel-Haenszel Test.
The degrees of freedom in the associated chi-sq test are the number of free parameters in the saturated model minus the number of free parameters in the joint independence model. \(df = (IJK-1) - ((IJ-1) + (K-1)) = (IJ-1)(K-1)\).
From the Boy Scouts study, the section on joint independence found that scouting and SES were jointly independent of delinquency, so it must be that either scouting is independent of delinquency, or SES is independent of delinquency, or both are independent of delinquency. The marginal test below establishes the independance of (scout, delinquent).
## delinquent
## scout no yes
## no 360 64
## yes 343 33##
## Pearson's Chi-squared test
##
## data: ct
## X-squared = 7, df = 1, p-value = 0.006The odds ratio 0.54 is a strong negative relationship between boy scout status and delinquency - boy scouts are 46% less likely (on the odds scale) to be delinquent than non-boy scouts.
ct <- xtabs(count ~ scout + delinquent, scout)
p <- prop.table(ct, margin = 1)
(or <- (p[2, 2] / p[2, 1]) / (p[1, 2] / p[1, 1]))## [1] 0.54However, one may ask is scouting and delinquency conditionally independent given SES? You can answer this question with three chi-square tests, one for each level of SES.
## delinquent
## scout no yes
## no 169 42
## yes 43 11##
## Pearson's Chi-squared test
##
## data: ct
## X-squared = 0.006, df = 1, p-value = 0.9## delinquent
## scout no yes
## no 132 20
## yes 104 14##
## Pearson's Chi-squared test
##
## data: ct
## X-squared = 0.1, df = 1, p-value = 0.8## delinquent
## scout no yes
## no 59 2
## yes 196 8## Warning in chisq.test(ct, correct = FALSE): Chi-squared approximation may be
## incorrect##
## Pearson's Chi-squared test
##
## data: ct
## X-squared = 0.05, df = 1, p-value = 0.8You can add up the three \(X^2\) values and run a chi-sq test with \(df = 3\). p = 0.984 indicating that the conditional independence model fits extremely well.
## [1] 0.98The other way of testing for independence is the Cochran-Mantel-Haenszel test.
The Cochran-Mantel-Haenszel (CMH) test statistic
\[M^2 = \frac{\left[ \sum_k{(n_{11k} - \mu_{11k})}\right]^2}{\sum_k{Var(n_{11k})}}\]
For the Boy Scouts study, the test is
##
## Mantel-Haenszel chi-squared test without continuity correction
##
## data: scout_ct
## Mantel-Haenszel X-squared = 0.008, df = 1, p-value = 0.9
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.6 1.6
## sample estimates:
## common odds ratio
## 0.98X2= 0.008 (df = 1, p-value = 0.9287) indicating the conditional independence model is a good fit for this data, so do not reject the null hypothesis.
2.4.2.5 Homogenous Associations
Whereas the conditional independence model, (AB, AC) requires the BC odds ratios at each level of A to equal 1, the homogeneous association model, (AB, AC, BC), requires the BC odds ratios at each level of A to be identical, but not necessarily equal to 1.
The Breslow-Day statistic is
\[X^2 = \sum_{ijk}{\frac{\left(O_{ijk} - E_{ijk}\right)^2}{E_{ijk}}}\]
##
## Breslow-Day test on Homogeneity of Odds Ratios
##
## data: scout_ct
## X-squared = 0.2, df = 2, p-value = 0.9