Choosing the right correlation for ordinal data: practical recommendations from a simulation study

1 Introduction

Correlations are a statistical measure that assesses the strength of the relationship between two quantitative variables [1, 2]. Commonly used across many fields, correlations are a mainstay of psychometric analysis and used heavily to evaluate instruments measuring health-related concepts such as clinical outcome assessments (COAs). Notably, COAs often are composed of multiple ordinal items that provide observed responses meant to represent the true variable that we really want to measure, such as pain or physical functioning. Such variables cannot be directly measured, which means that they are considered latent, but they can be approximated through data from the series of ordinal items. If we could measure the continuous latent variable, then the most well-known correlation coefficient, Pearson's r, would be the appropriate coefficient [2, 3]. However, ordinal data do not produce the continuous-normal data required to use Pearson's r. Although under some circumstances data transformation may be considered to satisfy requirements of statistical tests, data transformation is not always a viable option, and thus alternative correlation types must be considered.

Spearman's rank correlation coefficient, ρ, is a widely used nonparametric alternative [4] to Pearson's r. Spearman's ρ assesses how well the relationship between two variables can be described using a monotonic function, whether the relationship is linear or not, and can be applied to either continuous or to ordinal (rank) data regardless of normality [5].

Polychoric and polyserial correlations are additional measures of bivariate association designed for use when either both (polychoric) or one (polyserial) of the observed variables are ordinal measures of an underlying latent variable [6–12]. Such correlations assume that the ordinal variables are a result of partitioning the underlying, latent continuous variable(s) into discrete ordered categories. The polychoric correlation extends the tetrachoric correlation to polytomous items; the polyserial correlation is an extension of the biserial correlation [6]. Importantly, both assume that a continuous-normal latent (unobserved) variable underlies each ordinal item, and the polyserial also assumes the observed continuous item is normal, with their goal being to estimate the correlation that would be computed if the latent variable(s) had been observed.

Notably, the data structures suitable for polychoric or polyserial correlations overlap with those suitable for Spearman's correlation. Prior literature has demonstrated that the polychoric correlation performs better than Spearman when its assumptions are met [13, 14], but the polychoric correlation is biased when its assumptions are violated [15–21]. However, less work has been conducted investigating similar properties for the polyserial correlation, and few papers directly compare the Spearman correlation to either polychoric or polyserial correlation methods under varying conditions. Furthermore, two versions of polychoric and polyserial correlations have been developed, one full maximum likelihood (ML) and one computationally simpler two-step (two-step) version [6, 8, 9, 22], yet few publications compare these methods' accuracy under different conditions.

Although the normality of latent variables is commonly assumed, some literature calls this assumption into question [23–28], leading to the question of whether Spearman correlations would outperform polychoric or polyserial correlations under non-normal conditions. Spearman correlations are more widely-known and appear more commonly utilized compared to polychoric or polyserial correlations. Nevertheless, because correlations can be the basis for making critical decisions about COAs, such as the removal of items or the evaluation of validity, it is imperative to select the correlation that provides the most precise estimation for the data type.

1.1 Research questions

The primary research question was: Under what conditions is the Spearman correlation coefficient more accurate than polychoric or polyserial correlation coefficients? This overarching question was addressed via two sub-questions, corresponding to Experiments 1 and 2 respectively: (1) Under conditions for which polychoric or polyserial correlations were intended, how does the performance of the Spearman correlation compare? (2) If the normality assumption of polychoric and polyserial correlations is violated, does the Spearman correlation outperform them? An additional question was, given two-step methods are commonly used to calculate polychoric and polyserial correlations, are there conditions under which ML methods demonstrate noticeably superior performance?

2 Methods

We used simulated data to test the accuracy of different correlation types under a variety of data conditions.

2.1 Item type and skewness conditions

