Predicting genome-wide functional constraints with GPN-Star

A fundamental challenge in biology is understanding the functional significance of genetic variants. Despite tremendous advances in experimental technologies in genomics, determining which variants influence phenotypes or contribute to diseases remains difficult. Evolution provides valuable insights: deleterious mutations tend to be purged by natural selection, although tolerated or advantageous changes may accumulate. Millions of years of evolution therefore provide a genome-wide record of functional constraint for variant interpretation. This idea has deep roots in molecular evolution and comparative sequence analysis⁶. More recently, large comparative genomics studies have quantified evolutionary constraints across the genomes of hundreds of species^7,8.

Genomic language models (gLMs) have emerged as a promising approach for extracting evolutionary information directly from raw DNA sequences through self-supervised learning (ref. ¹ and references therein). By predicting masked nucleotides from sequence context, gLMs estimate per-site likelihoods that reflect evolutionary constraint without requiring labelled data. These likelihoods have been shown to be effective predictors of genome-wide variant effects^1,9. However, gLMs based on standard language modelling frameworks still underperform compared with much simpler classical phylogenetic models on certain variant interpretation tasks—particularly in complex eukaryotic genomes such as in humans and in distal regulatory elements such as enhancers²—even with massive model sizes^3,4.

A long-standing approach to modelling evolutionary data is to construct multiple sequence alignments (MSAs). By algorithmically aligning homologous loci across biological sequences, MSAs show site-specific patterns of conservation and facilitate the inference of evolutionary preferences for variants at specific positions. Protein MSAs have enabled highly successful models, including AlphaFold¹⁰, MSA Transformer¹¹ and EVE¹², and have recently experienced renewed interest as scaling single-sequence protein language models has shown diminishing returns^13,14,15. Extending this concept from proteins to genomes, whole-genome alignments (WGAs) of complete assemblies from tens to hundreds of species enable genome-wide studies of sequence evolution^16,17. Classical parametric phylogenetic models fitted on WGAs—such as PhastCons¹⁸ and PhyloP¹⁹—have long been staples in variant interpretation for humans and other species. Recent large-scale consortia have further expanded high-quality WGA resources^20,21.

Therefore, combining the gLM framework with WGA data presents a promising direction, and we recently demonstrated this potential with GPN-MSA²² trained on a WGA of 90 vertebrate species. Here we present GPN-Star (genomic pretrained network with species tree and alignment representations), a general gLM framework that more effectively exploits multispecies WGAs through a new model architecture with phylogeny-aware designs; see Supplementary Information for key differences between GPN-MSA and GPN-Star. This flexible framework enabled us to apply GPN-Star to three WGAs relative to the human genome, spanning vertebrate, mammal and primate evolutionary timescales. GPN-Star achieves state-of-the-art performance across coding and non-coding variant interpretation tasks. The vertebrate model consistently outperforms GPN-MSA trained on the same alignment, whereas the mammal and primate models excel particularly at predicting non-coding variant effects. We further evaluate its use in human genetics by analysing the heritability of complex traits. Importantly, GPN-Star trained on the primate WGA prioritizes variants with an unprecedented amount of heritability enrichment across more than a hundred complex traits²³. Also, we find a striking connection between the most informative evolutionary timescale for constraint prediction and the effective polygenicity of traits.

To demonstrate the generality of our approach, we apply GPN-Star to five model organisms with minimal tuning and show its effectiveness in assessing variant effects in these species.

GPN-Star learns functional constraints by leveraging evolutionary signals embedded in multispecies phylogeny and WGAs. Inspired by classical evolutionary models, it characterizes how genomic positions evolve across species. To achieve this, we introduce a specialized transformer architecture that integrates sequence context and evolutionary information (Fig. 1a).

GPN-Star uses an encoder-only architecture trained with a masked language modelling (MLM) objective. Each input consists of a WGA window and its species tree. The alignment is partitioned into target and source sequences and the model predicts masked nucleotides in the targets conditioned on both their sequence context and the evolutionary information in the sources. This is accomplished through a stack of encoder blocks, each comprising a sequence-wise self-attention module to encode intrasequence context, a phylogeny-informed cross-attention module to encode evolutionary context from the source sequences and a feed-forward network to integrate the information. To efficiently represent phylogenetic structure, we group closely related source species according to their evolutionary distances and pool their sequences into clade-level representations. The cross-attention module adaptively weights source sequence contributions using evolutionary distances derived from the species tree using an attention mechanism. The resulting embeddings for the target sequences are then passed to an MLM head for nucleotide probability prediction. Full architectural details are provided in Methods.

Another methodological innovation of this work is that it accounts for mutation rate variation in gLMs. Because GPN-Star is trained on genomic sequences observed in nature, its predictions reflect both mutation and selection. To disentangle these signals and isolate selective constraints, we normalized the prediction by the model of each variant against its background mutation effect estimated from high-confidence neutral sites in the genome (Methods). Unless otherwise specified, all analyses in this study use these calibrated scores to capture functional constraints with minimal mutation rate biases.

The GPN-Star framework is general and flexible, designed to work with any alignment data from any species, requiring minimal hyperparameter tuning to achieve robust performance (Methods). We first applied it to the human genome, training three separate GPN-Star models using the largest available vertebrate, mammalian and primate WGAs (Fig. 1b; Methods). Throughout this study, we denote these models as GPN-Star (V), GPN-Star (M) and GPN-Star (P), respectively. We have experimented with different model sizes, but we focus here on the largest models (200 million parameters), trained for several days on 8 NVIDIA A100 GPUs (Supplementary Table 1). Existing single-sequence gLMs such as Nucleotide Transformer (trained on 128 A100 GPUs for a month)³ and Evo 2 (trained on over 2,000 H100 GPUs for several months)⁴ require substantial computational resources to implicitly learn evolutionary constraints, in exchange for flexibility and broader applications. We show that GPN-Star achieves superior constraint prediction with a much smaller resource footprint by leveraging the explicit evolutionary context provided by alignment data.

