
==== Front
Brief Bioinform
Brief Bioinform
bib
Briefings in Bioinformatics
1467-5463
1477-4054
Oxford University Press

10.1093/bib/bbae191
bbae191
Problem Solving Protocol
AcademicSubjects/SCI01060
Genotypic–phenotypic landscape computation based on first principle and deep learning
https://orcid.org/0009-0009-0856-1869
Liu Yuexing Guangzhou Laboratory, Guangzhou, Guangdong Province 510005, China

Luo Yao National University of Singapore, 21 Lower Kent Ridge Road, 119077, Singapore

Lu Xin Guangzhou Laboratory, Guangzhou, Guangdong Province 510005, China

Gao Hao Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai 200030, China

He Ruikun Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai 200030, China

Zhang Xin Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai 200030, China

Zhang Xuguang Mengniu Institute of Nutrition Science, Shanghai 200126, China

Li Yixue Guangzhou Laboratory, Guangzhou, Guangdong Province 510005, China
Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai 200030, China
GZMU-GIBH Joint School of Life Sciences, The Guangdong-Hong Kong-Macau Joint Laboratory for Cell Fate Regulation and Diseases, Guangzhou Medical University, Guangzhou 511436, China
Key Laboratory of Systems Health Science of Zhejiang Province, School of Life Science, Hangzhou Institute for Advanced Study, University of Chinese Academy of Sciences, Hangzhou 310024, China
School of Life Sciences and Biotechnology, Shanghai Jiao Tong University, Shanghai 200240, China
Collaborative Innovation Center for Genetics and Development, Fudan University, Shanghai 200433, China
Shanghai Institute for Biomedical and Pharmaceutical Technologies, Shanghai 200032, China

Corresponding authors. Yixue Li, Guangzhou Laboratory, Guangzhou 510005, China. Tel.: 86-20-62689091; Fax: 86-20-84098081; E-mail: li_yixue@gzlab.ac.cn; Xuguang Zhang, Mengniu Institute of Nutrition Science, Shanghai 200126, China. Tel.: 86-21-54920000; Fax: 86-21-54920078; E-mail: fxgzhang@mengniu.cn; Xin Zhang, Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai 200031, China. Tel.: 86-21-54920000; Fax: 86-21-54920078; E-mail: xzhang@picb.ac.cn
Yuexing Liu, Yao Luo, Xin Lu, and Hao Gao contributed equally to this work.

5 2024
02 5 2024
02 5 2024
25 3 bbae19114 12 2023
03 3 2024
10 4 2024
© The Author(s) 2024. Published by Oxford University Press.
2024
https://creativecommons.org/licenses/by-nc/4.0/ This is an Open Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (https://creativecommons.org/licenses/by-nc/4.0/), which permits non-commercial re-use, distribution, and reproduction in any medium, provided the original work is properly cited. For commercial re-use, please contact journals.permissions@oup.com

Abstract

The relationship between genotype and fitness is fundamental to evolution, but quantitatively mapping genotypes to fitness has remained challenging. We propose the Phenotypic-Embedding theorem (P-E theorem) that bridges genotype–phenotype through an encoder–decoder deep learning framework. Inspired by this, we proposed a more general first principle for correlating genotype–phenotype, and the P-E theorem provides a computable basis for the application of first principle. As an application example of the P-E theorem, we developed the Co-attention based Transformer model to bridge Genotype and Fitness model, a Transformer-based pre-train foundation model with downstream supervised fine-tuning that can accurately simulate the neutral evolution of viruses and predict immune escape mutations. Accordingly, following the calculation path of the P-E theorem, we accurately obtained the basic reproduction number (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}) of SARS-CoV-2 from first principles, quantitatively linked immune escape to viral fitness and plotted the genotype-fitness landscape. The theoretical system we established provides a general and interpretable method to construct genotype–phenotype landscapes, providing a new paradigm for studying theoretical and computational biology.

genotype-fitness landscape
interpretability
deep learning
immune escape
SARS-CoV-2
the relative basic reproduction number (R0)
Strategic Priority Research Program of Chinese Academy of Sciences XDB38050200 R&D Program of Guangzhou Laboratory SRPG22–001 SRPG22–007
==== Body
pmcINTRODUCTION

Due to the large amount of scientific data left to us, the SARS-CoV-2 virus has become a model virus for basic research on the transmission mechanism of virus-related epidemics, despite the Covid-19 pandemic having ended. To detect the evolutionary mutations that drive virus ‘immune escape’, researchers have employed various high-throughput experimental techniques [1, 2]. However, these techniques have faced challenges in considering the possibility of epistasis and recombination across virus strains and in directly linking ‘immune escape’ mutations to virus fitness. Now, the Global Initiative on Sharing All Influenza Data (GISAID) (https://gisaid.org) has collected more than 15 million SARS-CoV-2 genome sequences with spatio-temporal information. This massive data set provides an opportunity to uncover virus epidemic–related fitness information. Such findings can then be used to provide evidence for further experimental verification.

Hie et al. have made groundbreaking strides in the field of virus sequence analysis [3] by utilizing a natural language model. Obermeyer et al. [4] developed a hierarchical Bayesian multinomial logistic regression model, PyR0, to estimate the per-lineage fitness of SARS-CoV-2. Maher et al. [5] developed statistical models incorporating epidemiological covariates to account for the effects of driver mutations and the associated fitness of different viruses.

Our work focuses on constructing a computable representation of the genotype-fitness landscape by deciphering the effects of epistasis and recombination, identifying and predicting viral mutations associated with immune escape and quantifying virus fitness. To do this, following the fitness definition in population genetics, we first developed the Phenotypic-Embedding theorem (P-E theorem), a theoretical framework that allows us theoretically to quantitatively calculate the fitness of a virus population by selecting appropriate embedding expressions for the virus genome under the Encoder-Decoder Seq2Seq framework. Then, on the basis of the P-E theorem, a more universal first principle of Genotype–Phenotype for landscape computing was proposed, and on the basis of the P-E theorem, we provide a computable basis for Genotype–Phenotype landscape computation under the deep learning framework. By calculating the transmissibility \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} of SARS-CoV-2, a macrobiological phenotype determined by the genotype of the virus, we give an example of constructing a real genotype-fitness landscape following first principles by using the P-E theorem. We then constructed a Transformer-based pre-train foundation model with downstream supervised fine-tuning (SFT) to calculate and predict virus immune escape mutations and the fitness without introducing and artificially coupling macro-epidemiological information, resulting in a precise mathematical representation of the virus genotype-fitness landscape in embedding space. In our work, the fitness, a macroscopic biological variable, can be accurately described mathematically as the expectation of a latent variable function according to a hidden state distribution. Our model can also be used as a generative model to predict the likely occurrence of ‘immune escape’ mutations. Our prospective and retrospective calculations have confirmed these results.

RESULTS

P-E theorem, first principles and constructing genotype-fitness landscape

Unraveling the evolutionary causes and consequences of the genotype-fitness landscape is a fundamental yet challenging task. The features of the genotype-fitness landscape profoundly influence the evolutionary process and the predictability of evolution. However, the impact of the genotype-fitness landscape concept on evolutionary biology has been limited due to the lack of empirical information about the topography of natural genotype-fitness landscapes and the precise and quantified relationship between genomic variation and fitness [6]. To address this issue, we propose a novel approach combining basic statistical theory and deep learning models to computationally quantify the genotype-fitness landscape of the SARS-CoV-2 virus based on its genome sequence.

