Every misskappa estimate carries an asymptotic covariance matrix, and (for the influence-function estimators) a per-subject influence matrix fit$psi satisfying vcov == crossprod(psi) / n^2. That is exactly what a Wald test of a linear hypothesis needs, so kappa_test() and alpha_test() can compare a coefficient across fits and return a standard htest. Only the way the covariance of the contrast is assembled changes between designs:
- one sample: the coefficient’s own variance
- independent samples: the per-fit variances add (block diagonal)
-
paired / same-subject: the cross-fit dependence is read from the stacked influence functions
fit$psi
One sample: is a coefficient equal to a fixed value?
Pass a single fit and a value. This tests whether the (linearly weighted) Conger coefficient for the Gwet (2014) data equals 0.5:
kappa_test(kappa(dat.gwet2014, estimator = "ipw"), coef = "Conger", value = 0.5)
#> Warning: rater pair(s) rater4-rater5 co-observed by only one subject; the
#> corresponding pairwise covariance is degenerate and the standard error
#> unreliable.
#>
#> One-sample Wald test that Conger = 0.5
#>
#> data: fit1
#> X-squared = 0.46985, df = 1, p-value = 0.4931
#> sample estimates:
#> fit1
#> 0.4233051Independent samples: two studies, different subjects
When the fits come from disjoint sets of subjects, set paired = FALSE. The variances add. Here the Conger coefficient is compared between the Gwet (2014) and Klein (2018) studies:
kappa_test(
gwet = kappa(dat.gwet2014, estimator = "ipw"),
klein = kappa(dat.klein2018, estimator = "ipw"),
coef = "Conger", paired = FALSE
)
#> Warning: rater pair(s) rater4-rater5 co-observed by only one subject; the
#> corresponding pairwise covariance is degenerate and the standard error
#> unreliable.
#>
#> Independent-sample test of equal Conger across 2 fits
#>
#> data: gwet, klein
#> X-squared = 0.0019646, df = 1, p-value = 0.9646
#> sample estimates:
#> gwet klein
#> 0.4233051 0.4301049The same works for coefficient alpha, and psych::bfi makes the point on real missing data: 2800 respondents answer the 25-item Big Five inventory and about 360 of them skip at least one item. We reverse-key the negatively-worded items and treat each five-item scale as a battery. The normal-theory FIML estimator handles the missing responses under ignorable missingness.
Is the Conscientiousness scale equally reliable for men and women? Different respondents, so the samples are independent:
g <- split(seq_len(nrow(bfi)), bfi$gender)
alpha_test(
men = alpha(scale_items("C")[g[["1"]], ], estimator = "nt_fiml"),
women = alpha(scale_items("C")[g[["2"]], ], estimator = "nt_fiml"),
paired = FALSE
)
#>
#> Independent-sample test of equal alpha across 2 fits
#>
#> data: men, women
#> X-squared = 0.0082484, df = 1, p-value = 0.9276
#> sample estimates:
#> men women
#> 0.7249358 0.7231688Paired: same subjects, two coefficients
When both coefficients are computed on the same subjects, the two estimates are correlated and ignoring that dependence gives the wrong standard error. paired = TRUE (the default) takes the dependence from the influence functions.
A clean illustration on real crowdsourced data: in dat.mcduff2019, MTurk judges rate smiling images on two attributes (smile presence and image positivity). Each image is rated on a 1-6 scale by only a handful of the 121 judges, and almost no two judges share images, so the judges are treated as exchangeable and each image is summarised by its category counts. Because both coefficients are estimated on the same items, the contrast is dependent.
to_counts <- function(attribute) {
m <- with(dat.mcduff2019,
tapply(get(attribute), list(item, judge), function(z) z[1]))
ratings_to_counts(matrix(as.numeric(m), nrow = nrow(m)))
}
smile <- kappa_counts(to_counts("rating_smile"), weight = "quadratic")
positive <- kappa_counts(to_counts("rating_positive"), weight = "quadratic")
kappa_test(smile = smile, positive = positive, coef = "Fleiss")
#>
#> Paired (dependent) test of equal Fleiss across 2 fits
#>
#> data: smile, positive
#> X-squared = 0.013445, df = 1, p-value = 0.9077
#> sample estimates:
#> smile positive
#> 0.6421397 0.6451303The headline is that this paired test runs on a genuinely crowdsourced design: the exchangeable count estimator supplies both the point estimates and the influence functions that encode their dependence.
Coefficient alpha has the same structure. Every bfi respondent answers all five scales, so two scales’ alphas are estimated on the same people and the comparison is paired. Are Neuroticism and Extraversion equally reliable?
alpha_test(N = alpha(scale_items("N"), estimator = "nt_fiml"),
E = alpha(scale_items("E"), estimator = "nt_fiml"),
paired = TRUE)
#>
#> Paired (dependent) test of equal alpha across 2 fits
#>
#> data: N, E
#> X-squared = 27.493, df = 1, p-value = 1.576e-07
#> sample estimates:
#> N E
#> 0.8138131 0.7615000G-way homogeneity: are several coefficients all equal?
With more than two fits the test becomes the joint hypothesis that they are all equal, with G - 1 degrees of freedom.
For coefficient alpha, all five bfi scales are answered by the same respondents, so testing that they are equally reliable is a paired G = 5 homogeneity test (df = 4):
scales <- c("A", "C", "E", "N", "O")
fits <- lapply(scales, function(L) alpha(scale_items(L), estimator = "nt_fiml"))
names(fits) <- scales
do.call(alpha_test, c(fits, list(paired = TRUE)))
#>
#> Paired (dependent) test of equal alpha across 5 fits
#>
#> data: A, C, E, N, O
#> X-squared = 274.18, df = 4, p-value < 2.2e-16
#> sample estimates:
#> A C E N O
#> 0.7021035 0.7256194 0.7615000 0.8138131 0.5999111The test decisively rejects: Neuroticism (0.81) is far more reliable than Openness (0.6).
The kappa analogue is “is the two-rater quadratic coefficient constant across all rater pairs?”. The Zapf (2016) data has four raters and hence six pairs, all on the same 50 subjects:
pairs <- combn(4, 2)
pair_fits <- apply(pairs, 2, function(rs)
kappa(dat.zapf2016[, rs], estimator = "pairwise"), simplify = FALSE)
names(pair_fits) <- apply(pairs, 2, paste0, collapse = "")
do.call(kappa_test, c(pair_fits, list(coef = "Conger", paired = TRUE)))
#>
#> Paired (dependent) test of equal Conger across 6 fits
#>
#> data: 12, 13, 14, 23, 24, 34
#> X-squared = 5.4178, df = 5, p-value = 0.367
#> sample estimates:
#> 12 13 14 23 24 34
#> 0.9309091 0.8811545 0.9012694 0.8786611 0.9403509 0.8466753Under the hood
All four cases are built only on the exported coef(), vcov(), and fit$psi. If you need an arbitrary linear hypothesis (general contrasts with non-zero margins across several fits), the internal joint_vcov() and wald_test() provide the general engine that kappa_test() and alpha_test() are convenience wrappers around.