Unlike previous gLMs trained across extremely broad evolutionary distances (for example, spanning prokaryotes to humans), GPN-Star focuses on narrower, more recent phylogenetic distances closely related to humans (Fig. 1b). As we demonstrate below, modelling longer evolutionary histories is not always helpful. Instead, capturing recent evolutionary constraints proves especially advantageous for interpreting certain classes of genetic variants (Fig. 1c).

A critical challenge in human genetics is understanding the impact of genetic variants on disease susceptibility—an essential step towards improving the diagnosis of genetic disorders, identifying drug targets and realizing the promise of precision medicine. Many computational models have been developed to predict variant pathogenicity in the human genome. Here we systematically evaluate the performance of GPN-Star across a comprehensive set of benchmarks.

Our evaluation of coding variants focused on missense variants, the most prevalent class. We first considered classifying pathogenic versus benign variants in ClinVar, a widely used clinical variant database containing expert-curated pathogenicity labels²⁴. We compared our models with major genome-wide variant effect predictors, including classical evolutionary methods (PhyloP and PhastCons) that also use WGA, the ensemble model CADD²⁵, recent single-sequence gLMs (Nucleotide Transformer 2.5B multispecies and Evo 2 40B) and GPN-MSA. GPN-Star (V) achieved the highest area under the precision-recall curve (AUPRC), matching the performance of the protein language model ESM-1b²⁶ while outperforming ESM-2 (ref. ²⁷) and ESM-3-open²⁸ (Fig. 2a). For ease of visualization, only the best-performing GPN-Star model is shown in Fig. 2; results for all three versions are available in Extended Data Fig. 1.

We next considered somatic cancer variants from the COSMIC database²⁹, constructing a benchmark to distinguish missense variants frequently observed in tumours from common missense variants in the general population (gnomAD³⁰). As shown in Fig. 2b, GPN-Star (V) substantially outperformed all competing models, demonstrating strong predictive capability for pathogenicity beyond germline variants.

Beyond clinical variants, we also assessed performance on deep mutational scanning (DMS) data, functional assays that measure the fitness effects of all possible missense mutations in a given protein. Across 31 human DMS datasets from ProteinGym³¹, GPN-Star (V) outperformed all genome-wide models, although it slightly lagged behind the protein-specific models (Fig. 2d).

We also compared GPN-Star with two recent missense variant effect predictors, AlphaMissense³² and PrimateAI-3D³³. These models were supervised on population allele frequency data, a highly informative signal for pathogenicity prediction, which is also used to define labels in databases such as ClinVar³⁴. Interestingly, although both AlphaMissense and PrimateAI-3D achieved top performance on ProteinGym and ClinVar, a simple post hoc adjustment of GPN-Star (V) predictions using gnomAD allele frequencies boosted its performance on ClinVar to surpass PrimateAI-3D and approach that of AlphaMissense (Fig. 2d and Supplementary Fig. 1).

We next considered tasks involving non-coding variants, which are known to be especially challenging for predictive modelling. Here we demonstrate that GPN-Star is a powerful tool for identifying pathogenic non-coding variants in the human genome.

In addition to previous genome-wide variant effect prediction (VEP) models, we included three prominent sequence-to-function models—Enformer³⁵, Borzoi³⁶ and AlphaGenome³⁷—which are trained on extensive functional genomics data and have been widely applied to non-coding variant interpretation. We evaluated the models on classifying non-coding pathogenic variants from the OMIM³⁸ and HGMD³⁹ databases, both of which contain expert-curated annotations of human disease-associated variants. GPN-Star (M) achieved the best performance on both benchmarks (Fig. 2e,f). Notably, the sequence-to-function models performed substantially worse than alignment-based gLMs and phylogenetic models on this task, consistent with recent findings^2,37. Stratified analysis indicated that GPN-Star consistently achieves top performance across variant types, with the largest advantage over sequence-to-function models in distal enhancers (Supplementary Figs. 2, 3).

Given the critical role of promoter regions in transcription initiation and gene regulation, many specialized models have been developed to predict the impact of promoter variants. We therefore evaluated GPN-Star on promoter variants in OMIM and compared its performance with other methods, including three promoter-specific models: PromoterAI⁴⁰, SpeciesLM⁴¹ and GPN-Promoter². As shown in Fig. 2h, GPN-Star (M) demonstrated superior predictive performance compared with all competing models. In particular, the margin of improvement over the promoter-specific models is substantial.

Genome-wide association studies (GWAS) have been instrumental in identifying variants that contribute to genetic disease susceptibility. To further assess the use of GPN-Star, we evaluated its performance in classifying putatively causal versus non-causal missense variants resulting from fine-mapping of GWAS variants across 65 traits from the UK Biobank⁴². Among all competing models, GPN-Star (M) achieved the highest predictive performance on these fine-mapped missense variants. Notably, despite leveraging population allele frequency information, both AlphaMissense and PrimateAI-3D were substantially outperformed by GPN-Star in this task (Fig. 2c).