Starting from the definitions of evolutionary biology [6, 7], we derive a formal expression for the fitness of a virus population, which is consistent with the mathematical expression given by Obermeyer et al. [4]. Under the Encoder-Decoder Seq2Seq framework that satisfies the universality theorem of neural networks [8], we introduce the Bayes theorem and combine it with Gaussian mixed model (GMM) expansion and Expectation-Maximum (EM) algorithm to deduce the P-E theorem under the condition of linear expansion according to the virus lineage [Box 1; Supplementary Note S1, Formulas (1)–(11)]: ‘An observable macrobiological phenotype can be calculated under the Encoder-Decoder framework if we can find a reasonable embedded representation of the related microscopic genotype’.

Box 1 The mathematical foundation and derivation of the P-E theorem.

According to the basic definition of population genetics [16], the average fitness \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overline{\omega}$\end{document} of a virus population can be expressed mathematically as follows:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overline{\omega}={\lambda}_1{\omega}_1+{\lambda}_2{\omega}_2+\dots +{\lambda}_k{\omega}_k={\sum}_{k=1}^K{\lambda}_k{\omega}_k=\boldsymbol{\lambda} \ast \boldsymbol{\omega}$\end{document}(1)
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\lambda =\left({\lambda}_1,{\lambda}_2,\cdots, {\lambda}_K\right)$\end{document}, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\lambda$\end{document}k is the occurrence frequency of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} virus lineage in the whole population, which is an observable measure related to amino acid substitution characteristics, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\boldsymbol{\omega} =\left({\omega}_1,{\omega}_2,\cdots, {\omega}_K\right)$\end{document}, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\mathrm{\omega}}_{\mathrm{k}}$\end{document} is the fitness of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} virus lineage and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\boldsymbol{x}=\left({x}_1,{x}_2,\cdots, {x}_i,\cdots, {x}_L\right)$\end{document} represents L variation sites. Thus, each \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\omega}_k$\end{document} can involve the epistatic interaction of multiple mutation sites in the input sequence [12] as follows:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\mathrm{\omega}}_k\left(\boldsymbol{x}\right)={a}^{(0)}+\sum{a}_i^{(1)}{x}_i^k+\sum{a}_{ij}^{(2)}{x}_i^k{x}_j^k+\sum{a}_{ij k}^{(3)}{x}_i^k{x}_j^k{x}_k^k+\dots +\sum{a}_{1,2,\cdots, L}^{(L)}{x}_1^k{x}_2^k\cdots{x}_L^k$\end{document}(2)
We define a function \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document} as a state probability distribution of the sequence’s variation states in real space. This variation causes changes in viral fitness. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document} can have an expanded form using a GMM.
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)={\sum}_{k=1}^K{\pi}_kN\left(\boldsymbol{x}|{\boldsymbol{\mu}}_k,{\boldsymbol{\varSigma}}_k\right)$\end{document} (3)
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\pi}_kN\left(\boldsymbol{x}|{\boldsymbol{\mu}}_k,{\boldsymbol{\varSigma}}_k\right)$\end{document} is the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} component in the mixture model, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\pi}_k$\end{document} is the mixture coefficient; \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{\mu}}_k$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{\varSigma}}_k$\end{document} are the mean and variance values of the relative fitness, and they can be obtained indirectly through the EM algorithm or deep neural network. Combining formulas (1) and (3), the formal representation of the average fitness of the virus population can be obtained as follows:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overline{\omega}=E\left(\omega \right)=\int \omega \left(\boldsymbol{x}\right)P\left(\boldsymbol{x}\right)d\boldsymbol{x}={\sum}_{k=1}^K{\pi}_k{\omega}_k$\end{document}(4)
Let \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\lambda}_k$\end{document} be the occurrence frequency of the dominant cluster in the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} virus lineage when sampling virus strains; we can redefine the (dominant) per-lineage fitness as
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\omega}_k=\frac{1}{n_k}{\sum}_{i=1}^{nk}{\omega}_i^k$\end{document}(5)
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\omega}_k$\end{document} is a mathematical expectation of the contribution of a single virus strain to virus lineage or virus population (if \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\lambda}_k={\pi}_k$\end{document}) fitness. Referring to Formulas (3) and (4),
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\overline{\omega}}_{\mid P(x)}=\int \omega \left(\boldsymbol{x}\right)P\left(\boldsymbol{x}\right)d\boldsymbol{x}$\end{document}(6)
The subscript P(·) represents the state probability distribution corresponding to virus population. As long as we find a suitable functional form of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\omega \left(\boldsymbol{x}\right)$\end{document} and a corresponding probability distribution \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document}, we can obtain the per-lineage fitness of viruses by calculating their mathematical expectations by Formula (6) and then get the average fitness of the virus population by Formula (4). We can now try to do this within the context of deep learning, despite the fact that it is typically difficult to obtain and calculate \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\omega \left(\boldsymbol{x}\right)$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document} directly from real-world data with noise. Considering the independent and identically distributed (IID) properties of state random variables of different samples in general, we let P  \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\Big(\boldsymbol{z}\left|\boldsymbol{x}\right)$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{z}\right)$\end{document} be the hidden-state probability distribution in encoder–decoder spaces under the Encoder-Decoder Seq2Seq framework and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document} and P\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left(\boldsymbol{x}|\boldsymbol{z}\right)$\end{document} be the prior and posterior probability related to the input and output of the model. The relationship among \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document}, P\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\Big(\boldsymbol{z}\left|\boldsymbol{x}\right)$\end{document}, P\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left(\boldsymbol{x}|\boldsymbol{z}\right)$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{z}\right)$\end{document} is given by Bayes’ theorem:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)P\left(\boldsymbol{z}|\boldsymbol{x}\right)=P\left(\boldsymbol{x}|\boldsymbol{z}\right)P\left(\boldsymbol{z}\right)$\end{document} (7)
The function relationship defined by Formula (1) still holds when \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document}P\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\Big(\boldsymbol{z}\left|\boldsymbol{x}\right)$\end{document} maps to P\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left(\boldsymbol{x}|\boldsymbol{z}\right)P\left(\boldsymbol{z}\right)$\end{document}. When the prior probability \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document} is unknown, obtaining the fitness of virus lineages in the encoder or decoder space is transformed into finding the hidden-state probability distributions P\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\Big(\boldsymbol{z}\left|\boldsymbol{x}\right)$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{z}\right)$\end{document}, the posterior probability distribution P\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left(\boldsymbol{x}|\boldsymbol{z}\right)$\end{document} and the counterpart \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overset{\sim }{\omega}\left(\boldsymbol{z}\right)$\end{document} of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\omega \left(\boldsymbol{x}\right)$\end{document}.
Under the Encoder-Decoder Seq2Seq framework together with Bayes’ theorem and Universality Theorem of multilayer feedforward networks, because \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\boldsymbol{x}\to{\boldsymbol{x}}^{\prime }=W\boldsymbol{z}+B$\end{document},
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\overline{\omega}}_{\mid P\left(\boldsymbol{x}\right)}=\int \omega \left(\boldsymbol{x}\right)P\left(\boldsymbol{x}\right)d\boldsymbol{x}=\int \omega \left(\boldsymbol{x}\right)P\left(\boldsymbol{x}\right)\left(P\left(\boldsymbol{x}|\boldsymbol{z}\right)/P\left(\boldsymbol{x}|\boldsymbol{z}\right)\right)d\boldsymbol{x}=\int \omega \left({\boldsymbol{x}}^{\prime}\right)P\left({\boldsymbol{x}}^{\prime}\right)d{\boldsymbol{x}}^{\prime }=\int \omega \left(W\boldsymbol{z}+B\right)P\left(W\boldsymbol{z}+B\right)\mid \frac{d{\boldsymbol{x}}^{\prime }}{d\boldsymbol{z}}\mid d\boldsymbol{z}=\int \overset{\sim }{\omega}\left(\boldsymbol{z}\right)P\left(\boldsymbol{z}\right)\mid W\mid d\boldsymbol{z}$\end{document}
and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overset{\sim }{\omega}\left(\boldsymbol{z}\right)=\omega \left(\boldsymbol{x}\left(\boldsymbol{z}\right)\right)\left(P\left(\boldsymbol{x}\left(\boldsymbol{z}\right)|\boldsymbol{z}\right)/P\right(\boldsymbol{z}\left|\boldsymbol{x}\left(\boldsymbol{z}\right)\Big)\right)$\end{document}, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\frac{dx}{dz}=W$\end{document} is the Jacobian matrix, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\mid W\mid$\end{document} is a determinant of Jacobian matrix that is the sum of combinations of products of elements in different positions of all the columns of matrix W. Considering the quasi-orthogonal property of base distribution in GMM expansion for ideal classification, let \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\lambda}_k={\left|W\right|}_k$\end{document} is a mapping of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\mid W\mid$\end{document} corresponding to the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} virus lineages, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\overset{\sim }{\omega}}_k$\end{document} is a hidden space state function corresponding to the fitness of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} virus lineage. We set \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overset{\sim }{\boldsymbol{\omega}}\left(\boldsymbol{z}\right)=\left({\overset{\sim }{\omega}}_1,{\overset{\sim }{\omega}}_2\cdots, {\overset{\sim }{\omega}}_K\right)$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overset{\sim }{\boldsymbol{\omega}}\left(\boldsymbol{z}\right)=\left({\overset{\sim }{\omega}}_1,{\overset{\sim }{\omega}}_2\cdots, {\overset{\sim }{\omega}}_K\right)$\end{document} are vectors, then \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left|W\right|\overset{\sim }{\omega}\left(\boldsymbol{z}\right)\to \boldsymbol{\lambda} \ast \overset{\sim }{\boldsymbol{\omega}}\left(\boldsymbol{z}\right)$\end{document}. Finally, we have
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\overline{\omega}}_{\mid P\left(\boldsymbol{x}\right)}=\int \omega \left(\boldsymbol{x}\right)P\left(\boldsymbol{x}\right)d\boldsymbol{x}=\int \left|W\right|\overset{\sim }{\omega}\left(\boldsymbol{z}\right)P\left(\boldsymbol{z}\right) dz=\int \left(\boldsymbol{\lambda} \ast \overset{\sim }{\boldsymbol{\omega}}\left(\boldsymbol{z}\right)\right)P\left(\boldsymbol{z}\right)d\boldsymbol{z}$\end{document}(8)
Obviously, Formula (8) is the corresponding formal representation of Formula (1) under the Encoder-Decoder Seq2Seq framework of deep learning; then, we get a simple representation:
\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\overline{\omega}}_{\mid P(x)}=\boldsymbol{\lambda} \ast{\overline{\overset{\sim }{\boldsymbol{\omega}}}}_{\mid P(z)}$\end{document}(9)
Because of the hierarchical structure of the population, lineage and cluster of the virus, and referring to Formulas (1), (3), (4), (6) and (8), and based on the basic principles of population genetics (14), the elements of the vector λ are the corresponding elements of mixture coefficient in formula (4). The element λ is an occurrence frequency of a dominant cluster in virus lineage when sampling virus strains and is a macroscopic parameter that links the hidden space with the actual space. Formula (9) gives a mathematical framework for representing phenotypes in embedded Spaces, leading to the P-E theorem: ‘An observable macro-biological phenotype can be computed under the Encoder-Decoder Seq2Seq framework if we can find a reasonable embedded representation of the related microscopic genotype’.
Genotype-fitness landscape: we can plot the genotype-fitness landscape in the embedded space based on the variation state of the viral genome sequence. Starting from P-E theorem, we can regard the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overset{\sim }{\omega}\left(\boldsymbol{z}\right)$\end{document} score as the ‘immune escape’ potential of the virus, then, as a two-dimensional hypersurface, the genotype-fitness landscape can plot in a three-dimensional mapping space with ‘\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overset{\sim }{\omega}\left(\boldsymbol{z}\right)$\end{document}’, ‘\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{z}\right)$\end{document}’ and ‘time’ as axes. The stable points on the hypersurface correspond to the genomic variation states that contribute the most to the ‘immune escape’ ability. We can obtain the fitness of the virus lineage by integrating the surface density of the specific region corresponding to the virus lineage on this hypersurface. Now, if we take the gradient of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\overset{\sim }{\omega}\left(\boldsymbol{z}\right)$\end{document} along the time and the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{z}\right)$\end{document} axis on a two-dimensional hypersurface, we get and define the immune escape force of the virus.	

