Most appropriate test for comparing multiple changes in proportions?

Question

I'd like to test the significance of changes in proportions of a number of features from the same population.

The data looks like this:

# Construct data
code_table_df <- data.frame(
    YEAR = rep(c(2019, 2020), each = 4),
    code = rep(rep(1:4), 2),
    value = c(61, 14, 15, 13, 50, 21,  3, 11)
)

# Create table
code_table <- xtabs(value ~ YEAR + code, data = code_table_df)
code_table

      code
YEAR    1  2  3  4
  2019 61 14 15 13
  2020 50 21  3 11

A regular Z test of proportions is not appropriate here as I do not want to test if the proportions for codes 1-4 are equal within YEAR. Nor do I want to test if the proportions of YEAR are equal within code.

What I want to check is, are the proportions of the individual codes changing accross years? That is, for codes 1-4, have the proportions changed between 2019 and 2020? Or, in the plot above, is the there a statistically significant difference in height between the dark grey and light grey bar, for each code. For code 3, this would be testing $p_{3,2019}=15/(61+14+15+13)=0.146$ is not equal to $p_{3,2020}=3/(50+21+3+11)=0.035$.

The first thought I had was to use the $\chi^2$ test to test for association between the 2 features, YEAR and code. This gives me

> chisq.test(code_table)

    Pearson's Chi-squared test

data:  code_table
X-squared = 9.016, df = 3, p-value = 0.02908

This indicates there is some evidence for association between the 2 features.

I could also do a proportion test (prop.test() in R) for each level of code, but is there a test that allows me to estimate the differences in proportions year on year for all codes with confidence intervals and also adjust/correct for the fact im performing multiple tests?

Answer 1

For the chi-squared test of homogeneity in R:

TBL = rbind(c(61,14,15,14), c(50,21,3,11));  TBL
     [,1] [,2] [,3] [,4]
[1,]   61   14   15   14
[2,]   50   21    3   11
chi.out = chisq.test(TBL)
chi.out

        Pearson's Chi-squared test

data:  TBL
X-squared = 9.0313, df = 3, p-value = 0.02888

Significant differences detected at 5% level. Look at Pearson residuals: Largest absolute values among them may point the way to significant differences for various codes. Here, Code 3 seems of particular interest.

chi.out$resi
            [,1]      [,2]      [,3]        [,4]
[1,] -0.01015505 -1.198408  1.618984  0.06562061
[2,]  0.01123284  1.325598 -1.790811 -0.07258510

I have found some accounts of the syntax for prop.test to be confusing. The version I use here works without difficulty. Counts are sufficiently small to warrant use of the (default) continuity correction. [I tend to turn it off (with parameter cor=F) for large counts--certainly, for 100 or more.]

prop.test(15, 18)

    1-sample proportions test with continuity correction

data:  15 out of 18, null probability 0.5
X-squared = 6.7222, df = 1, p-value = 0.009522
alternative hypothesis: true p is not equal to 0.5
95 percent confidence interval:
 0.5773525 0.9559302
sample estimates:
         p 
 0.8333333 

I wouldn't declare a significant difference from 1/2 unless the P-value is smaller than $.05/4 = 0.0125,$ but that level is met here.

Similarly, for any additional ad hoc tests.