pyreml: a general-purpose, GPU-accelerated mixed model solver in Python

Linear mixed models (LMM) are a popular modeling and inference framework, widely used in bioinformatics and other scientific fields. To obtain an unbiased estimation of the variance components, an LMM has to be fitted using restricted maximum likelihood (REML), a computationally intensive step. Most REML software is closed source, and currently there is no Python library offering complete LMM functionalities, such as the specification of spatial and genomic variance structures, multivariate analyses, random regression and factors, and heteroscedastic residuals. Here we propose pyreml, an LMM solver that fills this gap, comes with a user-centered interface and speeds up calculations by running on GPU. We show that pyreml fits a genomic best linear unbiased prediction (GBLUP) model more than one order of magnitude (>12\times) faster than the current state-of-the-art, closed source REML software. We also show that \pyreml can estimate a latent factor-analytic variance structure, a challenging feature that no current open source REML solution implements. Finally, we illustrate the use of pyreml to analyze longitudinal data, and to perform spatial analysis. Anchored by 500+ functional tests assessing the equivalence of its outputs to those from reference software (lme4, rrBLUP and BLUPF90), we demonstrate pyreml to be a performant, reliable and open-source solution available for academic, research and industrial LMM users.