We next evaluated the models on fine-mapped non-coding GWAS variants across 83 traits from the UK Biobank⁴², again assessing their ability to distinguish putatively causal versus non-causal variants. GPN-Star (M) continued to outperform all other models on this benchmark (Fig. 2g). Although sequence-to-function models (Enformer, Borzoi and AlphaGenome) exhibited moderate performance, Evo 2 showed relatively limited predictive value, as was previously observed². Stratified analysis showed robust performance of GPN-Star across variant types (Supplementary Fig. 4). For fine-mapped variants located in promoter regions, GPN-Star (M) again outperformed all models, including the promoter-specific models PromoterAI, SpeciesLM and GPN-Promoter (Fig. 2h).

Having observed the strong performance of GPN-Star in predicting pathogenic and fine-mapped variants, we explored its use in rare variant association testing (RVAT), an important but challenging task in statistical genetics. The great interest in RVAT stems from the fact that most large-effect variants tend to be rare, as they are subject to negative selection pressure⁴³. To overcome statistical power issues with low-frequency variants, RVAT is typically carried out at the gene level by aggregating variants in each gene⁴⁴. These testing procedures often rely on variant annotations that reflect functional importance to prioritize variants in the aggregation. DeepRVAT⁵ is a recent method that uses a deep set network to integrate such variant annotations for RVAT. By combining an expressive deep learning framework with a powerful set of variant annotations, it demonstrated improved statistical power and computational efficiency compared with previous methods in extensive evaluations on UK Biobank whole-exome sequencing (WES) data⁵.

We conducted RVAT experiments using DeepRVAT, enhancing the annotations used in the published version with all three GPN-Star predictions during both the training and association testing phases. Following the published benchmarking procedure⁵, we performed RVAT experiments using DeepRVAT on 34 quantitative traits in UK Biobank with the WES data from 161,822 unrelated individuals of European ancestry, retaining variants with minor allele frequency (MAF) <0.1% (Methods). To account for differences in stochastic model initialization, this procedure was run three times with different random seeds. As shown in Fig. 2j, adding the predictions from the GPN-Star models as variant annotations into DeepRVAT resulted in an increase in the number of discovered genes at a family-wise error rate <0.05 (on average across runs, 402 genes compared with 383 from the original method). Among the significant discoveries, there were also more gene–phenotype associations that replicated in conventional RVAT studies on UK Biobank with larger sample sizes^45,46 (Methods) (on average 353 compared with 338), indicating high robustness of the further discoveries. Notably, the original DeepRVAT already includes state-of-the-art missense variant effect predictors such as AlphaMissense and PrimateAI, as well as the sequence-to-function model DeepSEA⁴⁷, in the annotations. Nonetheless, GPN-Star seemed to offer complementary information and boosted the power of the tests.

To investigate the prediction of causal variants for complex trait GWAS, we turned to stratified linkage disequilibrium score regression (S-LDSC)⁴⁸. S-LDSC is a principled way to estimate the informativeness of an annotation for complex trait heritability while leveraging the signal from all single nucleotide polymorphisms (SNPs), including those not confidently fine-mapped or genome-wide significant. It also serves as a backbone for functionally informed fine-mapping⁴⁹ and polygenic risk scores⁵⁰. Because S-LDSC requires scoring about 10 million variants, it was feasible to evaluate only the most scalable models (or those with precomputed scores). We binarized model scores to select a specific fraction of common variants. In our main analysis, we used the top 0.1% most constrained variants. We ran S-LDSC separately for each model score while conditioning on 96 baseline features and meta-analysed the results across 106 independent traits²³.

Since the early work on S-LDSC, conservation scores have been found to be the most enriched annotation for complex trait heritability⁴⁸. More recently, primate-specific conservation—particularly PhastCons (P)—has emerged as the state-of-the-art for complex trait heritability enrichment^7,8. Notably, GPN-Star (P) substantially improves on this result, followed by GPN-Star (M) (Fig. 3a). These improvements are even more striking when examining the heritability coefficient \({\tau }^{\star }\) (Fig. 3a), which quantifies the unique contribution of an annotation to heritability after adjusting for baseline features (Methods). Although per-trait enrichment estimates are noisier, GPN-Star (P) or GPN-Star (M) consistently rank at the top (Supplementary Fig. 5). This advance is particularly notable given the long-standing role of complex trait heritability enrichment as a meaningful benchmark in human genetics.

The improvement of GPN-Star persists when restricting to top-ranked variants in coding and non-coding regions separately (Fig. 3b,c), with stronger gains in non-coding regions. We further performed an analysis of the best-performing model GPN-Star (P), stratifying the non-coding region by functional category, and found particularly strong improvements in enhancer regions (Extended Data Fig. 2). Moreover, across a range of binarization thresholds and all three evolutionary timescales, GPN-Star consistently outperforms previous conservation scores (Fig. 3d).

As the final analysis of complex trait heritability, we investigated the performance of tissue-specific scores. Sequence-to-function models such as Enformer³⁵ and Borzoi³⁶ do not provide a single variant effect score but instead predict changes in activity across thousands of functional genomics tracks representing different assays and tissues. In our comparisons, we used Enformer and Borzoi scores aggregated across all tracks (tissue agnostic) as well as in seven specific tissues, following procedures from ref. ⁵³ Given the low enrichment of tissue-agnostic Enformer and Borzoi scores when meta-analysed across all traits (Fig. 3a), we investigated the performance of tissue-specific scores. We meta-analysed these scores only across traits where we expected the tissue to be relevant (Supplementary Table 2). Tissue-specific Enformer and Borzoi scores consistently outperform their tissue-agnostic scores (Fig. 3g, Extended Data Fig. 4 and Supplementary Fig. 6).

