📊 statistics

Bayesian Estimation of the Univariate ACE Model With Ordinal Data

This paper introduces an efficient marginal likelihood approach implemented in Stan for Bayesian estimation of the univariate ACE model with ordinal data, which overcomes frequentist limitations regarding small samples and zero-frequency cells while decoupling computational cost from sample size.

Original authors: Mark Lai, Christopher Beam

Published 2026-07-31

📖 4 min read☕ Coffee break read

CC BY 4.0

Original authors: Mark Lai, Christopher Beam

Original paper licensed under CC BY 4.0 (https://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to solve a giant puzzle to figure out what makes people tick. Are we shaped mostly by the genes we inherit from our parents, the environment we grow up in, or the unique experiences that happen just to us? Scientists have been trying to crack this code for decades using "twin studies." They compare identical twins, who share 100% of their DNA, with fraternal twins, who share about 50%. By seeing how much more alike identical twins are compared to fraternal twins, researchers can estimate how much of a trait—like height, personality, or even how often someone feels sad—is due to genetics versus environment.

However, there's a catch. In real life, we rarely get perfect numbers. We often have to ask people to choose from a list of options, like "never," "sometimes," or "always." In the world of statistics, this is called "ordinal data." It's like trying to measure the exact temperature of a room using only a thermometer that says "cold," "warm," or "hot." Traditional math tools struggle with these fuzzy categories, especially when you don't have a huge crowd of people to study or when some answer choices are rarely picked (creating "empty" spots in the data). When the math gets stuck, the answers can be wild guesses or even impossible, like saying a trait is 120% genetic. This paper tackles the problem of how to get reliable answers from these messy, fuzzy surveys without needing a massive army of participants.

The authors, Mark Lai and Christopher Beam, have built a new, super-efficient way to solve this puzzle using a method called Bayesian estimation. Think of traditional math as trying to guess the weight of a mystery box by lifting it once and hoping for the best. If the box is light or the data is messy, you might guess wrong. The old way of doing this with twins often hits a wall when the data is sparse or has empty boxes. The authors' new approach is more like having a smart, curious detective who gathers clues, updates their theory with every new piece of evidence, and keeps a running list of all the possible weights the box could be, rather than just one guess.

The paper introduces a clever shortcut. Instead of trying to simulate every single person in a study one by one (which is slow and computationally heavy, like counting every grain of sand on a beach), their method looks at the "summary" of the data—the count of how many people picked each answer. It's like looking at a map of the beach to see where the sand is piled up, rather than walking the whole beach. This makes the math incredibly fast, taking less than a second even for large groups, whereas older methods could take minutes or hours.

Crucially, this new method uses "priors," which are like giving the detective a hint based on what we already know. For example, we know that the parts of a trait (genetics, shared environment, unique environment) must add up to 100% and can't be negative numbers. The new method builds these rules directly into the math. This prevents the "impossible" answers that plague older methods. In their tests, when they used a small dataset, the old methods often crashed or gave results like "negative environment," which makes no sense. The new method, however, stayed on track, providing a full picture of the uncertainty. It didn't just say "the answer is 60%"; it said, "the answer is likely between 35% and 84%, and here is the shape of that possibility."

The authors also showed that this method works well even when the data is tricky, like when some answer categories are empty. They ran simulations with small groups of twins (as few as 50 pairs) and found that while other methods often failed to find a valid answer, their new approach consistently found a solution that respected the rules of the universe (no negative percentages). They didn't just claim it worked; they proved it by comparing their results to the "gold standard" of older methods on a larger dataset, showing that their fast, new way gave answers just as accurate as the slow, heavy ones, but without the headaches.

In short, this paper offers a faster, smarter, and more reliable way to understand how much of us is nature and how much is nurture, even when the data is messy, small, or full of gaps. It turns a math problem that often breaks into one that is robust, flexible, and ready for the real world of behavioral science.

Technical Summary: Bayesian Estimation of the Univariate ACE Model With Ordinal Data

Problem Statement
The estimation of the classic univariate ACE (Additive genetic, Common environment, Unique environment) model with ordinal data presents significant challenges for frequentist methods. Common issues include difficulties with small sample sizes, the presence of zero-frequency cells in contingency tables, and the inability to reliably quantify uncertainty for derived variance components (heritability and environmental proportions). While polychoric correlations are theoretically appropriate for ordinal data, frequentist estimation often yields inadmissible solutions (e.g., negative variance components) or relies on asymptotic approximations that fail when sampling distributions are non-normal. Furthermore, traditional Bayesian approaches for this problem, which rely on data augmentation of underlying latent continuous variables, are computationally intensive and scale poorly with sample size, rendering them impractical for large datasets.

Methodology
The authors propose an efficient marginal likelihood approach for Bayesian estimation of the ACE model using Stan. Instead of augmenting data with latent variables, this method models the sufficient statistics of the contingency table directly.

Data Structure: For a univariate ordinal phenotype with $K$ categories, the data for monozygotic (MZ) and dizygotic (DZ) twins are summarized into $K \times K$ contingency tables. Due to the exchangeability of twin pairs, the authors utilize reduced frequency vectors containing $K(K+1)/2$ cells, combining symmetric off-diagonal counts.
Likelihood: The observed cell counts are modeled using a multinomial distribution. The cell probabilities are derived from the cumulative bivariate normal distribution defined by the polychoric correlations ( $\rho_{MZ}, \rho_{DZ}$ ) and shared threshold parameters. The authors implement an approximation of the cumulative bivariate normal probability using Owen's $T$ function to facilitate computation within Stan.
Parameterization: To ensure admissible solutions (non-negative variance components summing to one), the authors introduce a simplex parameterization. Instead of estimating correlations first and deriving ACE components, the model directly estimates the variance components ( $h^2, c^2, e^2$ ) as a simplex vector. These are then transformed back into the polychoric correlations ( $\rho_{MZ} = h^2 + c^2$ and $\rho_{DZ} = 0.5h^2 + c^2$ ) within the model block.
Priors: The study employs weakly informative priors, such as Beta distributions for transformed correlations and a Dirichlet distribution for the simplex ACE components, to regularize estimates and prevent extreme values (e.g., correlations near $\pm 1$ or infinite thresholds).

Key Contributions

Computational Efficiency: The marginal likelihood approach decouples computational cost from sample size, depending instead on the number of cells in the contingency table. Benchmarking against the blavaan package (which uses data augmentation) demonstrated that the proposed method is orders of magnitude faster (e.g., ~1 second vs. hundreds of seconds for datasets with up to 1,552 twin pairs).
Handling Inadmissible Solutions: By using the simplex parameterization, the method inherently constrains variance components to be non-negative and sum to one, eliminating the need for ad-hoc constraints or post-hoc corrections common in frequentist maximum likelihood estimation.
Robustness with Sparse Data: The approach effectively handles zero-frequency cells without requiring arbitrary adjustments (such as adding 0.5 to cells), providing stable estimates where frequentist methods often fail to converge or produce boundary estimates.
Comprehensive Uncertainty Quantification: Unlike frequentist methods that rely on the delta method or profile likelihoods (which may fail in small samples), the Bayesian framework provides full posterior distributions for all parameters, including derived quantities like heritability.

Results
The paper demonstrates the method using a subset of data from Jöreskog (2021) regarding BMI categories (normal, overweight, obese).

Convergence and Fit: The MCMC sampling converged successfully with $\hat{R}$ statistics of 1.00 and high effective sample sizes. Posterior predictive checks indicated that the model fit the data reasonably well, though some discrepancies were noted in specific cells for MZ twins.
Estimation Performance: In the small sample subset, the unconstrained Bayesian model produced inadmissible estimates (e.g., negative common environment variance). The simplex parameterization corrected this, yielding plausible estimates (e.g., $h^2 \approx 0.63$ , $c^2 \approx 0.12$ ) that aligned closely with estimates derived from the full dataset.
Comparison with Frequentist Methods: In a simulation study with small sample sizes ( $N=50$ and $100$ per zygosity), the frequentist approach (umx) failed to produce confidence intervals in 100% of replications and yielded inadmissible solutions in a majority of cases. While the unconstrained Bayesian model showed lower bias in some metrics, it also produced inadmissible solutions in a significant portion of replications. The simplex parameterization traded a slight increase in bias for a dramatic increase in the rate of admissible solutions (55–65% vs. 24–35% for frequentist methods), providing more reliable estimates for applied researchers.

Significance
The paper claims that this marginal likelihood approach offers a computationally efficient and theoretically sound alternative for twin modeling with ordinal data. Its primary significance lies in making Bayesian estimation of the ACE model feasible for large datasets and robust in small-sample scenarios where frequentist methods struggle. By integrating theoretical constraints (non-negativity and summation to one) directly into the parameterization via the simplex, the method ensures that results are always admissible. The authors emphasize that this framework allows researchers to incorporate prior knowledge effectively and obtain rich uncertainty quantification, which is particularly valuable when data are sparse or distributions are skewed. The paper concludes that while the current focus is univariate, the approach is directly extensible to other univariate models (AE, CE, ADE) and potentially to multivariate models via pairwise likelihoods.

More like this