We examined two conditions defining the types of items that would be correlated. The first involved correlating two ordinal items that were each produced from latent continuous variables (θ; “theta” in figures), corresponding to when polychoric correlations might be used. The second involved correlating an ordinal item, produced from a latent continuous variable, with a continuous item, corresponding to when polyserial correlations might be used.

We varied characteristics of the distributions of item pairs, such that either or both could be symmetric or skewed, they could be skewed in the same or opposite directions (pairs skewed in the same direction were always positively skewed), and the skewness severity ranged from mild to high. Because we examined both ordinal-ordinal and ordinal-continuous item pairs, at times we use the term “symmetric” for simplicity to encompass both ordinal and continuous items.

2.1.1 Experiment 1

Experiment 1 focuses on conditions that do not violate assumptions of polychoric or polyserial correlations (extending Coon [29] to polyserial correlations). When correlating two ordinal items (ordinal-ordinal; corresponding to polychoric correlations), the distribution conditions are: (1) both observed ordinal items symmetric, (2) one observed ordinal item symmetric and one skewed (i.e., non-symmetric), (3) both observed ordinal items skewed in the same direction, and (4) both observed ordinal items skewed in opposite directions. When correlating one ordinal to one continuous item (ordinal-continuous; corresponding to polyserial correlations), the conditions were: (1) both observed items symmetric (i.e., symmetric ordinal item and normally-distributed continuous item), and (2) the observed ordinal item skewed. A skewed observed continuous item would violate an assumption of the polyserial correlation, so this condition is considered in Experiment 2.

In all cases, the latent variables and the observed continuous variable were normal, and skewness in the ordinal items was produced by shifting the item response thresholds; that is, in Experiment 1, ordinal items with skewed (non-symmetric) distributions reflected changes to the process of partitioning a normal latent variable into an ordinal item. This is important to note given the assumptions of the correlations in question. That is, polychoric correlations specifically assume normality of the latent variables underlying the observed ordinal items, but do not assume symmetry in the observed ordinal items; similarly, polyserial correlations assume normality of the observed continuous item and of the latent variable underlying the observed ordinal item, but do not assume symmetry in the observed ordinal item.

2.1.2 Experiment 2

Experiment 2 focuses on conditions that violate assumptions of polychoric or polyserial correlations (extending Coon [29] to when assumptions are violated). For ordinal-ordinal correlations, the distribution type conditions are: (1) one latent continuous variable skewed, (2) both latent continuous variables skewed in the same direction, and (3) both latent continuous variables skewed in opposite directions. For ordinal-continuous correlations, the conditions are: (1) observed continuous item skewed, (2) latent variable (underlying the observed ordinal item) skewed, (3) both the observed continuous item and the observed ordinal item skewed in the same direction, (4) both the latent variable (underlying the observed ordinal item) and the observed continuous item skewed in the same direction, (5) both the observed continuous item and the observed ordinal item skewed in opposite directions, and (6) both the latent variable (underlying the observed ordinal item) and the observed continuous item skewed in opposite directions. In all cases, the ordinal items used symmetric response thresholds.

In summary, for all conditions of Experiment 2, at least one latent variable (or, for ordinal-continuous item pairs, at least one of either the latent variable or the observed continuous item) was non-normal, and therefore violated an assumption of polychoric or polyserial correlations.

2.2 Common conditions tested

For both experiments and for all item types and distributions, we tested the same conditions related to the strength of the population (or “true”) correlation and the number of response categories used by the ordinal item, using a fully-crossed design for these two variables. Population correlation magnitudes were r = 0.3, r = 0.5, and r = 0.8, and there were 3, 5, 7, 9, or 11 response categories. Additionally, when any item or variable was skewed, the skew severity was varied (mild, medium, moderately-high, and high skewness) within each of the distribution conditions.

2.3 Simulation methods

Each condition had 500 replications conducted. For all conditions, simulations began by producing two continuous variables drawn from a normal population with a mean of 0 and standard deviation (SD) of 1 (N[0,1]) for a sample size of 500. The relationship between the two variables was governed by the strength of their intended correlation coefficient.