On the other hand, GPN-Star trained solely on DNA sequences is intrinsically tissue agnostic. Motivated by this observation, we devised a simple approach to incorporate tissue specificity into GPN-Star by only considering top-scoring variants located near tissue-specifically expressed genes (SEG) (inspired by LDSC-SEG⁵⁴) or in tissue-specific candidate cis-regulatory elements (cCREs) (obtained from ENCODE v.4⁵³). Although this improved performance for most tissues—such as brain and blood/immune—it did not consistently surpass the performance of tissue-agnostic GPN-Star annotations, potentially because the present approach to incorporating tissue specificity remains suboptimal (Fig. 3g and Extended Data Fig. 4). Nevertheless, GPN-Star outperformed Enformer and Borzoi across all seven tissues. We also verified that the annotation by SEG and cCRE alone without model scores exhibited much lower heritability enrichment, confirming that the observed performance is still largely attributed to GPN-Star (Fig. 3g,h and Supplementary Fig. 6). Zooming in on example traits, brain-specific GPN-Star showed the highest enrichment for schizophrenia, blood/immune-specific GPN-Star performed best for lupus, and both tissue-agnostic and liver-specific GPN-Star yielded similar results for IGF1 (Fig. 3h). These findings highlight the importance of incorporating tissue-specific information in variant effect prediction.

We examined the representations learned by GPN-Star and found that its embeddings distinguish major genomic functional elements, with stronger separation for evolutionarily conserved regions (Supplementary Information; Supplementary Figs. 7–9). These results suggest that GPN-Star is aware of core functional elements of the genome when making predictions.

Site-independent models such as PhyloP cannot, by definition, capture dependencies among nucleotides. By contrast, alignment-free gLMs have been shown to learn dependencies among well-known interacting elements, for example, in a transcription-factor binding site (TFBS) motif or between splice donors and acceptors⁴¹. To further probe understanding in GPN-Star of genomic syntax, we analysed learned nucleotide dependencies by systematically mutating each position and quantifying the resulting probability changes at other positions in the sequence^41,55. This analysis showed biologically meaningful co-evolutionary signals, including dependencies among TFBSs, coding regions and splice sites, as well as evidence for primate-specific regulatory constraint. Representative examples are shown in Extended Data Fig. 5, with more case studies presented in Supplementary Figs. 10–13.

These results demonstrate that GPN-Star can leverage co-evolutionary signals embedded in genomic sequence context to learn meaningful nucleotide dependencies that align with known functional dependencies, representing a notable advance over traditional conservation scores. Finally, we used nucleotide dependency to examine a few cases in our disease variant benchmarks where GPN-Star correctly predicts them to be deleterious, whereas PhyloP and PhastCons predict neutral scores (Extended Data Figs. 6, 7 and 8). We found that these variants do not exhibit strong site-wise conservation but reside in functionally important elements with distinct sequence contexts shown by nucleotide dependencies. As schematically depicted in Extended Data Fig. 9, GPN-Star effectively leverages these contextual signals to outperform traditional conservation metrics. For example, in the most extreme case of perfectly conserved sites, traditional conservation scores are unable to distinguish different variant types, whereas GPN-Star tends to predict more deleterious scores for higher impact variants (Supplementary Fig. 14).

To more directly assess the connection between model predictions and evolutionary constraints in the genome, we leveraged allele frequency data from gnomAD v.3.1.2, which aggregates whole-genome sequencing samples from 76,156 human individuals³⁰. Allele frequencies in the human population serve as informative indicators of selective constraint: more deleterious alleles tend to have lower frequencies as a result of purifying selection^9,22,30.

In this evaluation, we focused on comparison with PhyloP and PhastCons, which also learn evolutionary constraints from WGA data. To assess how well each model captures the relationship between allele frequency and constraint, we considered the vertebrate, mammal and primate versions of each model and obtained predictions for all gnomAD v.3 variants on chromosome 22, the held-out chromosome not used in training GPN-Star. We then compared the mean minor allele frequencies in several quantile bins defined by different models. As shown in Fig. 4a, across all three evolutionary timescales, variants in lower GPN-Star quantile bins have consistently lower average allele frequencies compared with those in corresponding PhyloP and PhastCons bins, suggesting that GPN-Star more accurately captures selective constraints in the human genome.

We next performed a more quantitative evaluation focusing on the most deleterious tails of the model score distributions, a regime especially important in many human genetics applications. We quantified the enrichment of rare variants relative to common variants in the most constrained tail predicted by each model. Here we defined rare variants as singletons and common variants as those with allele frequency >5%. Because rare variants are, on average, more deleterious, a model with more accurate constraint predictions should show higher enrichment of singletons in its most constrained tail. As shown in Fig. 4b and Supplementary Fig. 15a, all three GPN-Star models yielded substantially higher enrichment of rare variants than did PhyloP, PhastCons or CADD. Among the GPN-Star models, the vertebrate model showed the strongest enrichment overall and outperformed GPN-MSA (also trained on vertebrate genomes). When stratifying variants by molecular consequences, GPN-Star again achieved the highest enrichment in every category (Fig. 4c). Notably, GPN-Star (V) performed the best for missense variants, whereas GPN-Star (M) led for the synonymous and non-coding categories, mirroring trends observed in previous benchmarks. Lastly, our procedure for accounting for local mutation rate variation substantially reduced correlations with mutation rate estimates while improving performance across most downstream benchmarks (Fig. 4d and Extended Data Fig. 10; Supplementary Information).

