Source-linked AI summary

Polygenic Modeling with Bayesian Sparse Linear Mixed Models

Xiang Zhou, Peter Carbonetto, Matthew Stephens

arXiv:1209.1341v2q-bio.QMq-bio.GNstat.APstat.ME

TL;DR

Polygenic modeling must make assumptions about genetic effects when markers greatly outnumber individuals, but the appropriate assumptions are often uncertain. This paper develops BSLMM, a data-adaptive hybrid of mixed models and sparse regression that matches or exceeds LMM accuracy for PVE estimation and outperforms both in phenotype prediction.

  • Problem

    Because genome-wide datasets contain far more markers than individuals, polygenic modeling requires assumptions about genetic effects whose accuracy is often uncertain.

  • Method

    BSLMM combines large individual genetic effects with many small effects while inferring their balance from the data.

  • Results

    BSLMM matches or exceeds LMM accuracy for PVE estimation across genetic architectures and outperforms both models and related models for phenotype prediction.

  • Takeaways & Limitations

    BSLMM provides a flexible, computationally tractable approach for polygenic modeling in GWAS datasets.

  • Takeaways & Limitations

    BSLMM effect estimates may require separate consideration when individual genetic effects, rather than relationships between phenotypes and genotypes, are the primary interest.

Abstract

from arXiv · show

Both linear mixed models (LMMs) and sparse regression models are widely used in genetics applications, including, recently, polygenic modeling in genome-wide association studies. These two approaches make very different assumptions, so are expected to perform well in different situations. However, in practice, for a given data set one typically does not know which assumptions will be more accurate. Motivated by this, we consider a hybrid of the two, which we refer to as a "Bayesian sparse linear mixed model" (BSLMM) that includes both these models as special cases. We address several key computational and statistical issues that arise when applying BSLMM, including appropriate prior specification for the hyper-parameters, and a novel Markov chain Monte Carlo algorithm for posterior inference. We apply BSLMM and compare it with other methods for two polygenic modeling applications: estimating the proportion of variance in phenotypes explained (PVE) by available genotypes, and phenotype (or breeding value) prediction. For PVE estimation, we demonstrate that BSLMM combines the advantages of both standard LMMs and sparse regression modeling. For phenotype prediction it considerably outperforms either of the other two methods, as well as several other large-scale regression methods previously suggested for this problem. Software implementing our method is freely available from http://stephenslab.uchicago.edu/software.html

Author Summary

Polygenic modeling seeks to clarify how genetic variation relates to variation in observed characteristics, including quantitative traits and disease susceptibility. Better models could ultimately inform clinical practice and breeding strategies.

  • Polygenic modeling aims to understand relationships between genetic variation and variation in observed characteristics.
  • Applications include quantitative traits such as human cholesterol levels and cattle milk production, as well as disease susceptibility.
  • Improved polygenic modeling could ultimately support changes in human clinical practice and better agricultural breeding or mating strategies.

Introduction

The paper introduces BSLMM as a hybrid of linear mixed models and sparse regression for polygenic modeling. It evaluates BSLMM for PVE estimation and phenotype prediction, finding robustness for PVE and superior prediction performance.

  • Motivation: LMMs and sparse regression make nearly opposite polygenic-modeling assumptions: LMMs treat every variant as causal with normally distributed effects, whereas sparse models assume relatively few causal variants.
  • Motivation: The proposed hybrid builds on related models previously used for breeding-value prediction, gene selection, phenotype prediction, and complex-trait mapping.
  • Motivation: The study focuses BSLMM on two genome-wide association applications: estimating PVE, also called chip heritability, and predicting phenotypes.
  • Motivation: BSLMM combines the advantages of LMMs and sparse regression, remains robust across PVE-estimation settings, and outperforms both models and related models in phenotype prediction.These findings come from simulations and applications to real data.

Methods