The P-E theorem provides a mathematical basis to establish a deep learning model to calculate virus fitness. This means that under an Encoder-Decoder framework with reasonable embedded representations, we can use only viral genomic sequence data to accurately compute viral population fitness as an observable macroscopic biological phenotype and construct the genotype-fitness landscape of viral populations (Figure 1A).

Figure 1 The P-E theorem, the framework of the CoT2G_F model and simulating the real evolution scenario of viruses. (A) Biological phenotype calculation processes. The process has five steps: (1) first principle, (2) the P-E theorem, (3) the generalized encoder-decoder deep learning model, (4) solving biological phenotypes in hidden space using the P-E theorem and (5) mapping the genotype–phenotype landscape. (B) CoT2G_F is a Transformer-based pre-train foundation model with downstream SFT, it is a standard encoder–decoder architecture that introduces a co-attention and continuous span masking mechanism. (C) and (D) show the co-attention combined pre-training and fine-tuning steps of the model (Supplementary Figure S4).

Inspired by the P-E theorem and the formal logic of fitness calculation in the proof of the P-E theorem, we propose a more general first principle of genotype–phenotype: ‘A macroscopically observable biology variable can express as a mathematical expectation of a microscopic state function according to a state probability distribution’ [Supplementary Note S4, Formulas (35) and (35')]. These ‘first principles’ provide a mathematical basis for precisely defining biological phenotypes, allowing us to deeply understand and solve computational problems of biological phenotypes, and at the same time, the P-E theorem provides a computable basis for the application of first principles. On the basis of this work, we further deduce how to start from the first principles of universality and combine the P-E theorem to achieve a more general quantitative calculation framework of biological phenotypes (Figure 1A). With the first principles and the P-E theorem, we can then design a suitable deep computational model to calculate the biological phenotypes of interest.

Modeling of virus evolution mechanisms based on virus sequences

In recent years, natural language models have made great strides in learning the composition rules of DNA and amino acid sequences and exploring the evolution of biological species [3, 9–11]. To effectively utilize the P-E theorem for fitness calculation, we need to devise a natural language model implemented within the Encoder-Decoder framework.

Different from Hie’s work [3], we construct a pre-trained + SFT foundation model, a Co-attention-based Transformer pre-train model CoT2G-F to bridge Genotype and Fitness (Figure 1): (i) introducing co-attention and self-attention mechanisms [12] to extract the spatio-temporal correlation and long-range site interaction information within and across virus sequences. These mechanisms enable the model to extract epistatic and recombination signals that promote the occurrence of immune escape mutations. (ii) Dividing training process into two stages: pre-training, simulating the neutral evolution of the virus by randomly masking one or more continuous bases at any position of the input virus sequence, learning the ‘grammar’ rules of sequence composition, and fine-tuning, in chronological order, construct the evolutionary map of the virus sequences dynamically, perform SFT according to the spatial–temporal correlation constraints and target sequences and reproduce the role of environmental selection, bringing the model closer to the virus evolution process.

Our approach stands out from other NLP models used in the study of viral sequence evolution primarily due to our optimized SFT strategy. Specifically, we first start from the pre-training model, combine the P-E theorem to determine the mutation that contributes the most to the fitness of the virus as the target sequence and then use the time-dependent dynamic fine-tuning mode to obtain the spatio-temporal selection signal of virus evolution that may not be directly attainable through other models (Figure 1C and D).

As of November 2022, 10.2 million SARS-CoV-2 spike proteins from the GISAID were selected and subjected to strict quality control to train the CoT2G-F model (Supplementary Figures S1 and S2; Supplementary Table S1). These data were used for pre-training, with further stages of fine-tuning to verify the model’s predictive ability. The CoT2G-F model can obtain the hidden state distribution of DNA bases or amino acids, and a sequence ‘embedding’ that satisfies this hidden state distribution is guaranteed to conform to the ‘grammar’ rules even in the presence of ‘semantic’ changes. This allows for the identification and inference of ‘immune escape’ mutations from the degree of ‘semantic’ change (Figure 1B–D; Supplementary Notes S2–S4). Subsequent work will focus on using the CoT2G-F model to determine the micro-state function of a latent variable corresponding to the virus fitness in the hidden space, starting from the P-E theorem [Supplementary Note S1, Formula (11)].

Identifying and predicting ‘immune escape’ mutations

Following the constrained semantic change search (CSCS) decision rule for semantics and syntax changes proposed by Hie [3], we tested the performance of our model CoT2G-F against two related mainstream models from the framework design in three technical dimensions. Our model demonstrated good performance. As shown in Figure 2A, our model outperformed the BiLSTM and Vanilla Transformer models, and the ability to identify and predict immune escape mutations was improved. By considering the spatio-temporal dynamics of virus evolution during training processes, the model was able to identify existing and predict future ‘immune escape’ mutations. We tested this using more than 2 million SARS-CoV-2 genomic sequence data from the UK as an example (Figure 2B–D; Supplementary Figure S3).

Figure 2 Modeling of virus evolution mechanisms based on virus sequences. (A) The figure shows a comparative experiment result among our method (CoT2G-F), Vanilla Transformer and Bi-LSTM. (B) The virus lineage with a high prevalence rate has a significantly higher CSC score (Supplementary Note S4), calculated by the UK-submitted spike protein sequence data for June 2020 to February 2022. The left side of (C) visualizes the semantic changes and grammaticality output by our model CoT2G-F with the horizontal axis representing the grammaticality and the vertical axis representing the semantic changes. The bottom side of (C) visualizes the semantic changes and grammaticality output by our model CoT2G-F with the horizontal axis representing the grammaticality and the vertical axis representing the semantic changes. The upper-right corner of scatter plot (C) means that the bigger the semantic changes and the more likely the virus strain is to immune escape. The upper side of (C) shows the predicted mutation sites on the Spike protein sequence of the Omicron B.A.2 virus lineage with high semantic changes and grammatical fitness, these mutation sites are real and marked at the bottom of the right figure. (D) Take the virus sequence data collected in the UK from August to October 2022 as input to identify and reproduce emergent ‘immune escape’ mutations. (E) Input sequences collected from the UK between March 2022 and June 2022 were selected to predict immune escape mutations.

Convergent evolution is a widespread phenomenon observed throughout the evolutionary process, including the later stages of viral epidemics. It refers to the independent emergence of similar traits or adaptations in different virus lineages, often in response to similar environmental pressures or functional requirements. This evolution typically occurs under strong environmental selection constraints [13]. Results showed that our model can predict and reproduce the convergent evolutionary mutations. For instance, 9 of the 18 convergent mutation sites in the representative lineages of pre-Omicron and Omicron [13–15] were predicted with the model (Supplementary Figures S4–S6).

Notably, our model was fine-tuning trained using data from the GISAID as of 31 March 2022, demonstrating its ability to give accurate predictions up to 6 months in advance without dynamic fine-tuning. Besides, our model can also be used as a generative model to generate possible immune escape mutations for screening and early warning study (Supplementary Figure S6; Supplementary Table S2).

Deciphering the intrinsic correlation between ‘semantic’ change-related ‘immune escape’ mutations and virus fitness

The immune escape ability of a virus lineage determines its fitness or the relative basic reproduction number (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}). The problem is establishing a quantitative relationship between the ‘immune escape’ mutation and the fitness of the virus conforming to the P-E theorem under the framework of the natural language model. If evolution by natural selection is to occur, there must be a genotypic change in the virus population across generations. This change dependent on the genotype; if the genotype does not change, there is no change in fitness. Therefore, ‘semantics’ and ‘grammar’ together determine the virus’s fitness. Under the theoretical framework of deep learning, the above property can be expressed mathematically as the convolution of ‘semantic’ and ‘grammatical’ terms of the CSCS criterion [3]. When multiplied by the coefficient λ, the final mathematical representation is consistent with the P-E theorem [Box 1; Supplementary Notes S1 and S4, Formulas (9), (11) and (31)]. Based on the above reasoning, we can argue that Obermeyer et al. [4] defined virus per-lineage fitness as precisely a special case of the P-E theorem. These results suggest that the term related to the ‘semantics’ change of the model CoT2G-F and the CSCS criterion [3] may be the micro-state function of the latent variable corresponding to the virus’s fitness that we hope to obtain.

Viewing the correlation between the trends of virus epidemics and the absolute mean distribution of conditional semantic change score

The P-E theorem proves that we can construct the function of a latent variable through the Encoder-Decoder Seq2Seq model or its extension model. The mathematical expectation of this function in terms of the hidden state distribution is the corresponding macroscopic observable variable that we hope to obtain. We can extend and define a conditional semantic change score (CSC score) from the CSCS score presented by Hie et al. [3]. The CSC score, as a function of latent variables describing ‘semantic’ changes, is used to measure the ‘immune escape’ ability related to virus sequence mutation. Thinking that the absolute mean of the CSC score is an index of the immune escape ability of a virus lineage, we want to see how this index changes across different virus lineages [Supplementary Notes S2–S4, Formula (27)]. We selected 2.7 million SARS-CoV-2 spike protein sequence data from the UK to calculate the absolute mean distribution of the CSC score according to the hidden-state probability distribution (Figure 3).

Figure 3 Absolute mean distribution plot of the semantic change parameter \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${CSC}_i^k$\end{document}. Take the virus sequence data collected in the UK after March 2022. (A) Using the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\ell}_1$\end{document} norm. (B) Using the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\ell}_2$\end{document} norm. The horizontal axis in the figures is the different virus lineages according to the emerging chronological order. The vertical axis is the absolute mean of the semantic change \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${CSC}_i^k$\end{document} of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} virus lineage, which is defined as the immune escape capacity of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} virus lineage. (A) and (B) both clearly show that with the progress of the global epidemic, from the Wuhan virus lineage to the Delta lineage and then to the Omicron lineage, the immune evasion ability continues to increase and even accelerates. (A) and (B) show the same absolute mean distribution trend of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${CSC}_i^k$\end{document}, but (B) changes more smoothly.