For skewed continuous variables, the normal distribution was transformed into a gamma distribution with a population SD of 1 (see Section 2.5 for further details). Skew severity was varied by changing the gamma shape parameter (α), specifically using the values 8 (mild), 5 (medium), 2 (moderately-high), and 1 (high); the gamma rate parameter (λ) was set using λ=α to maintain the population SD of 1 (based on the property of the gamma distribution, σ2= αλ2 where σ is the population standard deviation [30]). Based on the expected properties of the normal or gamma distributions, in the latter case using the formulae Skewness= 2α and Kurtosis= 6α ([30]; see Supplementary Methods for additional properties of the gamma distribution and its relationship to the normal distribution), the expected skewness for each distribution condition is 0 for the normal or symmetric conditions, 0.71 for mildly-skewed, 0.89 for medium skew, 1.41 for moderately-high skew, and two for high skew, and the kurtosis is 0 for the normal or symmetric conditions, 0.75 for mildly-skewed, 1.2 for medium skew, three for moderately-high skew, and six for high skew. Since the gamma distribution can only be positively skewed, producing a negatively skewed latent variable required multiplying the resulting variable by −1, in which case the true population correlation value also became negative; for simplicity, figures are shown relative to the absolute value of the population correlation value, i.e., the population correlation magnitude.

In order to ensure comparability between the conditions using normal and skewed continuous variables, a version of the simulation was run with the normal variables also gamma-transformed, taking advantage of the gamma distribution's property that it converges to the normal distribution as α asymptotically increases (e.g., Leemis and McQueston [31]; specifically, the parameters α = 999, λ = 31.61, and σ = 1 were used). Results of these comparisons are shown in the Supplementary Methods.

To produce an ordinal item, the values of the continuous variable were polychotomized (classified into a series of discrete, ordered values) using response thresholds. Item response thresholds were systematically chosen based on the properties of the corresponding population distribution. Symmetric item response thresholds were chosen by evenly splitting the central 2.5 σ's of N(0,1) based on the number of response categories. For the conditions where ordinal items were skewed with a normal underlying distribution (Experiment 1), skewed thresholds were created by (1) computing evenly split percentiles based on k+2 where k is the number of response categories (starting at 0.0001 and ending at 0.9999 to avoid infinite values), (2) calculating the values corresponding to those percentiles on the relevant gamma distribution using the gamma quantile function and the chosen α and λ parameters, (3) computing the differences between sequential values in step 2, (4) computing the corresponding values for the normal distribution by starting at μ−2.5σ (for α ≤ 2) or μ−3σ (for α ≥ 5) and incrementally adding each of the differences computed in step 3, and using these values as the response thresholds. The thresholds used for each condition are shown in Table 1.

Table 1. Item response thresholds by condition.

2.4 Analysis

2.4.1 Correlation types

We calculated three correlation types between the two items for each replication; Spearman was always calculated, along with polychoric two-step and ML for ordinal-ordinal conditions, or polyserial two-step and ML for ordinal-continuous conditions.

The process for calculating a Pearson product moment correlation is well-known and is simply the covariance of two variables divided by the product of the standard deviation of each variable. This calculation assumes a linear relationship between the variables, so the Spearman correlation circumvents that assumption by focusing on the relationship being monotonic. The Spearman correlation uses the same equation as the Pearson except that instead of the covariance and standard deviations of the two variables, it relies on the covariance and the standard deviations of the rank ordering of the data within the two variables. The full equations for these correlations can be found in Supplementary Methods.

