In Chapter 7 of the textbook Practical Methods for Design and Analysis of Complex Surveys, testing procedures for one-way and two-way tables are presented. The test statistics include a design-based Wald test statistic and various Rao-Scott adjusted test statistics. These test statistics are constructed to account for the sampling complexities, such as stratification, clustering and weighting, thus providing asymptotically valid testing procedures for complex surveys.
In Training Key 250, a design-based Wald test statistic of independence hypothesis for a two-way table is calculated and compared with an unadjusted Pearson test statistic, thus reproducing the results of Example 7.3. The Occupational Health Care Survey data set is used.
Example 7.3. In this Training Key we will study whether the variables PHYS (physical health hazards of work: None or Some) and PSYCH3 (overall psychic strain classified into three nearly equally sized classes) are associated or not. We use OHC Survey data as an example data set (Download ohc data). The general purpose is to show how the sampling design affects to the test of independence of two categorical variables.
We will demonstrate how to carry out the test of independence of health hazards of work and psychic strain in the OHC Survey setting. The test statistics are introduced in Lehtonen and Pahkinen (2004) Section 7.5, and some of the results of Example 7.3 will be reproduced.
NOTE: Due to the simplicity of the problem, we do not include the interactive analysis option in this training key.
A hypothesis of independence is stated as \(H_0: \: p_{jk} = p_{j+}p_{+k}\) with \(j=1,2\) and \(k=1,2,3\), or, analogously, \(H_0: \: p_{11} - p_{1+}p_{+1} = 0\) and \(p_{12} - p_{1+} p_{+2} = 0\).
Basically, two different testing approaches are available: the Rao-Scott adjusted SRS-based test statistics and the design-based Wald test statistics (for an overview see Section 7.5 of Lehtonen and Pahkinen 2004). Both approaches are illustrated in this Training Key.
The two-way table of PHYS and PSYCH3 is first displayed. The tables include cell frequencies and percentages and total percentages. The figures correspond to those in Table 7.4 of Lehtonen and Pahkinen (2004).
Load ohc data, packages (survey and tidyverse) and get a glimpse of the ohc data.
# Load ohc data and packages survey and tidyverse.
load("ohc.Rdata")
library(survey); library(tidyverse)
glimpse(ohc)## Rows: 7,841
## Columns: 15
## $ CLUSTER <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,...
## $ STRATUM <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,...
## $ SUBJECT <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,...
## $ ID <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, ...
## $ SEX <dbl> 1, 1, 2, 1, 2, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 1, 2, 2, 1, 2,...
## $ AGE <dbl> 38, 23, 23, 18, 20, 20, 54, 61, 45, 27, 34, 41, 32, 61, 29,...
## $ AGE2 <dbl> 1, 1, 1, 1, 1, 1, 2, 2, 2, 1, 1, 1, 1, 2, 1, 2, 2, 1, 2, 1,...
## $ PHYS <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 0, 0, 0, 0,...
## $ CHRON <dbl> 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 0,...
## $ PSYCH <dbl> -0.9977110, -0.4106113, -0.7097656, -0.4106113, -0.4264838,...
## $ PSYCH2 <dbl> 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0,...
## $ PSYCH3 <dbl> 1, 2, 1, 2, 2, 3, 2, 1, 3, 1, 3, 3, 3, 3, 3, 1, 2, 1, 1, 1,...
## $ NEWCLU <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,...
## $ NEWSTR <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,...
## $ AGE4 <dbl> 2, 1, 1, 1, 1, 1, 4, 4, 3, 1, 2, 3, 2, 4, 1, 4, 4, 3, 4, 3,...The sampling design of the OHC Survey is an example of stratified cluster sampling where both one-stage and two-stage sampling have been used. The pedacogical data set constructed for training purposes in this example includes a total of 250 clusters, i.e. industrial establishments, organized in five strata (stratums 2 to 6), and a total of 7841 persons. Stratification is based on the type of industry and cluster size (the number of salaried employees). Clusters having at least 10 employees are included in the OHC example data set. There is variable number of clusters per stratum in the design. The average cluster sample size is 11.2 employees. A more detailed description of the study design and sampling design of the OHC Survey are given in Section 6.1 of Lehtonen and Pahkinen (2004).
Relevant variables for our example for modelling are the following:
Next we will construct a ohc.design object. This ohc.design object will be used for all subsequent analysis commands: NOTE. nest = TRUE.
# Construct a ohc.design object.
# This ohc.design object will be used for all subsequent analysis commands:
# NOTE. nest = TRUE
ohc.design <- svydesign(ids = ~CLUSTER, data = ohc, strata = ~STRATUM, nest = TRUE)The two-way table of PHYS and PSYCH3 is first displayed. The tables include cell frequencies and percentages and total percentages.
paste("Cell frequencies:")## [1] "Cell frequencies:"svytable(~PHYS+PSYCH3, ohc.design) # cell frequencies## PSYCH3
## PHYS 1 2 3
## 0 1785 1716 1629
## 1 910 821 980paste("Row percentages:")## [1] "Row percentages:"prop.table(svytable(~PHYS+PSYCH3, ohc.design), 1) # row percentages## PSYCH3
## PHYS 1 2 3
## 0 0.3479532 0.3345029 0.3175439
## 1 0.3356695 0.3028403 0.3614902paste("Column percentages:")## [1] "Column percentages:"prop.table(svytable(~PHYS+PSYCH3, ohc.design), 2) # column percentages## PSYCH3
## PHYS 1 2 3
## 0 0.6623377 0.6763894 0.6243772
## 1 0.3376623 0.3236106 0.3756228paste("Total percentages:")## [1] "Total percentages:"svytable(~PHYS+PSYCH3, ohc.design, Ntotal=100) # percents## PSYCH3
## PHYS 1 2 3
## 0 22.76495 21.88496 20.77541
## 1 11.60566 10.47060 12.49841Design-based test results for the independence hypothesis are next displayed. Four test for independence are included (see documentation of survey::svychisq):
paste("Rao-Scott Rao-Scott adjusted F-test statistics with second-order correction:")## [1] "Rao-Scott Rao-Scott adjusted F-test statistics with second-order correction:"svychisq(~PHYS+PSYCH3, ohc.design) # test for independence. Note: The default (statistic="F") is the Rao-Scott second-order correction.##
## Pearson's X^2: Rao & Scott adjustment
##
## data: svychisq(~PHYS + PSYCH3, ohc.design)
## F = 7.1377, ndf = 1.9624, ddf = 480.7799, p-value = 0.0009558paste("Rao-Scott Rao-Scott adjusted Pearson chi-square test statistics:")## [1] "Rao-Scott Rao-Scott adjusted Pearson chi-square test statistics:"svychisq(~PHYS+PSYCH3, statistic="Chisq", ohc.design) # test for independence. Note: corrects the Pearson's chi-squared statistic ##
## Pearson's X^2: Rao & Scott adjustment
##
## data: svychisq(~PHYS + PSYCH3, statistic = "Chisq", ohc.design)
## X-squared = 16.569, df = 2, p-value = 0.0007946paste("Design-based Wald test:")## [1] "Design-based Wald test:"svychisq(~PHYS+PSYCH3, statistic="Wald", ohc.design) # test for independence. Note: test is that proposed by Koch et al (1975) and used by the SUDAAN software package. It is a Wald test based on the differences between the observed cells counts and those expected under independence.##
## Design-based Wald test of association
##
## data: svychisq(~PHYS + PSYCH3, statistic = "Wald", ohc.design)
## F = 6.614, ndf = 2, ddf = 245, p-value = 0.001594paste("Adjusted design-based Wald test")## [1] "Adjusted design-based Wald test"svychisq(~PHYS+PSYCH3, statistic="adjWald", ohc.design) # test for independence. Note: adjusted Wald statistics##
## Design-based Wald test of association
##
## data: svychisq(~PHYS + PSYCH3, statistic = "adjWald", ohc.design)
## F = 6.587, ndf = 2, ddf = 244, p-value = 0.001636The next tables provide results from the corresponding SRS-based tests, including the Pearson chi-square test and likelihood ratio test. Here, all the complexities (stratification and clustering) are ignored, clearly providing overly liberal testing results. Thus, these testing procedures give too large \(\chi^2\)-values and therefore too small p-values. This is because the clustering effect is not accounted for in the SRS-based tests of independence, contrary to the design-based testing procedures.
paste("Crosstabulation table:")## [1] "Crosstabulation table:"table(ohc$PHYS, ohc$PSYCH3)##
## 1 2 3
## 0 1785 1716 1629
## 1 910 821 980paste("Chi-square test and p-value:")## [1] "Chi-square test and p-value:"chisq.test(ohc$PHYS, ohc$PSYCH3)##
## Pearson's Chi-squared test
##
## data: ohc$PHYS and ohc$PSYCH3
## X-squared = 16.569, df = 2, p-value = 0.0002524tab <- xtabs(~PHYS + PSYCH3, data = ohc)
library(vcd)
paste("Likelihood ratio test:")## [1] "Likelihood ratio test:"assocstats(tab)## X^2 df P(> X^2)
## Likelihood Ratio 16.500 2 0.00026130
## Pearson 16.569 2 0.00025238
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.046
## Cramer's V : 0.046Pay attention to the following. The Rao-Scott adjusted SRS-based tests and the design-based Wald tests provide a reasonable testing procedure for an independence hypothesis of PHYS and PSYCH3 in the OHC survey setting.