Figure 3A and B clearly show that, from the Wuhan lineage to the Delta lineage and then to the Omicron lineage, the immune evasion ability of the emerging virus lineages continues to increase and even presents an accelerated trend. The ‘immune escape’ ability of the virus lineage is a calculable and quantifiable parameter directly related to the virus genome sequence mutation, which can be used to assess the cumulative collective effects of the immune escape mutations in an emerging virus lineage. The absolute mean of the CSC score accurately describes the changing trend of the ‘immune escape’ ability of virus lineages and thus reveals the epidemic potency of the virus. This result further suggests that the CSC score is the microstate function corresponding to the virus fitness indicated by the P-E theorem.

Expressing the fitness of a virus lineage as a convolution of ‘semantics’ and ‘grammar’ and the mathematical expectation of a semantic change function according to a hidden-state probability distribution

It has long been a goal in biology to create genotype-fitness landscapes by mapping DNA sequences to mutation combinations observed in phylogeny or evolutionary experiments. However, due to the consideration of the epistatic effects of mutations, when the mutation number is large, the model becomes complicated and difficult to calculate and verify experimentally [6, 7, 16]. Vaishnav et al. [17] argue that the complete fitness landscape defined by a fitness function maps every sequence (including mutations) in the sequence space to its associated fitness. However, no theoretical model has yet been able to provide a concise computational model and framework for the genotype-fitness landscape while fully considering epistasis and recombination potential. Under the framework of the P-E theorem, we attempt to solve this problem through our deep-learning natural language model CoT2G-F combined with CSC scores.