Both Pearson and Spearman correlations are closed-form solutions, whereas the polychoric and polyserial correlations require more complex iterative algorithms. For 2 observed ordinal variables, the polychoric correlation is defined as the product-moment correlation (ρ) between the corresponding underlying (i.e., latent, or unobserved) normal variables. For 1 observed ordinal variable and 1 observed continuous variable X, the polyserial correlation is defined as the product-moment correlation between the observed continuous variable and the normal latent variable that underlies the ordinal variable. More information about the mathematical equations that underlie the polychoric and polyserial correlations can be found in Drasgow [6]; the mathematical process and relevant equations are also summarized in Supplementary Methods. Importantly, in addition to the calculation of the correlations themselves, thresholds that represent the polychotomization process for going from a continuous latent variable to an ordinal observed variable, must also be estimated. Polychoric and polyserial correlations can either be estimated using joint maximum likelihood, where the correlations and the thresholds between each of the ordinal variable levels are estimated iteratively, or using a two-step procedure, where the thresholds are estimated in step one and then the correlation is estimated iteratively [32, 33].

It is notable that for some circumstances, data may be transformed to resemble a normal distribution, thus satisfying requirements of many statistical analyses, and there is robust literature on various methods for such transformations (e.g., Atkinson and Cox [34]). However, this paper is primarily intended for circumstances when data transformation is expected to be avoided where possible, such as when interacting with regulatory agencies [35], or when data transformation may not be appropriate, such as may be the case for some ordinal variables [36, 37]. At reviewer request, comparison of Spearman and polychoric or polyserial correlations to correlation after data transformation have been included in the Supplementary Results.

2.4.2 Bias and absolute error

For each condition, the percent bias and absolute percent error were calculated comparing the observed correlation value to the population correlation magnitude, according to the formulae given below.

Percent Bias= (Observed Correlation-Population Correlation)Population Correlation×100

Absolute Percent Error=|(Observed Correlation-Population Correlation)Population Correlation|×100

Primary results figures show percent bias for each condition. Because some correlation methods were observed to produce consistently positive bias while others were observed to produce consistently negative bias, to supplement these result figures and determine which methods produced the least overall magnitude of error, statistical models were fit to the absolute percent error.

As would be expected, the absolute percent error distribution was positively skewed, therefore these results were analyzed using a gamma regression model with an inverse link, and joint hypothesis tests were used to produce a single main effect or interaction result (non-sequential). Predictors included as factors were Correlation Method and the Number and Direction of Skewed Items, while predictors included as numeric were Skewness Severity (quantified as the skewness of the resulting distribution if continuous), Number of Response Categories, and Population Correlation magnitude. Because the fully saturated model would be quite large with higher-order interactions that are difficult to interpret, only up to three-way interactions were included in the models. Primary results figures show percent bias for each condition; figures showing absolute percent error results are presented in Supplementary Results.

2.5 Technical specifications

All simulations and analyses were conducted using R software version 4.5.0 [38]. Correlated normal continuous variables were simulated and transformed to gamma distributions using the faux package, specifically the rnorm_multi and norm2gamma functions, respectively [39]; the norm2gamma function relies on the normal distribution function (pnorm R function [38]) and the gamma quantile function (qgamma R function [38]; full details and code can be found at DeBruine [40]). When using the faux package, the argument “empirical” was set to false, indicating that the provided correlation value, mean, and SD should be treated as population values not empirical values.

Spearman correlations were calculated using the stats package and the cor.test function [38] and polychoric and polyserial correlations were calculated using the polycor package, the former with the polychor function and the latter with the polyserial function, using both the two-step and full ML methods [8, 41]. Observed skewness and kurtosis values were calculated using the moments package [42].

Statistical models were run using the stats package [38], joint hypothesis tests were conducted using the marginaleffects package [43], and post-hoc marginal means, contrasts, and trends were computed using the emmeans package [44], in R. Results are reported in the original metric.

R code used for this study can be found at https://github.com/ndayasse/Polychoric_Polyserial_Sims_2025.

3 Results

