Efficient mixed model approach for large-scale genome-wide association studies of ordinal categorical phenotypes

Abstract

In genome-wide association studies, ordinal categorical phenotypes are widely used to measure human behaviors, satisfaction, and preferences. However, because of the lack of analysis tools, methods designed for binary or quantitative traits are commonly used inappropriately to analyze categorical phenotypes. To accurately model the dependence of an ordinal categorical phenotype on covariates, we propose an efficient mixed model association test, proportional odds logistic mixed model (POLMM). POLMM is computationally efficient to analyze large datasets with hundreds of thousands of samples, can control type I error rates at a stringent significance level regardless of the phenotypic distribution, and is more powerful than alternative methods. In contrast, the standard linear mixed model approaches cannot control type I error rates for rare variants when the phenotypic distribution is unbalanced, although they performed well when testing common variants. We applied POLMM to 258 ordinal categorical phenotypes on array genotypes and imputed samples from 408,961 individuals in UK Biobank. In total, we identified 5,885 genome-wide significant variants, of which, 424 variants (7.2%) are rare variants with MAF < 0.01.

Keywords: GRM; GWAS; POLMM; PheWAS; UK Biobank; food and other preferences; genetic relationship matrix; genome-wide association studies; mixed model approach; ordinal categorical data; phenome-wide association studies; proportional odds logistic mixed model; saddlepoint approximation; unbalanced phenotypic distribution.

Figure 1

Empirical type I error rates of POLMM, BOLT-LMM, and fastGWA methods at a significance level 5$\times$10⁻⁸ We simulated 1,000 families with a total sample size $n=$ 10,000 and an ordinal categorical phenotype including four levels with sample sizes $n_1$, $n_2$, $n_3$, and $n_4$. From left to right, the plots consider four scenarios: balanced $\left(n_1:n_2:n_3:n_4=1:1:1:1\right)$, moderately unbalanced $\left(n_1:n_2:n_3:n_4=10:1:1:1\right)$, unbalanced $\left(n_1:n_2:n_3:n_4=30:1:1:1\right)$, and extremely unbalanced $\left(n_1:n_2:n_3:n_4=100:1:1:1\right)$. From top to bottom, the plots consider three variance components, tau, $\tau =$ 0.5, 1, and 2. We simulated common, low-frequency, and rare variants with MAFs of 0.3, 0.01, and 0.005, respectively. In total, 10⁹ replications were conducted in each scenario.

Figure 2

Empirical powers of POLMM, SAIGE, BOLT-LMM, and fastGWA methods at significance level 5$\times$10⁻⁸ We simulated 1,000 families with a total sample size $n=$ 10,000 and an ordinal categorical phenotype including four levels with sample sizes $n_1$, $n_2$, $n_3$, and $n_4$. From left to right, the plots consider four scenarios: balanced $\left(n_1:n_2:n_3:n_4=1:1:1:1\right)$, moderately unbalanced $\left(n_1:n_2:n_3:n_4=10:1:1:1\right)$, unbalanced $\left(n_1:n_2:n_3:n_4=30:1:1:1\right)$, and extremely unbalanced $\left(n_1:n_2:n_3:n_4=100:1:1:1\right)$. From top to bottom, the plots consider two MAFs of 0.3 and 0.01 to simulate common and low-frequency variants. We let the variance component $\tau =1$. For SAIGE, we dichotomize phenotype as 0 or 1 depending on whether the subject is in level 1 or not. For BOLT-LMM, the empirical powers were calculated on the basis of the empirical significance levels because it cannot control type I error rates for low-frequency variants.

Figure 3

Manhattan plots for UK Biobank data analysis The left panels show Manhattan plots based on BOLT-LMM, the middle panels show Manhattan plots based on FastPOLMM-NoSPA, and the right panels show Manhattan plots based on FastPOLMM. The redline represents the genome-wide significance level 5$\times {10}^{-8}$.