Based on the functional relationship between the CSC score and the latent variable of the CoT2G-F model [Supplementary Note S4, Formula (22)], we know that the CSC score is a relative measure of semantic change corresponding to the relative immune escape ability of the virus and the ‘semantic’ mutation corresponding to genomic mutation. Therefore, the CSC score measures the degree of the genomic mutation of emerging virus strains relative to wild-type virus strains. Furthermore, the CSC score and the ‘grammatical’ coincidence degree determined by the mutation state distribution—the hidden state distribution—together form a measure of the virus’s immune escape ability in the form of convolution [Supplementary Note S4, Formula (31)]. Referring to the P-E theorem and by the properties of the CSC score and the formal equivalence of Formulas (11) and (31), the CSC score should be the corresponding hidden space microstate function of the macroscopic observable variable \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}. Formula (30) is the concrete realization of the P-E theorem [Formula (11)] under the CoT2G-F model. Finally, we acquire the basic formula for calculating \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} [Supplementary Note S4, Formula (33)]. The correctness of the derived calculation formula for the per-lineage fitness is evident and verified (Figures 3A and B and 4). This result and the P-E theorem provide a mathematical basis for redefining macro-phenotypes such as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} and the resulting computability.

Figure 4 Relative fitness derived from the natural language model CoT2G-F versus date of lineage emergence. (A) This figure uses a time interval of 10 days. Referring to the Pango lineage designation and assignment, the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}, which is the fold increase in relative fitness of the virus lineages according to the Wuhan lineage, is plotted in different colors. The results for 5- and 20-day time intervals are shown in Supplementary Figure S6. All the results are almost consistent, reflecting the robustness of the model. (B) The transmissibility potential (i.e. the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\hat{R}}_0$\end{document}) of virus lineages of SARS-CoV-2 by using an observable actual number of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\hat{R}}_0$\end{document} for the Wuhan lineage of SARS-CoV-2. (C) Comparison of the R/RA calculated by the Co2TG-F with the values calculated with PyR0. (D) Comparison of the R/RA calculated by the Co2TG-F with the official results.

The genotype-fitness landscape

The calculation formula of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} has a novel meaning [Supplementary Note S4, Formulas (31)–(33)]: it is the convolution of the semantic change function and the syntactic state distribution function in the embedding space. Here, convolution means giving the average of the cumulative effects of the contribution of all sampled virus variants to the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} according to the reference sequence (we selected Wuhan-Hu-1 as a reference sequence). More importantly, the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} can be precisely rewritten as a mathematical expectation of a semantic change function according to a hidden-state probability distribution, giving the quantitative relationship between the ‘immune escape’ ability of the virus and the fitness of the virus lineage. Finally, we present a comprehensive and precise formulation of the virus genotype-fitness landscape. The consistency between these two different representations provides a solid mathematical basis and research paradigm for applying deep learning models to biological problems and gaining interpretability (Supplementary Notes S1–S4).

Studying the intrinsic relationship between statistical mechanics and deep learning has always been an intriguing subject [18–20]. We know that the core concept of statistical mechanics is: ‘The macroscopic physical observables can be characterized as the ensemble average of the corresponding microstate variable functions’. Inspired by this, the inference is naturally drawn: ‘Under the Encoder-Decoder framework, a macro-biological observable can be expressed as the mathematical expectation of a function of a latent variable according to a hidden state distribution in the decoder space’. We call this framework bridging genotype–phenotype landscapes is an important corollary of the P-E theorem and gives a novel interpretable application of deep learning theory to life science [Figure 5; Supplementary Notes S1 and S4, Formulas (11), (35) and (36)]. Finally, we can plot the genotype-fitness landscape in the embedded space based on the variation state of the viral genome sequence according to Equation (30). It is a two-dimensional hypersurface, and we can obtain the fitness of the virus lineage by integrating the surface density of the specific region corresponding to the virus lineage on this hypersurface, thereby defining the immune escape force of the virus (Figure 5B; Supplementary Note S4).

Figure 5 CoT2G-F framework for building genotype–phenotypic landscapes. (A) CoT2G-F is a natural language deep learning model that introduces co-attention and continuous span masking mechanism and takes the Transformer as a kernel. The model links the prior probability \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}\right)$\end{document} and hidden state distribution \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{z}|\boldsymbol{x}\right)$\end{document}to the prior probability \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{z}\right)$\end{document} and the posterior probability distribution \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P\left(\boldsymbol{x}|\boldsymbol{z}\right)$\end{document} in the encoder–decoder space via Bayes’ theorem, mapping viral protein sequences as latent variables reflecting semantic and grammatical changes in the embedding space. Then, refer to the semantic and grammatical changes to identify the virus sequence mutation. According to the hidden-state posterior probability distribution \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $P(z)$\end{document} in the decoder space, the model can express the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} as the mathematical expectation of a specific latent variables function related to viral sequence mutation according to a hidden-state probability distribution, building a genotype-fitness landscape. (B) Plot virus genotype-fitness landscape. As a two-dimensional hypersurface, the genotype-fitness landscape can plot in a three-dimensional mapping space with ‘semantics’, ‘syntax’ and time as axes. Starting from the sampled data, we fitted a quadric surface function and drew a sketch of the Genotype-Fitness landscape, in which only three virus lineages, including Omicron, were shown.

Inferring the R0 of the virus lineage

According to our theoretical framework, computing the fitness of each virus lineage and the fold increase in relative fitness to obtain the transmissibility potential, i.e. basic reproduction number (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}) of SARS-CoV-2 requires adequate sampling for each virus lineage in a specific region and time interval. In this region and time interval, the virus begins to emerge, sustain its spread, increase in number and eventually reach a plateau. Based on the biological definition of fitness, the selected sampling time interval needs to ensure the virus is transmitted for enough generations to reduce the violent fluctuation of the signal [4]. It is also crucial for determining the occurrence frequency of the dominant cluster of the different lineages and simultaneously calculating the absolute mean of the CSC score [Formula (30); Figures S3A and B; Supplementary Note S4]. Currently, the GISAIDwww.gisaid.org contains over 2 million genome sequences of SARS-CoV-2 from the UKhttp://gisaid.org, and the spatio-temporal distribution of the virus lineages reflects the major trend of the world epidemic of SARS-CoV-2, which can serve as a model system for a reliable and feasible G-P principle.

We computed the virus per-lineage fitness using 2.7 million SARS-CoV-2 spike proteins from the GISAIDwww.gisaid.org and adopted the Pango lineage designation and assignment. To calculate the factor \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\lambda}_k$\end{document} in the per-lineage fitness [Box 1; Supplementary Note S4, Formulas (9) and (30)], we first performed lineage assignment with Pangolin. We then selected the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} lineage and its nearest sublineages to form a population, sampling for a given time interval \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${t}_{bin}$\end{document}, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\lambda}_k$\end{document} was the occurrence frequency of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${k}_{th}$\end{document} lineage in this small population. The calculation result is shown in Figure 4, where the time interval we chose was 10 days (Supplementary Figure S6). Our model correctly inferred that the WHO classification variant Omicron lineages (Pango lineage BA.2.X) had very high relative fitness in a pandemic. It was approximately 10–12 times higher than the original Wuhan and lambda lineages {[95% confidence interval (CI) 16.78 to 16.84]; Figure 4A; Supplementary Figure S7; Supplementary Figure S8}, accurately predicting its rise in the spread regions.