GPN-Star requires only a multispecies WGA and the corresponding phylogeny to predict variant effects, making it readily applicable to species for which genome-wide variant interpretation has remained challenging. To demonstrate this, we trained GPN-Star models on five important model organisms: Mus musculus, Gallus gallus, Drosophila melanogaster, Caenorhabditis elegans and Arabidopsis thaliana. In each case, we showed that the learned constraints are highly informative for interpreting genetic variation in that species. For each of the five species, we collected a WGA dataset comprising aligned genomes from 18 to 135 species and applied training procedures similar to those used for the human models, with only minor adjustments (Supplementary Table 1).

Owing to the scarcity of large-scale curated variant datasets with functional annotations in non-human species, we used population genetic data for evaluation. Specifically, we gathered population-level variant data for each of the five species (Methods) and examined the enrichment of rare variants in the most deleterious tail of the prediction distributions. Across all five species, GPN-Star scores exhibited substantially higher rare variant enrichment compared with PhyloP and PhastCons (Fig. 5a), indicating more accurate predictions of genome-wide evolutionary constraints. Across different variant categories based on molecular consequences, GPN-Star consistently showed the highest enrichment in almost all cases (Supplementary Figs. 15b–f and 16). In A. thaliana, the performance gap between GPN-Star and other models was smaller compared with the results for the other species, which could be due to either the small alignment or the known lower quality of WGAs in plants⁵⁶.

For three of the species, we were able to gather curated pathogenic variants to carry out further evaluation. For the M. musculus model, we used pathogenic variants from the MMRdb database (Methods). GPN-Star outperformed both PhyloP and PhastCons in distinguishing these pathogenic variants from common variants in the population, both genome-wide (Fig. 5b) and in each variant category (Supplementary Fig. 16). Similarly, we collected experimentally validated lethal variants for D. melanogaster (from FlyBase⁵⁷) and for C. elegans (from an experimental study⁵⁸) and found that GPN-Star consistently outperformed the other models on both datasets (Fig. 5c,d and Supplementary Fig. 16). We also evaluated the mouse output head of Enformer³⁵ on the non-coding variants in the mouse benchmarks (Supplementary Fig. 17). It was markedly outperformed by GPN-Star, echoing the non-coding VEP results for humans.

Beyond variant prioritization, GPN-Star also enables exploration of learned functional elements and their co-evolution through nucleotide dependency analysis. As an example, we highlight learned nucleotide dependencies between TFBS at the well-characterized MSE enhancer locus in D. melanogaster (Fig. 5e). The coordinated activity of these transcription factors drives precise spatial patterning across the fly embryo⁵⁹. This showcases how GPN-Star can serve as an efficient and low-cost tool for investigating functional elements and their dependencies across genomes, complementing and guiding experimental studies.

GPN-Star is a phylogeny-informed genomic language modelling framework that consistently outperforms existing methods for predicting functional constraint and deleterious variants across the human genome. Notably, GPN-Star achieved particularly strong predictive performance for distal enhancer variants, suggesting that evolutionary constraint information derived from cross-species comparisons can provide valuable signals for modelling this historically challenging class of regulatory element⁶⁰. Beyond predictive accuracy, it learns biologically meaningful representations spanning functional elements from enhancers to TFBS and their co-evolutionary dependencies—all without any supervision. Classical WGA-based models, such as PhastCons¹⁸ and PhyloP¹⁹, have been indispensable tools to biologists and clinicians for over two decades. Our results indicate that GPN-Star complements and extends these tools, offering improved performance across diverse species and alignments.

A central insight from our study is that the evolutionary timescale represented in the training data strongly influences the constraint learned by gLMs (Figs. 1c and 2i), corroborating and extending previous observations for classical phylogenetic methods⁷. Models trained on deeper timescales better capture coding and other highly constrained regions, whereas shallower timescales better inform rapidly evolving regulatory elements, including those constrained on primate-specific regulation⁸. The evolutionary timescale should therefore be considered an important design consideration for future gLMs.

Like classical conservation scores, GPN-Star predictions quantify evolutionary constraint at specific timescales rather than general variant pathogenicity. Appropriate timescales should therefore be chosen according to the biological question. For integrating or selecting among models in predictive or discovery workflows, we recommend data-driven approaches that learn appropriate weights and transformations from task-specific data, as we demonstrated with DeepRVAT in our RVAT analysis.

GPN-Star achieves state-of-the-art performance in functional constraint prediction with substantially smaller model and context sizes compared with existing single-sequence gLMs. This effectiveness probably stems from the explicit homology information provided by the alignment, which is not readily available to single-sequence models. Although we observed modest performance gains from increased model sizes, improvements from expanding context sizes were relatively small (Supplementary Fig. 18). A plausible explanation is that WGAs are constructed relative to a reference species and often consist of small, highly fragmented synteny blocks. As a result, when context size increases, nucleotides in a window may not be contiguous in the actual genomes of non-reference species, introducing misleading context information to the model. Exploring how to leverage larger context sizes with WGA data is a promising direction for future research. Advancement in model architecture design or alignment data structure is probably necessary.

A practical limitation is that GPN-Star requires WGAs during inference. To facilitate its use, we have released genome-wide predictions through public repositories. Also, the reliance on a fixed alignment format makes the model less suitable for analysing indels and structural variants. On the other hand, single-sequence gLMs offer greater flexibility and applicability to a broad range of tasks and can naturally leverage large context and process indels and structural variants. We envision that future work can bridge these paradigms potentially through flexible homology retrieval mechanisms, allowing models to dynamically incorporate evolutionary information without relying on rigid alignment structures.