The methods develop BSLMM as a flexible hybrid of mixed-model and sparse-regression assumptions, then address prior specification and scalable posterior computation. The model separates genome-wide random effects from sparsity-inducing effects while retaining LMM and BVSR as special cases.

  • Motivation: The high-dimensional setting, with p much larger than n, requires modeling assumptions because genetic effects and PVE cannot otherwise be estimated accurately.Existing polygenic approaches arise from different assumptions about the genetic effects β, motivating a model that can accommodate uncertainty about those assumptions.
  • Effect-size assumptions: The two-normal mixture prior is more flexible than LMM and BVSR because setting π = 0 yields LMM, while σb = 0 yields BVSR.The model treats all variants as having small normally distributed effects, with proportion π receiving an additional normally distributed effect.
  • BSLMM model: BSLMM combines a standard linear mixed model with sparsity-inducing priors on regression coefficients, separating random effects from sparse effects.The random effects have known covariance matrix K, while the sparse effects emphasize variable-specific contributions through their prior.
  • Model interpretation: Although the mixture-prior linear model is equivalent to BSLMM for a genotype-derived K, the random effect u and coefficient vector β represent different parameterizations.The random effect captures combined small marker effects in BSLMM, whereas those effects are included directly in β in the simple linear model.
  • Prior specification: A uniform prior on log π allows the expected nonzero-marker proportion to range from 1/p to 1 and reflects uncertainty spanning orders of magnitude.The alternative uniform prior on π places greater weight on large numbers of nonzero sparse effects.
  • Computation: The new likelihood-evaluation algorithm replaces per-iteration cubic matrix inversion and determinant calculations with linear computational burden in n after an initial eigen-decomposition.This scaling is described as matching BVSR and enabling analysis of large studies.

Results

Across simulations and real data, BSLMM provided stable PVE estimation and generally matched or outperformed competing methods for phenotype and disease prediction. Its hybrid modeling also offered computational advantages in highly polygenic settings, though disease-prediction accuracy was not yet clinically practical.

  • PVE estimation: BSLMM showed stable PVE-estimation performance across genetic architectures, whereas BVSR excelled only for sparse models and performed poorly for more polygenic models.BVSR’s poor performance under highly polygenic models was attributed to strong downward bias in its PVE estimates.
  • Phenotype prediction: For phenotype prediction, BSLMM outperformed LMM in nearly all simulated settings and exceeded BVSR in more polygenic scenarios.BSLMM matched BVSR when causal-SNP counts were smaller, up to S = 1,000 in Scenario I, while its random-effect component captured polygenic signal.
  • PVE estimation: For five traits, BVSR produced consistently smaller PVE estimates, while LMM and BSLMM estimates were almost identical for two traits and similar for the others.Both LMM and BSLMM estimates explained over 50% of narrow-sense heritability for the five traits.
  • Disease prediction: BSLMM performed as well as or better than LMM and BVSR for disease prediction across all seven diseases.Its performance also appeared favorable relative to previously applied methods, although precise comparisons were difficult because validation strategies differed.
  • Trait prediction: For CD8 and MCH, BSLMM consistently outperformed all other methods, whereas all methods performed equally poorly for low-heritability BMI.Prediction accuracy was measured by RMSE against BSLMM, with similar results when measured by correlation.
  • Computational performance: In simulated data with S = 10,000 and PVE = 0.6, BSLMM could be over an order of magnitude faster than BVSR, while matching BVSR’s time for few causal SNPs.The speed advantage is consistent with BSLMM behaving similarly to LMM for highly polygenic traits.

Discussion

BSLMM combines sparse large effects with polygenic small effects while remaining tractable for moderately large datasets. The discussion emphasizes flexible, data-adaptive modeling and acknowledges substantial computational costs and uncertainty about whether richer effect-size distributions would improve performance.

  • Contributions: BSLMM jointly models a few large genetic effects and many small effects, inferring their balance from the data while remaining computationally tractable for moderately large datasets.The implementation can handle at least 10,000 individuals with 500,000 genetic markers.
  • Extensions: BSLMM can be adapted to binary phenotypes by treating case-control status as quantitative for PVE estimation and correcting observed-scale PVE to the latent-liability scale.The paper supplies an alternative derivation of the correction factor in Text S3.
  • Limitations: Despite computational innovations, posterior sampling remains demanding because of memory requirements and the large number of iterations needed for convergence.Additional computational resources may increase the feasible dataset size, but the burden remains substantial.
  • Modeling principles: The discussion recommends estimating hyper-parameters from data in an integrated procedure and specifying their priors through variance proportions such as PVE and PGE.Integrated estimation outperformed the two-step BSLMM-EB procedure for phenotype prediction.
  • Limitations: It remains unknown whether BSLMM’s two-normal mixture adequately captures actual effect-size distributions or whether richer mixtures would meaningfully improve performance.The authors note that more flexible distributions could improve modeling but make the accuracy–computational burden tradeoff difficult to predict.
  • Modeling principles: More flexible effect-size distributions generally improve phenotype prediction, while BSLMM balances flexibility, performance, and computational tractability relative to available alternatives.The discussion contrasts BSLMM’s two-normal mixture with the single-normal LMM and point-normal BVSR distributions.