In our study, we propose a simple calibration process to obtain the absolute basic reproduction number (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\hat{R}}_0$\end{document}) of virus lineages by using an observable actual number of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\hat{R}}_0$\end{document} for the Delta lineage of SARS-CoV-2. [Figure 4A; Formula (34); Supplementary Note S4]. Our calculated \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\hat{R}}_0$\end{document} align with the results obtained by Obermeyer et al. [4] (Figure 4C) and are also in agreement with the actual and official monitoring results (Figure 4D). This result validates our proposed P-E theorem, allowing for the precise computation of the relative fitness of the virus lineage and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}.

DISCUSSION

Applying AI and deep learning methods to biology has two key points. The first is how to make the model more in line with biological logic, and the second is to delve into the connotation of the model and provide a mathematical basis for interpretability. This study explores these two aspects and presents an initial research paradigm. In the first step, we propose a general theorem, the P-E theorem, which is a computational path of macrobiological phenotypes under the framework of deep learning. Then, in the second step, we constructed an encoder-decoder Seq2Seq deep learning model CoT2G-F that is more in line with the evolutionary biology scenario. The model can accurately predict the immune escape mutation of SARS-CoV-2 and calculate the immune escape ability of the virus lineage. In particular, because the CoT2G-F is a basic pre-trained foundation model for virus evolution research, the model can perform different downstream analysis tasks by various SFT.

In our work, the calculation of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} was realized by the model CoT2G-F under the guidance of the P-E theorem, which gives the first principle-based interpretability of the model. The \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} is a mathematical expectation of the ‘semantic’ change according to a hidden state distribution, which directly gives quantitative relationships among ‘semantics’ change and ‘grammatical’ fitness, ‘immune escape’ mutation and the fitness of virus lineages [Supplementary Note S4, Formulas (30) and (33); Supplementary Figure S9]. This result can be seen as the correspondence in the biology of the core concept in statistical mechanics that ‘macroscopic observable physical quantities are the ensemble average of the corresponding microscopic state variable functions’, and is a corollary of the P-E theorem, which provides a new valuable perspective and interpretability for applying deep learning to virus evolution. It reveals the intrinsic correlation between deep learning theory and statistical mechanics studied by many researchers [18–20].

The practice and results of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} calculation in the present work provide a feasible paradigm for the computational modelling of genotype–phenotype landscapes and a compelling commentary on the first principle of ‘selection imprinting recorded in viral genomic mutation’, which is a concrete and successful example of the P-E theorem. We believe that, following the research path proposed in the present work, many research examples that adhere to the P-E theorem will emerge in the future.

In our study, we employ the Transformer as a kernel in our model, incorporating co-attention, self-attention and continuous span masking mechanisms to mimic the spatial–temporal dynamical associations of virus evolution. This comprehensive approach that combines the main deep learning algorithms provides several noteworthy benefits. Firstly, it enables us to extract long-range correlation mutation information from both upstream and downstream regions of the viral sequence, thereby capturing crucial spatio-temporal features relevant to virus evolution. Secondly, our model can effectively simulate genetic drift and discerns which mutations are more likely to be advantageous and subsequently fixed through natural selection and exhibits significantly improved accuracy in the identification and prediction of potential future ‘immune escape’ mutations. The efficacy of our model is clearly demonstrated in Figure 2A–D, which provides compelling evidence of its predictive capabilities.

Importantly, we have constructed a formal three-dimensional embedding space to visualize the virus genotype-fitness landscape, utilizing the quantitative correlation between ‘immune escape’ mutations and viral lineage fitness. This pioneering approach represents the first precise mathematical representation of a virus’ genotype-fitness landscape, enabling us to understand the topological properties that drive virus evolution and immune escape. The development of this fundamental methodology establishes a revolutionary theoretical and technical framework for future research, with potential applications across diverse contexts.

The proposed P-E theorem may demonstrate its universality in biological research by effectively studying fundamental phenotypic states in biology, such as brain homeostasis, metabolic homeostasis, cell fate state, tumour microenvironment, survival capability and advanced aging state [21]. This theorem establishes a robust mathematical foundation for deterministically describing biological homeostasis and deciphering its stability mechanisms.

METHODS

Data access and preprocessing

Virial sequence and metadata were obtained from the Virus Genome Toolkits (ViGTK, https://www.biosino.org/ViGTK), which is a mirror site for SARS-CoV-2 genome data of the GISAID (https://gisaid.org) and the National Center for Biotechnology Information (NCBI), https://www.ncbi.nlm.nih.gov. We downloaded 16 395 521 genome sequences from ViGTK as of 26 November 2022. Records with missing time or location were discarded. After removing low-quality sequences of length less than 29 000 bp or with ambiguous bases more than 400, we obtained 14 959 290 genome sequences altogether (Supplementary Table S1). Lineage assignment was performed with Pangolin (version 4.1.2) [22] subsequently. To get the spike protein sequences, genome sequences were firstly aligned to SARS-CoV-2 reference genome (Wuhan-Hu-1, GeneBank accession number: NC_045512.2) [23] using the MAFFT7 (version 7.505) [24]. Then, coding sequences were extracted from the alignment and translated into protein sequences. We only considered spike protein sequences whose sequence length were between 1265 and 1275 amino acid residues. Additionally, spike proteins with ambiguous amino acids were excluded, which yielded a high-quality dataset of 10 184 273 sequences in total (Supplementary Figure S1). A total of 2 713 327 SARS-CoV-2 spike protein sequences that belong to the UK was subset from processed dataset. Tools for the data preprocessing are list in Supplementary Table S3.

Preparation of the dataset for pre-training and fine-tuning

After preprocessing, we bin time intervals into 1-month segments, resulting in 35 time-bins (Supplementary Figure S4). In each month, sequences were arranged by their evolution distance relative to Wuhan-Hu-1 strain of SARS-CoV-2 (Supplementary Figure S5). In order to find the nearest neighborhood for the sequence in t time, all sequences in t-1 month were aligned to this sequence. The spike protein sequence in t-1 month with minimum edit distance relative to Wuhan-Hu-1 strain were selected as the nearest neighborhood.

Co-attention-based transformer model bridging genotype and fitness

We built a pre-train foundation model CoT2G-F with downstream SFT. It is an Encoder-Decoder Seq2Seq architecture natural language model with a Transformer as the technical framework and introduced the co-attention and contiguous spans masking mechanisms. Pre-trained with 10.2 million SARS-CoV-2 Spike protein sequence data from GISAID, the data set with 2.7 million sequences from the UK region was selected to validate our model. Details are described in the Supplementary Text. Tools for the model building are listed in Supplementary Table S3.

Comparison with other major natural language models

Figure 2A shows a comparative experiment result with our method (CoT2G-F) and Vanilla Transformer without the introduction of the cross-attention mechanism, and Bi-LSTM, which is used in Hie [3] to capture the semantic and syntactic variations of protein sequences. It is worth noting that Hie (4) and our method (CoT2G-F) solve different problems: Hie [3] uses Bi-LSTM as an encoder to extract semantic variation and syntactic information and predict Immune escape, and our objective is to generate the protein sequence of immune escape in the next time slice, so we connected the decoder module in Bi-LSTM to generate protein sequence. In this figure, the comparison experiment is based on Bi-LSTM, Vanilla Transformer and our method (CoT2G-F). We compare the effect of the model from two dimensions; the first is the identification of sequence mutation sites. We hope that the model can accurately predict the mutated sites (i.e. high recall rate) without randomly predicting the mutated sites (i.e. high precise), so the f1 score combined with recall and precise is used to compare the identification of mutated sites by different models, denoted as ‘Site’ in this figure. Compared with Bi-LSTM and Vanilla transformer, our method (CoT2G-F) improved f1 score of mutation sites by 1.53% and 0.29%, respectively. The second dimension is to identify the amino acid mutated into the site based on the identification of the mutation site in the first dimension, denoted as ‘Amino acid’ in this figure, we find the effect of Bi-LSTM is not ideal when dealing with the task of generating long protein sequences and cannot accurately identify the mutated amino acids. However, the attention mechanism in Vanilla Transformer and our method (CoT2G-F) can capture the long-range dependencies of long protein sequences and can learn the mutation logic of amino acids relatively better. Compared with Vanilla Transformer, our method (CoT2G-F) improves the identification rate of mutated amino acids by nearly 1.64%, thanks to the information gain brought by the introduction of spatiotemporal information in our method (Figure 2A).

Model training based on the fundamental principles of evolutionary biology

CoT2G-F is an Encoder-Decoder Seq2Seq architecture natural language pre-training model (Supplementary Note S1) that introduces a co-attention and contiguous spans masking mechanisms and uses a Transformer as a basic framework. To make the model better reflect the real virus evolution scenario, we adopt a two-step training mode. The first step is the pre-training process. First, from the outbreak of the virus to the present, the time cross-section is divided into months, and the global SARS-CoV-2 genome sequences in 1 month are arranged according to the evolutionary distance to form a two-dimensional genome sequence matrix \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{X}}_{\boldsymbol{t}}=\left({\boldsymbol{x}}_{\left(t,1\right)},\cdots, {\boldsymbol{x}}_{\left(t,p\right)},\cdots, {\boldsymbol{x}}_{\left(t,N\right)}\right)$\end{document}; the row of the submatrix \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{x}}_{\boldsymbol{t},\boldsymbol{p}}$\end{document} represents the genome sequence of a virus around spatial position p and t is the corresponding time of a certain time cross-section. During pre-training, the co-attention module of CoT2G-F starts from the first row of the matrix \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{X}}_{\boldsymbol{t}}$\end{document}, and sequentially takes three spatially connected Spike protein sequences \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left({\boldsymbol{x}}_{\left(t,p-1\right)},{\boldsymbol{x}}_{\left(t,p\right)},{\boldsymbol{x}}_{\left(t,p+1\right)}\right)$\end{document}, get \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\Big({\boldsymbol{x}}_{\left(t,p-1\right)},{x}_{\left(t,p\right)}$\end{document}) and\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left({\boldsymbol{x}}_{\left(t,p\right)},{\boldsymbol{x}}_{\left(t,p+1\right)}\right)$\end{document} pairwise co-attention signals, and then the Transformer module randomly masks any amino acid sites of the Spike protein with a proportion under 15%, which is completely different from Obermeyer et al. [4]. We adopt such a pre-training process that can not only replicate the ‘neutral evolution’ source of virus variation but also consider the long-range associations between sites in the viral genome sequence composition rules, as well as possible recombination effects. Therefore, our pre-training model that reflects the ‘semantics’ and ‘grammar’ of viral sequence composition is more reasonable. After completing the first step of training, we can already calculate the ‘semantic’ changes of the Spike protein of each input virus genome sequence as defined in Supplementary Note S4 and select the target sequences according to the CSCS defined in Obermeyer et al. [4] for the second step training.