All presented figures depict descriptive statistics. Bar plots are used to indicate the mean proportion of the simulated samples endorsing each possible response option (i.e., histograms with proportions as the y-axis). Here these are shown only for conditions with five response categories; relative frequency histograms for all conditions can be found in the Supplementary Results. Other result figures show the mean percent bias (points; error bars indicating the 95% bootstrapped confidence interval were originally included but are small enough to be difficult to plot visibly); other than under fully symmetric conditions, these are shown here for population correlations of 0.3 and 0.8, with the population correlation 0.5 condition shown in Supplementary Results. Plots depicting the observed skewness and kurtosis of the simulated data can be found in the Supplementary Results.

3.1 Experiment 1 results

Figure 1 shows histograms of the ordinal-ordinal data structures in Experiment 1, specifically for the conditions with five response categories (relative frequency histograms for all conditions can be found in the Supplementary Results), and Figures 2, 3 show the percent bias results of observed compared to population correlation magnitudes; Figures 4–6 show the corresponding information for ordinal-continuous data.

3.1.1 Ordinal-ordinal data

For ordinal-ordinal data, looking first at the symmetric items, almost all examined variables have significant effects; almost all tested main effects and interactions are statistically significant, with the only non-significant term being the two-way interaction between Number of Response Categories and Population Correlation magnitude (Table 2). Importantly, the Spearman correlations result in overall greater error than either of the polychoric methods, and the difference between polychoric methods is negligible. Spearman correlations also demonstrate a greater effect on estimation error from increasing the Number of Response Categories than either of the polychoric methods (supported by marginal trends; Spearman: −1.19 [−1.23, −1.16]; polychoric two-Step: −0.17 [−0.20, −0.14]; polychoric ML: −0.17 [−0.19, −0.14]).

Table 2. Ordinal-ordinal absolute percent error results when all items are symmetric (Experiment 1).

Estimates are shown for numeric predictors and for specific marginal means and contrasts.

CI, confidence interval; DF, degrees of freedom; Correl., correlation; No. Resp. Cat., number of response categories; Pop. Correl. Mag., population correlation magnitude.

For the remainder of the conditions with normally-distributed latent variables, all examined variables have significant effects on estimation error (all tested main effects and interactions are statistically significant; Table 3). Again, the Spearman correlations result in greater error overall than either of the polychoric and the difference between polychoric methods is negligible. Spearman correlations also demonstrate a larger effect of Skew Severity and of the Number of Response Categories on estimation error compared to either of the polychoric correlations (see Table 4 estimated marginal trends).

Table 3. Ordinal-ordinal absolute percent error results when items are not symmetric but have normally distributed latent variables (Experiment 1).

Estimates are shown for numeric predictors and for specific marginal means and contrasts.

CI, confidence interval; DF, degrees of freedom; Correl., correlation; No., number; Vars., variables; No. Resp. Cat., number of response categories; Pop. Correl. Mag., population correlation magnitude.

Table 4. Ordinal-ordinal absolute percent error marginal trends when items are not symmetric but have normally distributed latent variables (Experiment 1).

3.1.2 Ordinal-continuous data

For ordinal-continuous data, all examined variables have significant effects on estimation error (all tested main effects and interactions are statistically significant; Tables 5, 6). Similarly to the ordinal-ordinal data, the Spearman correlations result in overall greater error than either of the polyserial methods. Spearman correlations also show a larger effect of the Number of Response Categories on estimation error compared to either polyserial method (symmetric items' marginal trends: Spearman −0.80 [−0.83, −0.76]; polyserial two-Step −0.08 [−0.11, −0.06]; polyserial ML −0.07 [−0.10, −0.05]; remaining conditions' marginal trends: Spearman −1.49 [−1.52, −1.46]; polyserial two-Step −0.14 [−0.16, −0.13]; polyserial ML −0.14 [−0.16, −0.13]). The differences between polyserial correlation methods are consistently small; for symmetric items, the difference does not reach statistical significance, and for not-symmetric items, the two-step method demonstrates greater absolute percent error than the ML, but by less than one percentage point.