Text S1 Detailed Methods … Details of BSLMM

The methods use four GWAS datasets for PVE estimation and phenotype prediction, simulation scenarios with controlled causal architectures, and prediction metrics that account for genotype correlations. BSLMM analyses center genotypes and relatedness and specify hyper-parameter expressions and priors under these assumptions.

  • GWAS Datasets: Four GWAS datasets support simulation, PVE estimation, and phenotype prediction across human and mouse populations.The datasets include Australian height, blood lipids, WTCCC disease phenotypes, and three mouse phenotypes.
  • Simulations: Simulations use human height genotypes and 20 phenotype replicates per scenario, varying causal-SNP counts, effect-size distributions, and true PVE values.Scenario I uses normal causal effects with 10, 100, 1000, or 10000 causal SNPs and PVE 0.2 or 0.6; Scenario II uses a two-normal mixture.
  • MSPE and RPG: Prediction accuracy is assessed primarily with MSPE and relative predictive gain, with MSPE extended to account for correlations among neighboring SNPs.Neighboring-SNP covariance is approximated from sample covariance for markers satisfying |i −j| ≤20 because linkage disequilibrium decays with distance.
  • MSPE and RPG: For LMM and BSLMM predictions, estimates of β rather than ˜β incorporate both sparse effects β and random effects u.The baseline MSPE uses the phenotype mean, denoted MSPE0.
  • Correlation: Correlation between predicted and true future observations provides a second prediction-accuracy measure, also using β estimates for LMM and BSLMM.The correlation expression depends on ˆβ, β, and τ.
  • Centering X and K: Genotypes are coded as 0, 1, or 2 reference-allele copies and centered column-wise; imputed genotypes use posterior mean values.The genotype matrix is denoted X, with centered entries Xij = (xij −¯xj).
  • Centering X and K: The relatedness matrix is centered so random effects have mean zero, either automatically for K ∝XXT with centered X or through projection by MKM.The transformed matrix remains positive-semidefinite when the original relatedness matrix is positive-semidefinite; the h and ρ derivations assume centered genotypes and relatedness.
  • Definition and derivation of expressions for h and ρ: The h and ρ expressions approximate expected PVE and PGE by ratios of expectations conditional on BSLMM hyper-parameters, with independent priors on h, ρ, and π.The induced priors have polynomial tails and make marker-effect variance σ2a decrease as the proportion of affected markers π increases.

PVE Estimation with LMM, BVSR and BSLMM · Phenotype Prediction with LMM, BVSR and BSLMM