Our results on pathogenicity prediction and complex trait heritability constitute a systematic evaluation of GPN-Star in the context of present knowledge of human genetic variation. However, the true value of the model depends on its capacity to enable biological discovery. We expect GPN-Star to hold great potential in human genetics applications by advancing our understanding of causal variants and human traits. Initial experiments demonstrate improved RVAT, and we expect even greater gains for non-coding analyses, in which variant prioritization remains challenging. GPN-Star may also improve functionally informed fine-mapping, polygenic risk prediction and phenotype prediction.

Beyond identifying evolutionary constraints and deleterious variants, an important future direction is to study human-specific adaptations. The present GPN-Star models focus on capturing evolutionary signals across species, which lack the resolution for studying selection in human populations. We anticipate this to be an important direction for future work and would require learning from even more localized evolutionary data, such as human population genomes and archaic human genomes.

It is remarkable how much can be learned from unlabelled DNA sequences alone. Evolutionary and functional genomics provide complementary views of genetic variation. Although evolutionary modelling outperformed present sequence-to-function approaches in our pathogenicity and heritability analyses, functional genomics remains essential for studying molecular phenotypes and tissue-specific regulation. Integrating these complementary sources of information into multimodal gLMs represents a promising direction.

Large-scale efforts to sequence and align genomes across the tree of life are accelerating^21,61. GPN-Star is well positioned to take advantage of this growth, as it scales readily to new alignments and may gain further capacity from increasingly diverse training data. In this age of rapidly expanding genomic data, we expect GPN-Star to be a valuable tool for leveraging these data to advance our understanding of genetic variation.