The second step of training is the dynamic fine-tuning process. The occurrence of virus ‘immune escape’ mutations are the result of selection. Pre-training can also indirectly sense the evolution of a virus over time. It is certain that if the evolution direction of the virus over time is highlighted in the training, it will be able to get a stronger ‘selection’ signal and can more accurately identify and predict ‘immune escape’ mutations. According to the time arrow, we take the three adjacent virus genome sequence matrices \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{X}}_{\boldsymbol{t}-\mathbf{1}},{\boldsymbol{X}}_{\boldsymbol{t}},{\boldsymbol{X}}_{\boldsymbol{t}+\mathbf{1}},$\end{document}the five selected sequences in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{X}}_{\boldsymbol{t}+\mathbf{1}}$\end{document} with the largest CSCS scores used as the training targets, and then the three spatially correlated sequences \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left({\boldsymbol{x}}_{\left(t,p-1\right)},{\boldsymbol{x}}_{\left(t,p\right)},{\boldsymbol{x}}_{\left(t,p+1\right)}\right)$\end{document} are selected according to the same pre-training mode in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{X}}_{\boldsymbol{t}}$\end{document}, and finally in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{X}}_{\boldsymbol{t}-\mathbf{1}}$\end{document} select a sequence \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{x}}_{\left(t-1,p\right)}$\end{document} with the closest evolution distance to the sequence \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{x}}_{\left(t,p\right)}$\end{document} in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{X}}_{\boldsymbol{t}}$\end{document}, when calculating the pairwise co-attention mechanism between \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\left({\boldsymbol{x}}_{\left(t,p-1\right)},{\boldsymbol{x}}_{\left(t,p\right)},{\boldsymbol{x}}_{\left(t,p+1\right)}\right)$\end{document}, while computing the co-attention mechanism between \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{x}}_{\left(t-1,p\right)}$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\boldsymbol{x}}_{\left(t,p\right)}$\end{document}. The sequence input of the dynamic fine-tuning process is sequentially expanded in spatio-temporal order, from the outbreak to the present. Finally, we obtained a trained model CoT2G-F that can accurately reflect the ‘mutation + selection’ mechanism during virus evolution, can predict and identify ‘immune escape’ mutations and then compute \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} (Supplementary Note S4).

‘Immune escape’ mutations prediction and R0 inference

As in Hie et al. [3], our model CoT2G-F can construct a semantic and grammatical ‘embedded’ representation for given virus genome sequences and even has significant advantages to extract the upstream and downstream long-range associated variation information of the sequence itself and obtain the epistatic effects of ‘immune escape’ mutations and recombination signals of virus strains. The model identifies and predicts virus ‘immune escape’ mutational patterns by the learned semantic and grammatical rules. We further prove that the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} or a fold increase in relative fitness of a virus lineage over the Wuhan lineage as a macro-observable variable is a mathematical expectation of a latent variable function according to a micro-hidden-state probability distribution (Supplementary Note S4).

When computing \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}, we find that this microscopic latent variable function is what we define as the CSC score. In our model, the CSC score counts the semantic change of the embedded virus sequence, reflecting the immune escape ability of the virus strain (Supplementary Note S4). We applied CoT2G-F to all publicly available SARS-CoV-2 spike protein sequences following the Pango lineage designation and assignment to identify and predict ‘immune escape’ mutations and calculate the per-lineage fitness of the virus without introducing additional epidemiological parameters and then obtain the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} of the virus lineage. The introduced conditional semantic change score \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $CSC={\left\Vert \hat{\boldsymbol{z}}-\hat{\boldsymbol{z}}\left[{\overset{\sim }{\boldsymbol{x}}}_{\boldsymbol{i}}\right]\right\Vert}_{{\mathbf{\ell}}_{\boldsymbol{x}}}$\end{document}, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\left\Vert \cdotp \right\Vert}_{\ell_{\mathrm{x}}}$\end{document} is the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${\ell}_x$\end{document} norm, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $x=1,2$\end{document}(Supplementary Note S1), in which \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\hat{\boldsymbol{z}}=\frac{1}{N}{\sum}_{i=1}^N{\hat{\boldsymbol{z}}}_{\boldsymbol{i}}$\end{document}. Considering the history of virus evolution and the reference of the relative fitness, we selected the virus sequence Wuhan-Hu-1 as a reference corresponding to the latent variable \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} $\hat{\boldsymbol{z}}$\end{document} to align and normalize the CSC score (Supplementary Note S4).

Key Points

We present a precise formulation of the fundamental principles governing the calculation of biological phenotypes. Specifically, macroscopically observable biological variables (biological phenotypes) can be mathematically represented by the expectations of their corresponding microstate functions within their microstate distribution.