The section describes how PVE is estimated under LMM, BVSR, and BSLMM, and how phenotypes are predicted from sparse and random effects. It also specifies accuracy measures and parameter-estimation procedures for simulated and real data.

  • PVE Estimation with LMM, BVSR and BSLMM: PVE is estimated from posterior means for BVSR and BSLMM, whereas LMM uses a REML-based approximation involving the variance component.The LMM estimate generalizes a prior formula and applies to any centered relatedness matrix K.
  • PVE Estimation with LMM, BVSR and BSLMM: LMM PVE standard errors are approximated from the second derivative of the restricted likelihood and the delta method.The derivative is evaluated efficiently using recursions.
  • PVE Estimation with LMM, BVSR and BSLMM: LMM-Bayes posterior-mean PVE estimates were nearly identical to LMM estimates, so only LMM results are presented.
  • Real Data: For real-data prediction, BSLMM and BVSR use posterior means of sparse and random effects from MCMC, while LMM uses BLUP for its random effects.Predictions for test observations use the conditional mean of future random effects given the training data.
  • Simulated Data: Simulated-data prediction accuracy is assessed using RPG and correlation, requiring estimates of β in the simple linear model.
  • Simulated Data: BSLMM combines small and sparse effects because β = α + ˜β, with LMM and BVSR recovered when π = 0 and σ2_b = 0, respectively.
  • Simulated Data: For simulated predictions, LMM obtains β through its α estimate, BVSR uses the posterior mean of ˜β, and BSLMM adds Rao-Blackwellised α and MCMC-based ˜β estimates.

Other Methods

The study compares BSLMM with LMM, BVSR, LMM-Bayes, Bayesian Lasso, BayesA-Flex, BayesCπ, and an empirical-Bayes BSLMM approximation. These methods differ in priors, fixed variance components, software, and MCMC settings.

  • Other competing methods: The comparison includes LMM fitted with GEMMA, Bayesian Lasso fitted with BLR, and BayesCπ fitted with GenSel.Bayesian Lasso uses a double-exponential coefficient prior, whereas BayesCπ is fit using online software.
  • BSLMM special cases: BVSR is implemented by fixing ρ = 1 in BSLMM, while LMM-Bayes is implemented by fixing ρ = 0.The BSLMM implementation gives BVSR slightly better results and runs faster than piMASS in all considered examples.
  • BayesA-Flex: BayesA-Flex estimates its scaling variance parameter, unlike BayesA, and this modification greatly improves prediction performance in the evaluated data set.BayesA fixes σ2 to a prespecified value, whereas BayesA-Flex places a prior on σ2 and estimates it from the data; BayesA itself gives poor results.
  • BSLMM-EB: BSLMM-EB fixes σ2_b to the null-model REML estimate, reducing computation but potentially discouraging large-effect SNPs and underestimating their effects.The approximation avoids updating σ2_b during each MCMC iteration and may reduce prediction performance, especially with large-effect SNPs.
  • Computation: MCMC-based methods use 2.1 million iterations with 0.1 million burn-in, except BayesCπ, which uses 1.1 million iterations with the same burn-in.BayesCπ uses fewer iterations because of web-server restrictions.

Text S2 Detailed MCMC Strategy for BSLMM

The BSLMM sampler uses Metropolis-Hastings on hyper-parameters with analytically integrated likelihoods and conditional updates for latent quantities. Eigen decompositions and bounded model size make computation feasible, while PVE and PGE are obtained from posterior samples.

  • Posterior sampling: BSLMM uses Metropolis-Hastings to sample (h, ρ, π, γ), integrating out ˜β, u, and τ analytically to evaluate the marginal likelihood.The sampler operates on (0,1) × (0,1) × (0,1) × {0,1}p.
  • Likelihood computation: An eigen decomposition of the relatedness matrix transforms phenotypes and genotypes, enabling efficient determinant and matrix-inverse calculations during likelihood evaluation.The decomposition is performed at the beginning of MCMC using K = UDUT.
  • Proposal distributions: The sampler combines rank-based add, remove, and switch proposals for γ with random-walk updates for h, ρ, and log(π), plus small-world proposals for longer-range moves.The add, remove, and switch moves have probabilities 0.4, 0.4, and 0.2; small-world proposals occur with probability 0.33 and compound 1–20 local moves.
  • Conditional updates: For each sampled hyper-parameter state, BSLMM samples τ, ˜β, and u from conditional distributions, sampling UT u elementwise rather than directly from an n-dimensional normal.The transformed conditional distribution has a diagonal covariance matrix.
  • Derived quantities: Posterior samples of ˜β, u, and τ yield samples of PVE and PGE, while the special case K = XXT /p supports a Rao-Blackwellised posterior mean for α.The α approximation avoids sampling α at every iteration and requires one p-dimensional matrix-vector multiplication at the end.
  • Computational considerations: MCMC computation is linear in the number of individuals and quadratic in |γ| per iteration, so |γ| is capped at 300 for simulations and human data and 600 for mouse data.The calculation of Ω is the most time-consuming component when |γ| is large.