Table 5. Ordinal-continuous absolute percent error results when all items are symmetric (Experiment 1).

Estimates are shown for numeric predictors and for specific marginal means and contrasts.

CI, confidence interval; DF, degrees of freedom; Correl., correlation; No. Resp. Cat., number of response categories; Pop. Correl. Mag., population correlation magnitude.

Table 6. Ordinal-continuous absolute percent error results when ordinal item is not symmetric but latent variable and continuous item are normally distributed (Experiment 1).

Estimates are shown for numeric predictors and for specific marginal means and contrasts.

CI, confidence interval; DF, degrees of freedom; Correl., correlation; No. Resp. Cat., number of response categories; Pop. Correl. Mag., population correlation magnitude.

3.2 Experiment 2 results

Figure 7 shows histograms of the ordinal-ordinal data structures in Experiment 2, specifically for the conditions with five response categories (relative frequency histograms for all conditions can be found in the Supplementary Results), and Figure 8 shows the percent bias results of observed compared to population correlation magnitudes; Figures 9 and 10 show the corresponding information for ordinal-continuous data.

3.2.1 Ordinal-ordinal data

For the ordinal-ordinal data, almost all examined variables have significant effects on estimation error; all tested main effects and interactions are statistically significant except for the two-way interactions between the Number and Direction of Skewed Items and Population Correlation magnitude, and between Skew Severity and Population Correlation magnitude (Table 7). Again, the Spearman correlations result in greater error overall than either polychoric method, and the difference between polychoric methods is negligible. Spearman correlations also show larger effects on estimation error from both Skew Severity and the Number of Response Categories (see Table 8 estimated marginal trends) than do the polychoric methods.

Table 7. Ordinal-ordinal absolute percent error results when latent variables are not all normally distributed (Experiment 2).

Estimates are shown for numeric predictors and for specific marginal means and contrasts.

CI, confidence interval; DF, degrees of freedom; Correl., correlation; No., number; Vars., variables; Dir., direction; No. Resp. Cat., number of response categories; Pop. Correl. Mag., population correlation magnitude.

Table 8. Ordinal-ordinal absolute percent error marginal trends when latent variables are not all normally distributed (Experiment 2).

CI, confidence interval; DF, degrees of freedom.

3.2.2 Ordinal-continuous data

For the ordinal-continuous data, almost all examined variables have significant effects on estimation error; all tested main effects and interactions are statistically significant except for the main effect of Skew Severity (Table 9; see also Table 10 estimated marginal trends). As with prior conditions, the Spearman correlations result in greater error overall than the polyserial methods. Additionally, the Number of Response Categories consistently, and Skew Severity sometimes, produce a larger effect on the estimation error of Spearman correlations compared to the polyserial methods. (Specifically, for most levels of the Number and Direction of Skewed Items, Spearman correlations are more affected by Skew Severity. When both observed items are skewed, in either direction, all three correlation methods are similarly affected by Skew Severity. When only the continuous item is skewed, the polyserial methods are more affected by Skew Severity.) Although the difference between the polyserial correlation methods is less than one percentage point, the two-step method demonstrates slightly greater absolute error than the ML; this pattern is similar to that observed for not-symmetric ordinal-continuous item pairs with normal latent variables in Experiment 1.

Table 9. Ordinal-continuous absolute percent error results when observed continuous item and latent variable are not both normally distributed (Experiment 2).

Estimates are shown for numeric predictors and for specific marginal means and contrasts.

CI, confidence interval; DF, degrees of freedom; Correl., correlation; No., number; Vars., variables; Dir., direction; No. Resp. Cat., number of response categories; Pop. Correl. Mag., population correlation magnitude.

Table 10. Ordinal-continuous absolute percent error marginal trends when observed continuous item and latent variable are not both normally distributed (Experiment 2).

CI, confidence interval; DF, degrees of freedom.

4 Discussion

4.1 Summary of results

