Test Bivariate Normality for Two Ordinal Variables

Likelihood-ratio chi-square test for bivariate normality, the assumption underlying polychoric correlation, following Joreskog (2005). The observed cell frequencies of the contingency table are compared with the expected frequencies under a bivariate normal threshold model.

Usage

test_bivariate_normality(x, y = NULL, alpha = 0.05, max_levels = 10L)

Arguments

x, y

Integer or factor vectors representing ordinal (or binary) variables. Alternatively, a pre-computed contingency table (matrix or table object) can be passed as x with y = NULL.

alpha

Significance level for the test. Default: 0.05.

max_levels

Integer. Cap on the number of categories per variable when building the contingency table from numeric vectors. A numeric variable with more than max_levels distinct values (e.g. a continuous or count variable passed in for a polyserial-style screen) is binned into max_levels quantile groups before the table is formed. This keeps the maximum-likelihood threshold model tractable: an unbinned continuous variable would induce one threshold per distinct value, making the ML fit prohibitively slow and the chi-square reference distribution meaningless. Genuine ordinal/binary inputs (at most max_levels distinct values) are never altered. Default: 10.

Value

An object of class "smartcor_normtest", a list with elements:

statistic

The likelihood-ratio chi-square statistic.

df

Degrees of freedom.

p.value

The p-value from pchisq(statistic, df).

rho

The two-step estimate of the polychoric correlation.

row_thresholds

Estimated thresholds for the first variable.

col_thresholds

Estimated thresholds for the second variable.

n

Total number of observations.

m1

Number of categories for the first variable.

m2

Number of categories for the second variable.

alpha

The significance level used.

conclusion

One of "reject", "fail_to_reject", or "saturated".

Details

The bivariate normal threshold model is fitted with the two-step estimator (Olsson 1979) via polycor::polychor(): the thresholds come from the marginal proportions and the correlation is then estimated by one-dimensional maximum likelihood given those thresholds. The likelihood-ratio statistic \(G^2 = 2 \sum_{ij} O_{ij} \log(O_{ij} / E_{ij})\) compares the observed counts \(O_{ij}\) with the expected counts \(E_{ij}\) under the fitted model (cells with zero observed count are skipped). The degrees of freedom equal \(m_1 \times m_2 - m_1 - m_2\), where \(m_1\) and \(m_2\) are the numbers of categories of the two variables. For a \(2 \times 2\) table the model is saturated (df = 0) and no test is possible. The Python smartcor package computes the same statistic, so the two implementations agree.

Examples

csv = system.file("extdata", "gss_2024_casestudy.csv", package = "smartcor")
gss = read.csv(csv)

# Check the latent-normality assumption for two ordinal variables
test_bivariate_normality(gss$degree, gss$happy)
#> 
#> ── Bivariate Normality Test (LR Chi-Square) ──
#> 
#> Table dimensions: 5 x 3
#> N = 3000
#> Polychoric rho: -0.1032
#> LR chi-square: 23.9731 (df = 7, p = 0.001152)
#> Conclusion: bivariate normality rejected at alpha = 0.05. Consider Kendall's
#> tau.