Text S3 The Probit BSLMM and Binary Traits

The section extends BSLMM to binary traits through a probit-linked latent-variable model and develops posterior sampling for its parameters. It also reports that ordinary BSLMM outperformed probit BSLMM on three tested traits and derives a liability-scale PVE correction for case-control data.

  • Probit BSLMM: Probit BSLMM models binary outcomes with a probit link and represents them using latent normal variables thresholded at zero.The latent model uses zi = µ + xi ˜β + ui + ǫi with ǫi ∼ N(0, 1), assigning yi = 1 when zi > 0.
  • Posterior inference: Posterior inference adds truncated-normal sampling for latent variables z to an MCMC strategy otherwise similar to the quantitative-trait BSLMM.The latent variables are sampled left-truncated at zero for yi = 1 and right-truncated at zero for yi = 0, with τ = 1.
  • Performance comparison: For all three tested traits, treating binary values as quantitative with BSLMM performed better than directly modeling them with probit BSLMM.Performance was compared using Brier score differences, where positive values indicate worse performance than BSLMM.
  • PVE correction: A Taylor-series derivation provides a correction factor, ka(1−ka)z2p, for converting observed-scale PVE estimates to the liability scale.The derivation assumes a probit population model and approximates the expected binary label using a linear mixed model.

Figures

The figures compare LMM, BVSR, and BSLMM for PVE estimation and phenotype prediction across simulations, human disease data, and mouse traits using multiple performance metrics.

  • Main comparisons: Figures 1–4 evaluate PVE estimation and prediction across simulation scenarios and WTCCC and heterogeneous-stock mouse datasets.Simulation analyses vary causal-SNP architecture; empirical analyses cover seven diseases and three mouse traits.
  • PVE estimation: Figure 1 compares LMM, BVSR, and BSLMM PVE estimates across two simulation scenarios with true PVE values of 0.2 and 0.6.RMSE and estimated-PVE boxplots are based on 20 replicates per case.
  • Simulation prediction: Figures 2 and S1–S3 compare simulation prediction using Relative Predictive Gain or correlation, with causal SNPs either included or excluded.The simulations use true PVE=0.6 and 20 replicates; scenarios vary causal or medium/small-effect SNP counts.
  • WTCCC prediction: Figures 3 and S4 compare prediction for seven WTCCC diseases using AUC and Brier score, respectively.Higher AUC indicates better performance, whereas lower Brier score indicates better performance; both figures order diseases by BSLMM performance.
  • Mouse-trait prediction: Figures 4 and S5 compare several models with BSLMM for three heterogeneous-stock mouse traits using RMSE and correlation differences.Each train/test split uses 20 replicates; positive RMSE difference indicates worse performance than BSLMM, while positive correlation difference indicates better performance than BSLMM.

CD8 MCH BMI

Figure S6 compares BSLMM with probit BSLMM for three binary traits in heterogeneous stock mouse data. Performance is evaluated using Brier-score differences across two train/test splits with 20 replicates each.

  • CD8 MCH BMI: Figure S6 compares prediction performance between BSLMM and probit BSLMM for three binary traits in heterogeneous stock mouse data.Positive Brier-score differences indicate worse performance than BSLMM.
  • CD8 MCH BMI: The comparison uses two train/test splitting schemes, each evaluated over 20 replicates.
  • CD8 MCH BMI: The mean Brier scores of BSLMM across the six cases are 0.185, 0.205, 0.201, 0.236, 0.245 and 0.249.

Tables · Selected References

Table 1 summarizes some effect size distributions proposed for polygenic modeling.

  • Tables: Table 1 summarizes proposed effect size distributions for polygenic modeling.The table provides an overview of selected distributional assumptions used in polygenic modeling.
  • Tables: The table focuses on effect size distributions rather than broader polygenic modeling methods.
  • Tables: Table 1 presents a selection of distributions that have been proposed in the polygenic modeling literature.
Loading 1209.1341v2…