When all latent variables are normal (i.e., the normality assumption was not violated; Experiment 1), for ordinal-ordinal conditions, the polychoric correlation coefficient consistently outperforms Spearman correlations. Under these conditions, Spearman correlation coefficients tend to underestimate the population correlation magnitude; in contrast, when polychoric correlations are biased, they tend to overestimate the population correlation magnitude. The results demonstrate that the advantage of polychoric correlations is greater when there are fewer response categories, higher population correlation magnitudes, and/or greater skew severity. The two-step and full ML polychoric methods demonstrate negligible differences.

For ordinal-continuous conditions with normality assumptions met, the polyserial correlation coefficient demonstrates a similar advantage over Spearman as the polychoric, and Spearman again tends to underestimate the population correlation magnitude. As with ordinal-ordinal data, the advantage of polyserial correlations is greater with fewer response categories, higher population correlation magnitudes, and/or greater skew severity. Although for symmetric items the two-step and ML polyserial methods demonstrate negligible differences, when skewness resulted from the ordinal item's thresholds, the polyserial ML method showed overall slightly lower error than the two-step; this difference was slightly more pronounced at higher population correlation magnitudes and with more response categories.

When the underlying latent or observed continuous variables are skewed, violating assumptions of polychoric and polyserial correlations (Experiment 2), the pattern of results remains similar for ordinal-ordinal correlations but changes somewhat for ordinal-continuous correlations. Polychoric correlations maintain a clear advantage over Spearman with negligible differences observed between polychoric ML and two-step methods.

For ordinal-continuous data that violate normality assumptions, polyserial correlations tend to generally maintain their advantage over Spearman correlations. The degree of bias was similar between Spearman and polyserial correlations when there were fewer response categories, only the continuous item was skewed, and that item had high skew severity. However, Spearman correlations demonstrated less absolute error than either polyserial method when there were more response categories, a low population correlation magnitude, and either only the continuous item was skewed or skew severity was low. Under other conditions, Spearman correlations were less accurate than polyserial correlations. Notably, when Spearman correlations had lower absolute error than polyserial, the difference never exceeded 5 percentage points, however, when polyserial correlations had lower absolute error than Spearman, the difference commonly exceeded 5 or even 10 percentage points; under some conditions, this difference exceeded 50 percentage points (specifically, when there were fewer response categories, the latent variable and continuous item were highly skewed, regardless of direction, and the population correlation magnitude was higher).

Comparing the polyserial ML and two-step methods when their assumptions are violated, the ML method is more likely to underestimate, and the two-step method to overestimate, the population correlation; the exception is when only the continuous item is skewed and all three correlations tend to underestimate the population correlation. Generally, the polyserial ML produces estimates with less error than the two-step, with some exceptions. The polyserial two-step results in similar or lower error than the ML with either (a) higher skew severity and lower population correlation magnitudes, or (b) fewer response categories, higher skew severity, and higher population correlation magnitudes. Of note, the difference in the absolute percent error of the two polyserial correlation methods never exceeded 5 percentage points.

Generally, the Spearman correlation estimations demonstrate a consistent decrease in absolute error with more response categories, while the estimation error of polychoric and polyserial correlations are overall less affected by the number of response categories. The biggest impact of the number of response categories on the estimation error of polychoric and polyserial correlations is observed for three response categories, mainly when items are skewed in opposite directions. As would be expected, all correlation methods tend to produce greater estimation error with increasing skew severity, though this effect is more pronounced for Spearman correlations and particularly for ordinal-ordinal data as opposed to ordinal-continuous.

Across conditions, when Spearman correlations demonstrate biased estimates, they tend to underestimate the population correlation magnitude. Although polychoric and polyserial correlations have generally less estimation error under the conditions investigated here, when they do demonstrate biased estimates, polychoric (ML and two-step) and polyserial two-step correlations tend to overestimate, and polyserial ML correlations usually underestimate, the population correlation magnitude.

4.2 Return to research questions