We introduce and substantiate the genotypic–phenotype embedding theorem (P-E theorem). This theorem, established for the first time, demonstrates that the foundational principles of biological phenotype calculation can be accurately computed within the framework of deep learning.

Employing a deep learning model, we successfully calculate the fitness (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document}) of SARS-CoV-2, an observable macro-phenotype. We provide a concrete example of biological phenotype computation.

The deep learning model we devised can effectively simulate the genetic drift and environmental selection involved in virus evolution. Furthermore, it accurately predicts the immune escape mutations of SARS-CoV-2.

Supplementary Material

Supplementary_data_bbae191

ACKNOWLEDGEMENTS

We thank Guoping Zhao, Li Jin, Jiarui Wu, Peiji Zeng, Ruibin Liu and Zhijian Lang for their deep and helpful discussions. We thank Tao Zeng and Junwei Liu for their assistance with the manuscript. We are also particularly grateful to Brian Hie for his constructive responses to our questions about the details of his related excellent models.

AUTHOR CONTRIBUTIONS

Conceptualization: Yixue Li. Methodology: Yao Luo, Yuexing Liu, Xin Lu, Hao Gao, Ruikun He, Xin. Zhang, Xuguan Zhang and Yixue. Li. Investigation: Yao Luo, Yuexing Liu, Xin Lu, Hao Gao and Ruikun He. Resources: Yixue Li. Data curation: Yuexing Liu, Yao. Luo, Xin Lu, Hao Gao and Ruikun He. Writing, original draft: Yixue Li, Yuexing Liu, Yao Luo, Xin Lu, Hao Gao. Writing, review and editing: Yixue Li, Yuexing Liu. Supervision: Yixue Li. Funding acquisition: Yixue Li.

FUNDING

This work was supported by the Strategic Priority Research Program of Chinese Academy of Sciences [XDB38050200 to Y.L.] and the R&D Program of Guangzhou Laboratory [SRPG22–001 to Y.L., SRPG22–007 to Y.L.].

DATA AVAILABILITY

GISAID sequence data are publicly available at https://gisaid.org. We used the following publicly available datasets for model training: SARS-CoV-2 spike protein sequences from the GISAID (https://gisaid.org); training and validation datasets are deposited to https://zenodo.org/record/7388491. The data set sampled from the UK region according to the Pango lineage designation and assignment, used for the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{upgreek} \usepackage{mathrsfs} \setlength{\oddsidemargin}{-69pt} \begin{document} ${R}_0$\end{document} calculation is available at https://zenodo.org/record/7388491.

CODE AVAILABILITY

The source code for data preprocessing and modelling are available at https://github.com/martyLY/CoT2G-F.

Author Biographies

Yuexing Liu is a postdoctoral fellow at Guangzhou Laboratory, Guangzhou, China.

Yao Luo is a master student at the National University of Singapore, Singapore.

Xin Lu is a PhD student at Guangzhou Laboratory, Guangzhou, China.

Hao Gao is a software engineer at Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai, China.

Ruikun He, PhD, is a senior scientist at Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai, China.

Xin Zhang, PhD, is a technical expert at Shanghai Institute of Nutrition and Health, Chinese Academy of Sciences, Shanghai, China.

Xuguang Zhang, PhD, is a professor at Mengniu Institute of Nutrition Science, Shanghai, China.

Yixue Li, PhD, is a principal investigator at Guangzhou Laboratory, Guangzhou, China.
==== Refs
References

1. Mathew  D, Giles  JR, Baxter  AE, et al.  Deep immune profiling of COVID-19 patients reveals distinct immunotypes with therapeutic implications. Science  2020;369 :eabc8511.32669297
2. Starr  TN, Greaney  AJ, Hilton  SK, et al.  Deep mutational scanning of SARS-CoV-2 receptor binding domain reveals constraints on folding and ACE2 binding. Cell  2020;182 :1295–1310.e20.32841599
3. Hie  B, Zhong  ED, Berger  B, et al.  Learning the language of viral evolution and escape. Science  2021;371 :284–8.33446556
4. Obermeyer  F, Jankowiak  M, Barkas  N, et al.  Analysis of 6.4 million SARS-CoV-2 genomes identifies mutations associated with fitness. Science  2022;376 :1327–32.35608456
5. Maher  MC, Bartha  I, Weaver  S, et al.  Predicting the mutational drivers of future SARS-CoV-2 variants of concern. Sci Transl Med  2022;14 :eabk3445.35014856
6. de  Visser  JAGM, Krug  J. Empirical fitness landscapes and the predictability of evolution. Nat Rev Genet  2014;15 :480–90.24913663
7. Fragata  I, Blanckaert  A, Dias Louro  MA, et al.  Evolution in the light of fitness landscape theory. Trends Ecol Evol  2019;34 :69–82.30583805
8. Hornik  K, Stinchcombe  M, White  H. Multilayer feedforward networks are universal approximators. Neural Netw  1989;2 :359–66.
9. Alley  EC, Khimulya  G, Biswas  S, et al.  Unified rational protein engineering with sequence-based deep representation learning. Nat Methods  2019;16 :1315–22.31636460
10. Bepler  T, Berger  B. Learning protein sequence embeddings using information from structure. arXiv preprint 2019. arXiv:1902.08661. https://arxiv.org/abs/1902.08661.
11. Rao  R, Bhattacharya  N, Thomas  N, et al.  Evaluating protein transfer learning with TAPE. Advances in Neural Information Processing Systems  2019;32 :9689–701.33390682
12. Vaswani  A, Shazeer  N, Parmar  N, et al.  Attention is all you need. Advances in Neural Information Processing Systems  2017;30 :5998–6008.
13. Cao  Y, Jian  F, Wang  J, et al.  Imprinted SARS-CoV-2 humoral immunity induces convergent omicron RBD evolution. Nature  2022;614 :521–9.36535326
14. Focosi  D, Quiroga  R, McConnell  S, et al.  Convergent evolution in SARS-CoV-2 spike creates a variant soup from which new COVID-19 waves emerge. Int J Mol Sci  2023;24 :2264.36768588
15. Ito  J, Suzuki  R, Uriu  K, et al.  Convergent evolution of SARS-CoV-2 omicron subvariants leading to the emergence of BQ.1.1 variant. Nat Commun  2023;14 :2671.37169744
16. Kondrashov  DA, Kondrashov  FA. Topological features of rugged fitness landscapes in sequence space. Trends Genet  2015;31 :24–33.25438718
17. Vaishnav  ED, de  Boer  CG, Molinet  J, et al.  The evolution, evolvability and engineering of gene regulatory DNA. Nature  2022;603 :455–63.35264797
18. Bahri  Y, Kadmon  J, Pennington  J, et al.  Statistical mechanics of deep learning. Annual Review of Condensed Matter Physics  2020;11 :501–28.
19. Carleo  G, Cirac  I, Cranmer  K, et al.  Machine learning and the physical sciences. Rev Mod Phys  2019;91 :045002.
20. Liu  Y, Yao  X. Ensemble learning via negative correlation. Neural Netw  1999;12 :1399–404.12662623
21. Wang  F, Liu  J, Gao  F, et al.  Exploring multi-omics latent embedding spaces for characterizing tumor heterogeneity and tumoral fitness effects. bioRixv preprint, 2023. 2023.02.09.527693. 10.1101/2023.02.09.527693.
22. Rambaut  A, Holmes  EC, O’Toole  Á, et al.  A dynamic nomenclature proposal for SARS-CoV-2 lineages to assist genomic epidemiology. Nat Microbiol  2020;5 :1403–7.32669681
23. Wu  F, Zhao  S, Yu  B, et al.  A new coronavirus associated with human respiratory disease in China. Nature  2020;579 :265–9.32015508
24. Katoh  K, Standley  DM. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol  2013;30 :772–80.23329690