Supplementary Information provides details of training data, model architecture, training setup, inference, mutation rate calibration, embedding analysis, nucleotide dependency analysis, evaluation benchmarks, heritability enrichment analysis and RVAT.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Code to train the model, run inference and reproduce the main analyses is available at GitHub (https://github.com/songlab-cal/gpn) and Zenodo (https://doi.org/10.5281/zenodo.21501177)⁶³.

1. Benegas, G., Ye, C., Albors, C., Li, J. C. & Song, Y. S. Genomic language models: opportunities and challenges. Trends Genet. 41, 286–302 (2025).

   Article  CAS  PubMed  Google Scholar 

2. Benegas, G., Eraslan, G. & Song, Y. S. Benchmarking DNA sequence models for causal regulatory variant prediction in human genetics. Preprint at bioRxiv https://doi.org/10.1101/2025.02.11.637758 (2025).

3. Dalla-Torre, H. et al. Nucleotide transformer: building and evaluating robust foundation models for human genomics. Nat. Methods 22, 287–297 (2025).

   Article  CAS  PubMed  Google Scholar 

4. Brixi, G. et al. Genome modelling and design across all domains of life with Evo 2. Nature 652, 1349–1361 (2026).

   Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

5. Clarke, B. et al. Integration of variant annotations using deep set networks boosts rare variant association testing. Nat. Genet. 56, 2271–2280 (2024).

   Article  CAS  PubMed  PubMed Central  Google Scholar 

6. Dayhoff, M. O., Schwartz, R. M. & Orcutt, B. C. A model of evolutionary change in proteins. Atlas Protein Seq. Struct. 5, 345–352 (1978).

7. Sullivan, P. F. et al. Leveraging base-pair mammalian constraint to understand genetic variation and human disease. Science 380, eabn2937 (2023).

   Article  CAS  PubMed  PubMed Central  Google Scholar 

8. Kuderna, L. F. et al. Identification of constrained sequence elements across 239 primate genomes. Nature 625, 735–742 (2024).

   Article  ADS  CAS  PubMed  Google Scholar 

9. Benegas, G., Batra, S. S. & Song, Y. S. DNA language models are powerful predictors of genome-wide variant effects. Proc. Natl Acad. Sci. USA 120, e2311219120 (2023).

   Article  CAS  PubMed  PubMed Central  Google Scholar 

10. Jumper, J. et al. Highly accurate protein structure prediction with AlphaFold. Nature 596, 583–589 (2021).

    Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

11. Rao, R. M. et al. MSA Transformer. In Proc. 38th Int. Conf. on Machine Learning vol. 139 (eds Meila, M. & Zhang, T.) 8844–8856 (PMLR, 2021).

12. Frazer, J. et al. Disease variant prediction with deep generative models of evolutionary data. Nature 599, 91–95 (2021).

    Article  ADS  CAS  PubMed  Google Scholar 

13. Truong, T. Jr & Bepler, T. PoET: a generative model of protein families as sequences-of-sequences. Adv. Neural Info. Process. Syst. 36, 77379–77415 (2023).

    Google Scholar 

14. Yang, K. K. et al. The Dayhoff Atlas: scaling sequence diversity for improved protein generation. Preprint at bioRxiv https://doi.org/10.1101/2025.07.21.665991 (2025).

15. Akiyama, Y. et al. Expanding the scope of protein language modeling to protein-protein interactions with MSA pairformer. Cell 189, 4964–4979 (2026).

    Article  CAS  PubMed  Google Scholar 

16. Blanchette, M. et al. Aligning multiple genomic sequences with the threaded blockset aligner. Genome Res. 14, 708–715 (2004).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

17. Armstrong, J. et al. Progressive Cactus is a multiple-genome aligner for the thousand-genome era. Nature 587, 246–251 (2020).

    Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

18. Siepel, A. et al. Evolutionarily conserved elements in vertebrate, insect, worm, and yeast genomes. Genome Res. 15, 1034–1050 (2005).

    Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

19. Pollard, K. S., Hubisz, M. J., Rosenbloom, K. R. & Siepel, A. Detection of nonneutral substitution rates on mammalian phylogenies. Genome Res. 20, 110–121 (2010).

    Article  CAS  PubMed  Google Scholar 

20. Christmas, M. J. et al. Evolutionary constraint and innovation across hundreds of placental mammals. Science 380, eabn3943 (2023).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

21. Rhie, A. et al. Towards complete and error-free genome assemblies of all vertebrate species. Nature 592, 737–746 (2021).

    Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

22. Benegas, G., Albors, C., Aw, A. J., Ye, C. & Song, Y. S. A DNA language model based on multispecies alignment predicts the effects of genome-wide variants. Nat. Biotechnol. 43, 1960–1965 (2025).

    Article  CAS  PubMed  Google Scholar 

23. Kim, A. et al. Identifying independent causal cell types for human diseases and risk variants. Cell Genom. https://doi.org/10.1016/j.xgen.2026.101325 (2026).

24. Landrum, M. J. et al. ClinVar: public archive of relationships among sequence variation and human phenotype. Nucleic Acids Res. 42, D980–D985 (2014).

    Article  CAS  PubMed  Google Scholar 

25. Rentzsch, P., Witten, D., Cooper, G. M., Shendure, J. & Kircher, M. CADD: predicting the deleteriousness of variants throughout the human genome. Nucleic Acids Res. 47, D886–D894 (2019).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

26. Rives, A. et al. Biological structure and function emerge from scaling unsupervised learning to 250 million protein sequences. Proc. Natl Acad. Sci. USA 118, e2016239118 (2021).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

27. Lin, Z. et al. Evolutionary-scale prediction of atomic-level protein structure with a language model. Science 379, 1123–1130 (2023).

    Article  ADS  MathSciNet  CAS  PubMed  Google Scholar 

28. Hayes, T. et al. Simulating 500 million years of evolution with a language model. Science 387, 850858 (2025).

    Article  ADS  Google Scholar 

29. Tate, J. G. et al. COSMIC: the catalogue of somatic mutations in cancer. Nucleic Acids Res. 47, D941–D947 (2019).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

30. Chen, S. et al. A genomic mutational constraint map using variation in 76,156 human genomes. Nature 625, 92–100 (2024).

    Article  ADS  CAS  PubMed  Google Scholar 

31. Notin, P. et al. ProteinGym: large-scale benchmarks for protein fitness prediction and design. Adv. Neural Info. Process. Syst. 36, 64331–64379 (2023).

    Article  Google Scholar 

32. Cheng, J. et al. Accurate proteome-wide missense variant effect prediction with AlphaMissense. Science 381, eadg7492 (2023).

    Article  CAS  PubMed  Google Scholar 

33. Gao, H. et al. The landscape of tolerated genetic variation in humans and primates. Science 380, eabn8153 (2023).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

34. Ghosh, R. et al. Updated recommendation for the benign stand-alone ACMG/AMP criterion. Hum. Mutat. 39, 1525–1530 (2018).

    Article  PubMed  PubMed Central  Google Scholar 

35. Avsec, Ž et al. Effective gene expression prediction from sequence by integrating long-range interactions. Nat. Methods 18, 1196–1203 (2021).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

36. Linder, J., Srivastava, D., Yuan, H., Agarwal, V. & Kelley, D. R. Predicting RNA-seq coverage from DNA sequence as a unifying model of gene regulation. Nat. Genet. 57, 949–961 (2025).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

37. Avsec, Ž et al. Advancing regulatory variant effect prediction with AlphaGenome. Nature 649, 1206–1218 (2026).

    Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

38. Amberger, J. S., Bocchini, C. A., Schiettecatte, F., Scott, A. F. & Hamosh, A. OMIM.org: Online Mendelian Inheritance in Man (OMIM), an online catalog of human genes and genetic disorders. Nucleic Acids Res. 43, D789–D798 (2015).

    Article  PubMed  Google Scholar 

39. Stenson, P. D. et al. The Human Gene Mutation Database (HGMD): optimizing its use in a clinical diagnostic or research setting. Hum. Genet. 139, 1197–1207 (2020).

    Article  PubMed  PubMed Central  Google Scholar 

40. Jaganathan, K. et al. Predicting expression-altering promoter mutations with deep learning. Science 389, eads7373 (2025).

    Article  CAS  PubMed  Google Scholar 

41. Tomaz da Silva, P. et al. Nucleotide dependency analysis of genomic language models detects functional elements. Nat. Genet. 57, 2589–2602 (2025).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

42. Kanai, M. et al. Insights from complex trait fine-mapping across diverse populations. Preprint at medRxiv https://doi.org/10.1101/2021.09.03.21262975 (2021).

43. Bomba, L., Walter, K. & Soranzo, N. The impact of rare and low-frequency genetic variants in common disease. Genome Biol. 18, 77 (2017).

    Article  PubMed  PubMed Central  Google Scholar 

44. Lee, S., Abecasis, G. R., Boehnke, M. & Lin, X. Rare-variant association analysis: study designs and statistical tests. Am. J. Hum. Genet. 95, 5–23 (2014).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

45. Backman, J. D. et al. Exome sequencing and analysis of 454,787 UK Biobank participants. Nature 599, 628–634 (2021).

    Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

46. Karczewski, K. J. et al. Systematic single-variant and gene-based association testing of thousands of phenotypes in 394,841 UK Biobank exomes. Cell Genom. 2, 100168 (2022).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

47. Zhou, J. & Troyanskaya, O. G. Predicting effects of noncoding variants with deep learning–based sequence model. Nat. Methods 12, 931–934 (2015).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

48. Finucane, H. K. et al. Partitioning heritability by functional annotation using genome-wide association summary statistics. Nat. Genet. 47, 1228–1235 (2015).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

49. Weissbrod, O. et al. Functionally informed fine-mapping and polygenic localization of complex trait heritability. Nat. Genet. 52, 1355–1363 (2020).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

50. Márquez-Luna, C. et al. Incorporating functional priors improves polygenic prediction accuracy in UK Biobank and 23andMe data sets. Nat. Commun. 12, 6052 (2021).

    Article  ADS  PubMed  PubMed Central  Google Scholar 

51. O’Connor, L. J. & Sella, G. Principled measures and estimates of trait polygenicity. Preprint at bioRxiv https://doi.org/10.1101/2025.07.10.664154 (2025).

52. Karollus, A., Mauermeier, T. & Gagneur, J. Current sequence-based models capture gene expression determinants in promoters but mostly ignore distal enhancers. Genome Biol. 24, 56 (2023).

    Article  PubMed  PubMed Central  Google Scholar 

53. Fabiha, T. et al. A consensus variant-to-function score to functionally prioritize variants for disease. Preprint at bioRxiv https://doi.org/10.1101/2024.11.07.622307 (2024).

54. Finucane, H. K. et al. Heritability enrichment of specifically expressed genes identifies disease-relevant tissues and cell types. Nat. Genet. 50, 621–629 (2018).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

55. Zhang, Z. et al. Protein language models learn evolutionary statistics of interacting sequence motifs. Proc. Natl Acad. Sci. USA 121, e2406285121 (2024).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

56. Song, B., Buckler, E. S. & Stitzer, M. C. New whole-genome alignment tools are needed for tapping into plant diversity. Trends Plant Sci. 29, 355–369 (2024).

    Article  CAS  PubMed  Google Scholar 

57. Öztürk-Çolak, A. et al. FlyBase: updates to the Drosophila genes and genomes database. Genetics 227, iyad211 (2024).

    Article  PubMed  PubMed Central  Google Scholar 

58. Qin, Z. et al. Genomic identification and functional characterization of essential genes in Caenorhabditis elegans. Genes Genomes Genet. 8, 981–997 (2018).

    Article  CAS  Google Scholar 

59. Small, S., Blair, A. & Levine, M. Regulation of even-skipped stripe 2 in the Drosophila embryo. EMBO J. 11, 4047–4057 (1992).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

60. Wong, E. S. et al. Deep conservation of the enhancer regulatory code in animals. Science 370, eaax8137 (2020).

    Article  ADS  CAS  PubMed  Google Scholar 

61. Lewin, H. A. et al. Earth BioGenome project: sequencing life for the future of life. Proc. Natl Acad. Sci. USA 115, 4325–4333 (2018).

    Article  ADS  CAS  PubMed  PubMed Central  Google Scholar 

62. Seplyarskiy, V. et al. A mutation rate model at the basepair resolution identifies the mutagenic effect of polymerase III transcription. Nat. Genet. 55, 2235–2242 (2023).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

63. Ye, C., Benegas, G., Albors, C., Li, J. C. & Song, Y. S. GPN-Star model source code. Zenodo https://doi.org/10.5281/zenodo.21501177 (2026).

64. Verbeek, M. M. et al. Mutations in the cyclic adenosine monophosphate response element of the tyrosine hydroxylase gene. Ann. Neurol. 62, 422–426 (2007).

    Article  CAS  PubMed  Google Scholar 

65. Kircher, M. et al. Saturation mutagenesis of twenty disease-associated regulatory elements at single base-pair resolution. Nat. Commun. 10, 3583 (2019).

    Article  ADS  PubMed  PubMed Central  Google Scholar 

66. Bennett, M. K., Ngo, T. T., Athanikar, J. N., Rosenfeld, J. M. & Osborne, T. F. Co-stimulation of promoter for low density lipoprotein receptor gene by sterol regulatory element-binding protein and Sp1 is specifically disrupted by the yin yang 1 protein. J. Biol. Chem. 274, 13025–13032 (1999).

    Article  CAS  PubMed  Google Scholar 

Download references

Author notes

1. These authors contributed equally: Chengzhong Ye, Gonzalo Benegas

Authors and Affiliations

1. Department of Statistics, University of California, Berkeley, CA, USA

   Chengzhong Ye & Yun S. Song

2. Computer Science Division, University of California, Berkeley, CA, USA

   Gonzalo Benegas, Carlos Albors, Jianan Canal Li, Sebastian Prillo & Yun S. Song

3. The Jackson Laboratory, Bar Harbor, ME, USA

   Peter D. Fields

4. Division of Computational Genomics and Systems Genetics, German Cancer Research Center (DKFZ), Heidelberg, Germany

   Brian Clarke

5. Center for Computational Biology, University of California, Berkeley, CA, USA

   Yun S. Song

6. Innovative Genomics Institute, University of California, Berkeley, CA, USA

   Yun S. Song

Contributions

C.Y., G.B. and Y.S.S. conceptualized the study. C.Y. developed and implemented the GPN-Star model. C.Y. and G.B. designed benchmarks, tested the method and analysed data, with contributions from C.A., J.C.L., S.P., P.D.F., B.C. and Y.S.S. P.D.F. contributed to collecting benchmarking data for M. musculus. B.C. performed the DeepRVAT analysis incorporating GPN-Star scores. Y.S.S. supervised the project. C.Y. and G.B. wrote the initial draft of the paper and Y.S.S. edited it. All authors reviewed the paper.

Corresponding author

Correspondence to Yun S. Song.

Competing interests

The authors declare no competing interests.

Peer review information

Nature thanks Miquel Anglada-Girotto, Mafalda Dias, Thomas Pierrot and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available.

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

Reprints and permissions