In answer to our first research question, the Spearman correlation demonstrates poorer performance than the polychoric or polyserial correlations under the conditions for which the latter two were designed (i.e., underlying normal distributions for ordinal items or observed normal distributions for continuous variables).

In answer to our second research question, for most conditions that violate polychoric or polyserial correlation assumptions, the polychoric or polyserial correlation still outperforms the Spearman correlation. While Spearman occasionally outperforms one or both of the polyserial correlations, the absolute percent error of the polyserial estimates is never more than 5 percentage points greater than Spearman estimates; in contrast, under some conditions, the absolute percent error of the Spearman estimates is as much as 50 percentage points greater than that of the polyserial estimates.

When assumptions of polychoric and polyserial correlations are met, performance of the full ML and two-step methods for these correlation types are highly similar. When their normality assumptions are violated, the performance of the polychoric two-step remains highly similar to its ML counterpart, but the performance of the polyserial two-step sometimes differs from the polyserial ML. Where their performance differed, the polyserial ML correlation usually, but not always, demonstrated better accuracy than the polyserial two-step.

4.3 Recommendations

Table 11 provides a summary of recommendations based on the simulation results presented here.

Table 11. Summary of correlation method recommendations.

We note that although some literature has endeavored to create a generalized version of the polychoric correlation method that is robust to non-normality [26, 45–48], testing of these methods has revealed mixed results [45, 48]. Additionally, estimation of a generalized polychoric correlation does not appear to be available in standard statistical software packages, and a generalized version of the polyserial correlation does not appear to have been developed. It is therefore difficult to suggest such a method for common statistical or psychometric application.

4.4 Limitations

As with all simulations, our tested conditions were not exhaustive and did not consider every possible violation of the assumptions of polychoric and polyserial correlations. In lieu of testing alternative distributions [16, 20, 21], we tested varying degrees of skewness using the gamma distribution, chosen due to the flexibility of this distribution and its ability to produce different degrees of skewness. Additionally, although every effort was made to produce polychotomous distributions with similar skewness to the corresponding continuous gamma distributions, the nature of ordinal items makes such a task difficult, and the observed skewness of the simulated items' distributions demonstrated more variability and overall trended toward slightly higher observed skewness for the ordinal items than for the continuous variables (see and ). Thus, it is possible that some of the differences observed between Experiments 1 and 2 in the correlations' estimation error could stem from these differences in observed skewness.

5 Conclusions

Without question, when conditions and assumptions are met for either polychoric or polyserial correlations, their accuracy is superior to Spearman correlations, and differences between two-step and ML methods for polychoric or polyserial correlations are negligible. When assumptions are violated by skewness of the latent variables or of the observed continuous item (the latter in the case of ordinal-continuous data), polychoric correlations maintain their superior accuracy compared to Spearman, and polyserial correlations are usually more accurate than Spearman or else relatively similar in the degree of error. Across the conditions tested here, polychoric correlations are always recommended over Spearman, and polyserial correlations are almost always recommended with limited exceptions. We strongly recommend that COA developers relying on correlational methods for evaluating their instrument's measurement properties carefully choose the most appropriate correlation type based on what is known and observed about the data, and what assumptions are testable, to ensure that proper and defensible conclusions are drawn about their COA items and scores.

Generative AI statement

The author(s) declared that Generative AI was not used in the creation of this manuscript.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fams.2026.1823819/full#supplementary-material

Supplementary material Supplementary material at the source

Footnotes

Edited by: Miguel Fernández, University of Valladolid, Spain

Reviewed by: Fuxia Cheng, Illinois State University, United States

Daniel Brooks Wright, University of Nevada, Las Vegas, NV, United States

Footnotes

2-step, two-step; θ, latent or unobserved variable; COA, clinical outcome assessment; IQR, interquartile range; ML, maximum likelihood; N[0,1], normal distribution with population mean of 0 and population standard deviation of 1; SD, standard deviation.