Introduction

The importance of big data science has led to the accumulation of huge amounts of data, recently, a generative AI model with large language model (LLM) using various types of big data has further been emerging across various fields, including healthcare, omics research, and industry^1,2. In response to the continuous emergence of new SARS-CoV-2 variants, our study introduces an innovative model and platform that integrate big data analysis, protein structure prediction using AlphaFold, and AI learning to predict highly transmissible novel SARS-CoV-2 variants. The validation of predictions made using AI is crucial, so we present the proofs through in-silico and in-vitro experiments for highly infectious variants.

The COVID-19 pandemic-related deaths have decreased, but the infectivity of the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) persists. This may be because of neutralizing antibodies from previous infections or numerous mutations in the virus, from the Alpha variant to current Omicron sublineages^3,4.

Genomic databases, such as Nextstrain and the Global Initiative on Sharing All Influenza Data (GISAID), Our World in Data, containing epidemiological data, have been instrumental in collecting valuable data during the COVID-19 pandemic^5-7. Reports on the current global efforts of genomic surveillance strategies and sequencing show varied levels of data accumulation⁸. Nonetheless, the increasing amount of available SARS-CoV-2 data has made various analytical approaches possible.

SARS-CoV-2 mutations can alter gene function. For example, the SARS-CoV-2 furin gene plays an important role in the cleavage of the spike protein, and mutations in the gene significantly affect SARS-CoV-2 fusion with the host cell membrane^9,10. Mutations in the genes have a direct impact on its protein structure, influencing the pathway used by the virus to infect host cell¹¹. Mutations in the spike protein leads to differences in the protein structure and binding affinity to ACE2, the receptor through which the virus penetrates the host cell^12,13.

Deep learning methods predict mutations and protein structures, and their application in research has led to enhanced genomic analysis with large datasets^14,15. This approach has proved valuable in predicting disease-related mutations and identifying genes of interest¹⁶.

Recently, deep learning has been utilized to analyze mutations in SARS-CoV-2^17,18. Prediction models, such as AlphaFold2, have been used to assess protein 3D structures¹⁹. This has led to the detection of structural changes caused by mutations in Delta (B.1.617.2) and Omicron (B.1.1.529) ²⁰. Alphafold2 has also been used to analyze protein genotypes and phenotypes in the RBD region²¹.

Prior approaches may not be intuitive without expertise in the field. Also, structure-based predictions of mutations cannot fully account for all aspects of SARS-CoV-2 infection. In this study, we analyzed SARS-CoV-2 mutations using artificial intelligence (AI) methods and leveraged large datasets to make various predictions. We discovered amino acid polarity changes and substitutions then evaluated the infectivity from significant mutations in the receptor binding motif (RBM) region. Hydrophobic properties in the amino acid sequence affects protein folding²². Coronavirus hydrophobicity has significant effects on amino acid properties and protein folding²³. Therefore, we specifically examined the effect of protein structures from hydrophobicity to hydrophilic and alkaline properties.

Compared with prior approaches, our evaluation involved a comprehensive analysis of epidemiological, and genomic data, capitalizing on the availability of large SARS-CoV-2 datasets. We extracted properties from the VOCs and each sublineage at the nucleotide and amino acid levels. We developed an evaluation model called Amino acid Property Eigen Selection Score (APESS) to calculate SARS-CoV-2 sequences. Various methods were utilized to validate our findings and we used machine learning to make further predictions.

Finally, we present our findings and evaluation models through Artificial Intelligence Analytics Toolkit for predicting virus mutations in protEin (AIVE), a web-based platform that integrates SARS-CoV-2 data, visualizes APESS score distribution, offers 3D protein predictions, and increases accessibility for researchers. Overall, our research provides a model for pandemic preparedness and the study of infectious diseases.

Results

Discovery of significant properties in the amino acid sequence

We identified consecutive hydrophilic amino acids within the SARS-CoV-2 spike protein in the following lineages: Wuhan-Hu-1, Alpha, Beta, and Omicron (BA.1, BA.2, BA.2.75, BA.4/BA.5, XBB, BQ.1). A series of hydrophilic amino acid regions were observed in the RBM region (amino acid sequence 437–508). In specific SARS-CoV-2 lineages, amino acid substitutions were observed in these regions (Fig. 1A and Supplementary Fig. 1). For example, the Delta and Omicron lineages contained a substitution of T478K of alkaline amino acids.

Figure 1.

Each position showed differences in the number of polarity changes; specifically, changes from hydrophilic to hydrophobic or positively charged amino acids. For most positions, we did not detect significant differences in polarity between the lineages investigated (Supplementary Tables 1-3).

As the virus progressed from the Wuhan-Hu-1 (40 in NN* and 31 in PP*) to Omicron (average 45.5 in NN* and average 27 in PP*), the number of polarity changes increased in NN* while decreasing in PP*. For these positions, polarity changes occurred in the lineages as they evolved to the Omicron variant (Fig. 1B). More polarity changes occurred in Delta, XXB, Omicron BQ.1, and Omicron BA.4/BA.5 than in the other VOCs.

We investigated the amino acid substitutions in the VOCs and compared them with the Wuhan-Hu-1 sequence as the reference. We discovered a two-fold increase in amino acid substitutions to lysine (K) and arginine (R) in later SARS-CoV-2 lineages. Meanwhile, glutamine (Q), phenylalanine (F), and glutamic acid (E) levels decreased by half in the VOCs. For phenylalanine (F), hydrophobic residue mutations occurred in areas irrelevant to the regions of consecutive hydrophilic amino acids (Fig. 1C and Supplementary Table 4 and Supplementary Fig. 2).

Various coronaviruses, including MERS-CoV, SARS-CoV-1, and SARS-CoV-2, did not show significant differences in polarity across positions (Supplementary Fig. 3A).

We investigated the polarity (hydrophilic, hydrophobic, alkaline, and acidic) of amino acids in Wuhan-Hu-1, Alpha, Beta, Omicron (BA.1, BA.2, BA.2.75, BA.4/BA/5, XBB, BQ.1), and the number of amino acids in the RBM sequence. Most amino acids were either hydrophilic or hydrophobic across lineages, with a slight increase in alkaline amino acid levels in the Omicron variants than in the others (Supplementary Fig. 3B).

We investigated the nucleotide sequences of the SARS-CoV-2 spike protein gene. We analyzed the transitions and transversions of 7,335,614 samples from the GISAID viral sequence. The rate of change from G to U was higher than that of U to G, and the rate from C to U was higher than that of U to C. Both C to G and G to C rates were extremely low (Fig. 2C and Supplementary Fig. 3C and Supplementary Table 5) ^24,25.

Figure 2.

Identification of properties in amino acid substitutions

We studied each mutation in the RBM region of lineages and sublineages using verified mutation data: Alpha, Beta, Delta, Omicron (BA.2.75), Omicron (BA.4/BA.5), Omicron (XBB), and Omicron (BQ.1).

From the viral sequences submitted to GISAID, we extracted data (n = 7,335,614) of the sublineages defined on the outbreak.info platform. The most noteworthy mutations were N440, L452, S477, T478, E484, F486, N501, and Y505 (Fig. 2A). Only seven individuals possessed mutations at D467 out of the 7,335,614 investigated (Supplementary Fig. 4A).

Importantly, the molecular work results and amino acid substitutions from the viral sequences were the same. Almost no mutations were found at the D467 position of the RBM sequence. In vitro experiments using a luciferase assay and viral entry experiments demonstrated that infectivity associated with D467 mutagenesis (D467P and D467I) was lower than the wild type (spike protein D614G) (Fig. 2C and Supplementary Fig. 5).

We investigated the mutation rates and amino acid substitutions in VOCs and variants under monitoring (VUMs). In the earlier Alpha, Beta, and Delta lineages, various amino acid changes occurred at multiple locations. In later lineages, we observed a reduction in the diversity of amino acid changes due to substitutions (Fig. 2B). L452R and T478K were the most numerous non-synonymous mutations while more synonymous mutations developed.

Differences in molecular cell levels that affect infectivity and severity in the evolution of SARS-CoV-2

We investigated the molecular cell level differences that affect the severity and infectivity of SARS-CoV-2 in its evolutionary process. From 1/23/2020 to 12/31/2022, we examined the number of infections and deaths using data from OWID and compared the periods where lineages occurred. Based on epidemiological data, we established three main periods: the major outbreak period of Delta, Omicron BA.5, and Omicron BQ. During the Delta period, there was an increase in the ratio of deaths to infections, for the Omicron BA.5 period, both infections and deaths increased, and finally in the Omicron BQ period, both decreased (Fig. 3A). We attributed these trends to mutations in the major lineages and further investigation of sublineages through viral sequence information revealed stabilization per lineage (Supplementary Tables 6 and 7).

Figure 3.

We investigated the epidemiological relationship between symptoms and disease severity in patients by examining mutations in the RBM region for major SARS-CoV-2 lineages and their sublineages. We found that these mutations became fixed over time (Fig. 3B and Supplementary Tables 8 and 9). In the Delta variant, a significant increase in symptomatic infections occurred for mutations L452R (odds ratio 4.346, 95% CI 2.378-8.191) and T478K (odds ratio 3.116, 95% CI 1.595-6.292). For most patients with Omicron BA.5 and BQ.1, mutations were asymptomatic. However, among patients with Omicron BQ.1 mutations, only those with K447T were symptomatic. The Delta, Omicron BA.5, and Omicron BA.1 variants more frequently produced milder outcomes with more mutations. Consequently, we showed that as SARS-CoV-2 evolves to Omicron, the number of asymptomatic patients with mild outcomes increased, with infection symptoms becoming less severe.

Since we cannot determine the severity and presence of symptoms in SARS-CoV-2 from mutations alone, we analyzed gene expression within host cells of infected patients. We utilized GSE235262 from the NCBI GEO database to analyze the gene expression of controls (uninfected), Alpha (B.1.1.7), Delta (B.1.617.2), and Omicron (B.1.529, BA.2, BA.4, and BA.5) variants during waves of COVID-19. In particular, we focused on the expression differences between Alpha (B.1.1.7) and Omicron (BA.4/BA.5), and between Delta (B.1.617.2) and Omicron (BA.4/BA.5), conducting enrichment analyses on genes that showed significant differences (Supplementary Table 10). Compared to Omicron (BA.4/BA.5), most genes in Delta (B.1.617.2) were up regulated in molecular and virus infection pathways (p<0.05). Notably, in the mTOR pathway, there was a correlation with receptor genes leading to the PI3K-AKT pathway, lysosome, ribosome biogenesis, autophagy, and lipid biosynthesis, In the virus infection pathway, endocytosis and metabolic pathways were related (Fig. 3C). Compared to Omicron (BA.4/BA.5), Delta is more associated with various pathways within the host cell, leading to inflammation, replication, and severity.

We believe that Omicron’s characteristics are not solely due to its molecular cell level features within the host cell but also stem from its evolution, affecting its binding affinity to the host’s ACE2. Therefore, we utilized pDockQ and HADDOCK to predict the binding affinity between SARS-CoV-2 and the ACE2 receptor affected by mutations^26,27. Based on the pDockQ results, with Wuhan-Hu-1 as the standard, the SARS-CoV-2 variants showed strong affinity in descending order of Delta (0.577), Omicron XBB (0.575), Beta (0.569), Omicron BQ.1 (0.564), and Omicron BA.4/BA.5 (0.560) (Fig. 3D and Supplementary Fig. 6). Compared to the binding affinity of Wuhan-Hu-1 (0.506), the lineages with higher binding affinity being associated with infectivity (Supplementary Tables 11-13 and Supplementary Fig. 7).

The results showed the following descending order of binding affinities for the ACE2 receptor: Delta, XBB, Beta, BQ.1, BA.4/BA.5, BA.1/BA.2, Alpha, and Wuhan-Hu-1. While Beta had a high affinity score, it has had a low occurrence rate since its emergence, and it did not receive a high score from our evaluation model.

To check for interspecies infection in humans, bats, and pangolins, we measured the bond affinity between the virus and host. We also measured the bond affinity of SARS-CoV-1. The bond strength of SARS-CoV-2 to the Homo sapiens ACE2 receptor is weaker than that of SARS-CoV-1²⁸. Although the bond strength with the same host is higher, it weakens with a different host. This suggests that viruses that evolve within a human host are likely to be the most infectious to humans. As SARS-CoV-2 evolved to Omicron, the binding affinity increased for the virus.

Overall, the Delta variant was found to be more associated with the host molecular pathway and severity. Meanwhile, the Omicron variant showed a higher interaction between the ACE2 and RBM region. The Omicron variant is more infectious.

Development of APESS, an evaluation scoring model and the evaluation of lineages

We discovered consecutive amino acids along with amino acid substitutions to lysine (K) and arginine (R). Based on these findings, we developed APESS, an evaluation model, to analyze viral sequences. APESS was calculated from four previous calculations: sub-clustering of protein structure (SCPS), polarity change score (PCS), mutation rate (MR), and biochemical properties eigen score (BPES) (Fig. 4A and Supplementary Tables 14-17 and Supplementary Fig. 8). APESS evaluates the infectivity of the SARS-CoV-2 lineage. Lower APESS scores indicated that the lineage is less infectious, whereas higher APESS scores indicated that it may be more infectious.

Figure 4.

We calculated the APESS scores of 7,335,614 viral sequences encompassing VOCs and their sublineages (Fig. 4B). The Alpha and Beta variants received lower APESS scores, whereas the Delta and Omicron variants (BA.5, BQ.1, and XBB) received higher scores.

We generated an APESS distribution graph for the VOC sublineages. Significant variants, including Delta and Omicron, showed APESS values higher than 1.62, indicating high infectivity (Fig. 4A and Supplementary Fig. 9). The APESS scores of the sublineages were similar to those of the VOCs. Some sublineages had lower APESS scores because of the wide variety of amino acid substitutions in the RBM region (Fig. 4B). The sublineages with the highest APESS scores possessed amino acid substitutions to lysine (K) and arginine (R) at S477.

Based on cumulative prevalence results from investigations on 21 March 2023 (outbreak.info), 2% cumulative prevalence of the S477K mutation was confirmed in the Omicron variants BA.5.1 and BA.5.2. Severe infections were associated with amino acid substitutions to lysine (K) and arginine (R) at the S477 position.

Utilizing a gaussian mixture model (GMM), we divided the APESS scores of 30,000 randomly selected sublineages into four components. Four groups were identified, G1 (Alpha and Beta: a centroid at 0.0289), G2 (BA.1, BA.2.75, and XBB: a centroid at 1.9474), G3 (Delta: a centroid at 1.6383). G4 (BA.2, BA.4/5 and BQ.1: a centroid at 2.0802). We also examined the distribution for all 7,335,614 lineages and found that the results were the same as the 30,000 sampled sublineages. We observed that most of Alpha and Beta were included in the component with a centroid of 0.0289, while the component with a centroid of 2.0627 included various Omicron variants.

Creation of candidate lineages and protein structure predictions

Candidate lineages were created by performing amino acid substitutions at lysine (K) and arginine (R) to evaluate their association with infectivity. APESS scores of these sequences were higher than of the SARS-CoV-2 variants. In particular, S477K and S477R amino acid substitutions were linked to increased infectivity. In contrast, amino acid substitutions of asparagine (N), serine (S), tyrosine (Y), and glycine (G) resulted in APESS scores similar to or lower than of the SARS-CoV-2 variants (Fig. 5A).

Figure 5.

Of the mutated sequences, amino acid substitutions of K and R showed the highest APESS and pDockQ scores (Fig. 5B).

We used deep learning and machine learning methods to determine the probability of amino acid substitutions at specific locations in the lineages before and after Omicron variant emergence. Models such as random forest, light-gradient boosting machine (LightGBM), extreme Gradient Boosting (XGBoost), and ensemble methods were used, and consistent results were obtained (Supplementary Fig. 10 and Supplementary Tables 18-22). Specifically, for N460K, we found that the probability increased as the lineage progressed to Omicron. For Q493R, the probability was relatively constant for Omicron lineages (Fig. 5C and Supplementary Table 23).

To validate the binding kinetics of Q493R and N460K variant RBDs to human ACE2 at the macromolecular level, we performed surface plasmon resonance (SPR). In line with our computational prediction, both variants show approximately a three-fold increase in binding affinity compared to the wild-type RBD (Supplementary Figs. 11 and 12).

We further verified the infectivity of the predicted mutations using luciferase and viral entry assays (Supplementary Fig. 5). N437R mutation led to a two-fold increase in luciferase activity and viral entry compared to the wild-type (spike protein D614G). In contrast, N460K and Q493R mutations led to a decrease in luciferase activity and viral entry compared to the wild type (spike protein D614G) (Fig. 5D and Supplementary Fig. 11). However, viral entry of N460K was 10-fold higher than that of Q493R. N460K mutations occurred at a low rate but had a higher probability of infection compared to other mutations.

Therefore, our results imply that the maintenance of the proliferative rate of virus was due to new mutations, high docking scores, and an increased presence of Omicron variants (Supplementary Figs. 13-22).

Based on our findings, we developed AIVE (https://ai-ve.org/), a web-based platform enabling rapid and precise prediction of protein structure and calculation of SARS-CoV-2 infectivity. AIVE facilitates the analysis of virus sequences entered by users and provides visualization and analysis report, allowing it to be used in environments without GPU installation (Fig. 6). For example, we used a customized sequence wherein N460K substitution was introduced into the Wuhan-Hu-1 sequence (Wuhan-Hu-1+N460K) as an input. The output generated the following four results:

Figure 6.

First, protein structure prediction results showed that Wuhan-Hu-1+N460K had increased folding compared to Wuhan-Hu-1. Meanwhile, XBB, which contains N460K substitution, demonstrated a more stabilized protein structure compared to Wuhan-Hu-1+N460K due to the alpha-helix in the 493-495 region and the beta-sheet in the 465-467 and 505-507 regions (Panel ① of Fig. 6. Second, for polarity changes, there was a decrease in consecutive hydrophilic amino acids in Wuhan-Hu-1+N460K due to the mutation (Panel ② of Fig. 6). Third, calculations of the APESS sub-scores showed a difference in scores for XBB compared to Wuhan Hu-1+N460K owing to mutations (Panel ③ of Fig. 6). Fourth, the APESS distribution was displayed as “apess”, which reflects the score at each position in the amino acid sequence. The apess score at N460K for both Wuhan-Hu-1+N460K and XBB (with N460K) was 0.112 (Panel ④-1 of Fig. 6). Furthermore, the distribution of the total “APESS” score, which comprehensively evaluates each position of the sequencing results, was used to determine the infectivity sections. Usage guidelines of AIVE can be found on supplementary information (Supplementary Fig. 23).

An apess score at N460K for XBB (with N460K) was 1.967. XBB (with N460K) mutation was located in the infectivity region with a statistically significant association (> 95%, pink color), whereas APESS score increased for Wuhan Hu-1+N460K compared to Wuhan Hu-1 but was not included in the high infectivity section containing XBB (Panel ④-2 of Fig. 6). Therefore, our research findings enable rapid protein structure prediction and APESS via meticulous and systematic structuring of data. Visualization of these analyses and the evaluation of infectivity are available on AIVE.

Discussion

In this study, we aimed to comprehensively analyze the viral sequence and protein structure of SARS-CoV-2 and investigated its association with epidemiological data. Specifically, we analyzed the impact of these factors on infection and predicted the occurrence of new mutations. Our approach involved a multistep analysis. First, we identified specific amino acid substitutions within the viral genome, focusing on their potential impact on SARS-CoV-2 protein structure. Second, we explored the impact of these structural changes on the interaction between SARS-CoV-2 and the host. Third, we assessed the potential effects of these changes on virus infectivity. This systematic analysis allowed us to gain useful insights into the behavior and evolution of the virus.

Hydrophobic protein structure plays an important role in protein stability and folding^29,30. This affects structural changes from hydrophilic to hydrophobic or positively charged residues. Flaviviruses, including Zika virus, Japanese encephalitis virus, and Dengue virus, are enveloped viruses that use envelope proteins to infiltrate the host. Infection is facilitated by E-protein folding, which causes the virus to fuse with its host³¹. Recently, SARS-CoV-2 treatments targeted the folding of the nonstructural protein NSP3 of the virus³². Research on avian coronaviruses has shown that the introduction of a hydrophobic domain into the infectious bronchitis virus E protein affects the infectivity of the virus³³. Amino acid substitutions such as L411F, F472S, D510S, and I529T have been reported in MERS-CoV^34,35. In addition to D510S, consecutive hydrophilic amino acids were observed in MERS-CoV, which we believe contributed to its low infectivity (Supplementary Fig. 24).

Our analysis primarily focused on the RBM (amino acids 437–508) region of the spike protein of SARS-CoV-2 that directly interacts with the ACE2 receptor, allowing the virus to infiltrate host cells and discovered mutations in VOCs and VUMs. In the Alpha, Beta, and Delta lineages, specific amino acid substitutions occurred in various locations. For Omicron, less diversity was observed for amino acid substitutions (Figs. 2B and 3A).

Amino acid substitutions became fixed as SARS-CoV-2 lineages progressed. Examination of epidemiological data revealed four distinct periods. During the first and third periods, there was less diversity in amino acid substitutions, leading to a decrease in both infections and deaths (Fig. 3A). In contrast, during the second period, there was more diversity in amino acid substitutions for the sublineages, leading to an increased number of cases and deaths. Thus, we established an association between the fixation of SARS-CoV-2 amino acid substitutions and infectivity.

The N437R mutation led to increased viral infectivity and received a high APESS score but was barely observed in patients. We detected only three instances of N437R amino acid substitutions, all with the AAT codon, where two or more codon positions must undergo alterations. Such changes were only in V445P, a mutation found in XBB and its sublineages, but it occurred at a miniscule rate of 0.14%. Although not yet prominent, N437R is expected to be associated with high infectivity and remains a primary candidate for close monitoring in future.

Frequency prediction of N460K through machine learning showed a 9-fold increase compared to BA.1. And the mutation showed decreased viral entry in vitro and low binding affinity in silico. The N460K mutation affects docking but combination with other mutations could increase infectivity. Ito et al. confirmed that N460K mutation forms a hydrogen bond to N-linked glycan on N90 of human ACE2³⁶.

The D467 position plays a key role in a salt-bridge interaction within the Delta variant³⁷. Due to this and mutations at the position being extremely rare suggests that D467 contributes to the structural stability of the virus. Specifically, amino acid substitutions D467P and D467I disrupt the salt bridge formed by R457 and R454 (Supplementary Fig. 4). Salt bridges play an important role in thermostabilization³⁸ and represent electrostatic interactions between positively charged and negatively charged residues. Alterations in positively charged residues can disrupt these salt-bridge interactions with negatively charged residues of SARS-CoV-2. Overall, mutations at this site impair the ability of SARS-CoV-2 to effectively infect the host. In silico results showed amino acid substitutions at D467 caused the alpha helix to present in the wild type to change into a linear structure (Supplementary Fig. 25).

XBB predicted by ML exhibited a decrease in frequency, displaying a low slope and signifying slow disappearance (Fig. 5C). Protein measurement proved to be difficult in vitro. Notably, there was an approximate 2-fold increase in binding affinity against the wild type in SPR analysis. For the Q493R mutation itself, the docking is low but combined with other mutations, it showed increased docking. We believe that the mutation itself negatively affects the spike protein structure but is important in binding with the host overall. Recent structural research which identified RBD-ACE2 complex structures of BQ1.1 (N460K) and B1.1.529 (Q493R) supports our results³⁹.

Our research shows that consecutive hydrophobic amino acids in the RBM region and specific amino acid substitutions affect not only infection but also protein structure. We assessed the infectivity of SARS-CoV-2 lineages using protein structure prediction. We observed that the emergence of various mutations corresponded to changes in binding affinity. As SARS-CoV-2 lineages progressed to Omicron, their binding affinity became stronger, resulting in increased infectivity, notably in the Delta and Omicron variants. Conversely, the Beta variant showed high binding affinity but not high infectivity. Our assumption is that for consecutive amino acids and amino acid substitutions, lysine (K) and arginine (R) have less weight; therefore, Beta is unlikely to develop into a significant variant. We speculate that Beta variant emerged early in the pandemic, potentially limiting its ability to infect a large population, while widespread vaccination efforts may have contributed to decreased infectivity.

During the evolutionary process of SARS-CoV-2, the increase in infectivity and severity of variants from Alpha to Delta was attributed to abnormalities in the molecular signaling systems of the host. In contrast, Omicron showed lower severity but increased infectivity. We believe this is because of structural changes that allow the variant to bind well with the ACE2 receptor (Fig. 3 and Supplementary Fig. 7). In addition, epidemiological data reveals a decrease in deaths, reproduction rate is maintained, and Omicron has high docking which confirms the high infectivity. We created mutagenesis sequences created from the backbone of Wuhan-Hu-1 and VOCs with amino acid substitutions to lysine (K) and arginine (R). Sequences from the backbone of Omicron XBB were predicted to have high infectivity and still liable to occur. We can therefore expect consistent infections caused by SARS-CoV-2 going forward.

We comprehensively evaluated genomic and epidemiological data. This multifaceted approach, along with the usage of APESS, machine learning, in vitro assays, and AIVE ensured a more comprehensive understanding of SARS-CoV-2 behavior and evolution. In addition, by applying GMM to the APESS score of a new sample, we can predict its membership within a specific component (Fig. 4A). This prediction allows us to make informed conjectures about its potential strains and infectivity.

Based on our findings, we predicted ∼440 mutations with a probability of occurring in the near future. One year after generating the predicted datasets, BA.2.86, EG.5.1, and HK.3 were included in the datasets were reported to the highest frequencies with 43%, 28%, and 12% at Nextstrain (https://nextstrain.org/ncov/gisaid/global/6m; 21 December 2023). APESS values of them also were high at 2.051, 1.982, and 1.986, respectively. These mutations were provided in AIVE and can be considered as important targets to prevent a new wave of infections.

AIVE features these analyses and can be easily utilized without expert knowledge of its algorithm and the analysis features can be used without GPU setup or the need for software library installations. Furthermore, AIVE offers a user-friendly interface with the flexibility to input and append experiment data for analysis. Additionally, ongoing and past analyses can be accessed after processing, providing many options for management. The time required for analysis in the local server or Google Colab can be verified. The local server run on RTX8000 takes 50 minutes, offers a large database, and is free to use (Fig. 6). After completing the job, user can examine and manage the analysis results permanently with their account generated on AIVE.

Prior research of SARS-CoV-2 included analysis of epidemiological data, molecular work, and AI applications. We adopted a comprehensive approach, utilizing real-world data and multi-faceted validation. Our research reveals that SARS-CoV-2 increased in infectivity over time, illuminating significant trends in viral infections. Our discoveries, evaluation model, and the AIVE platform will serve as the foundation for preparedness against further developments in future pandemics.

Materials and Methods

Data collection

We used the World Health Organization tracking of SARS-CoV-2 variants to verify the currently circulating VOCs. In the tracking SARS-CoV-2 variants page, variants are divided into currently circulating variants of concern (VOCs), currently circulating variants of interest (VOIs), and currently circulating variants under monitoring (VUMs). In our study, the primary currently circulating VOCs were Omicron, including BA.1, BA.2, BA.3, BA.4, and BA.5 (https://www.who.int/activities/tracking-SARS-CoV-2-variants).

To obtain the precise lineage and mutation data, cov-lineages.org (cov-lineages.org) and outbreak.info (outbreak.info) databases were used. For the sublineage analysis of the viral sequences collected through GISAID, we used accurately defined sublineage information from the outbreak.info platform.

Using WHO VOCs as standards, Alpha (B.1.1.7), Beta (B.1.351), Delta (B.1.617.2), and Omicron (BA.1, BA.2, and BA.4/BA.5) were selected as the main lineages. From WHO’s currently circulating VUMs as of 17 May 2023, we included Omicron BA.2.75 nicknamed “Centaurus,” BQ.1, and XBB.

We downloaded the mutation data of the currently reported 7,335,614 viral sequences from GISAID and utilized the data for the following purposes: first, the epidemiological data were categorized based on the patient status column: symptomatic or asymptomatic, and severity: mild or severe (including hospitalized and deceased). Second, for VOCs, including sublineages, we categorized lineage, date, and mutations in the RBM region. Third, to utilize the learning data, we categorized clades, lineages, and mutations in the RBM region. Our epidemiological data on COVID-19 used our world’s data (OWID, https://ourworldindata.org/). For confirmed cases and deaths, we used 7-day rolling average worldwide data. For vaccine doses, cumulative vaccine doses were used.

Detection of amino acid features from SARS-CoV-2 sequences

We defined four properties: hydrophilicity (polar, P), hydrophobicity (non-polar, N), alkalinity (basic, B), and acidity (A). We divided two consecutive positions of amino acids with these properties into NN*, NP*, NA*, NB*, PN*, PP*, PA*, PB*, AN*, AP*, AA*, AB*, BN*, BP*, BA*, and BB*, with the most common combination being NN*, NP*, PN*, and PP*.

We verified consecutive amino acids in the SARS-CoV-2 viral sequence spike protein. First, we categorized alanine (A), valine (V), leucine (L), glycine (G), isoleucine (I), methionine (M), tryptophan (W), phenylalanine (F), and proline (P) as hydrophobic amino acids. Second, we categorized serine (S), cysteine (C), asparagine (N), glutamine (Q), threonine (T), and tyrosine (Y) as hydrophilic amino acids. Third, lysine (K), arginine (R), and histidine (H) were categorized as alkaline amino acids (B: positively charged). Finally, aspartic acid (D) and glutamic acid were categorized as acidic amino acids.

We performed the following analysis to observe the polarity changes due to specific amino acid substitutions. Based on polarity and considering the possibility of mutations occurring, we classified each adjacent amino acid for each mutation and divided them into 64 patterns of three amino acids grouped together. We determined the pattern changes caused by the mutations. The analysis methods are available publicly in the web-based platform AIVE (https://ai–ve.org) and GitHub (https://github.com/Honglab-Research/AIVE).

Protein 3D structure evaluation

We utilized Alphafold2 to measure protein structure folding and used pDockQ and HADDock to measure the binding affinity between SARS-CoV-2 and the ACE2 receptor⁴⁰. We used AlphaFold2 to predict changes in the 3D structure of the RBM at amino acid positions 437–508 in the SARS-CoV-2 viral sequence. We performed predictions for Wuhan-Hu-1, Alpha, Beta, Delta, and Omicron (BA.1, BA.2, BA2.75, BA.4, BA.5, XBB, and BQ.1).

Through UniProt (P0DTC2 · SPIKE_SARS2), we used the viral sequence in the RBM region ‘NSNNLDSKVGGNYNYLYRLFRKSNLKPFERDISTEIYQAGSTPCNGVEGFNCYFPLQSYGFQP TNGVGYQPY.’ We performed an analysis using mutation data per lineage (https://www.uniprot.org/uniprotkb/P0DTC2/entry).

We used AlphaFold2 and through structure predictions, we observed structural similarity, distance and folding between amino acids. The predicted local distance difference test (pLDDT), Predicted Aligned Error (PAE), and AlphaFold2 results were used for analysis. Additionally, protein data bank (PDB) files and another AlphaFold2 result, were used to visualize the 3D structures using ChimeraX⁴¹. We measured the binding affinity between the human ACE2 receptor and the SARS-CoV-2 RBM region. We obtained the human ACE2 receptor sequence from UniProt (UniProt: Q9BYF1).

To increase the confidence of binding affinity analysis, first, we used pDockQ to score binding affinity. Second, we used HADDOCK 2.2 for analysis and comprehensively evaluated the HADDOCK score, root mean square deviation (RMSD), van der Waals energy, electrostatic energy, and desolvation energy (Supplementary Fig. 7). AIVE can predict the structural features of SARS-CoV-2 in different regions and the binding affinity between the human ACE2 receptor and the RBM region of SARS-CoV-2.

Analysis of the association between epidemiology and mutations

From 7,335,614 viral sequences sourced from GISAID, we identified mutations commonly discovered in VOCs and VUMs: L440K, K444T, L452R, N460K, S477N, T478K, E484A, F486V, Q498R, and Y501H. We measured the odds ratios for mutations commonly discovered in VOCs and VUMs: L440K, K444T, L452R, N460K, S477N, T478K, E484A, F486V, Q498R, and Y501H, and epidemiological indicators (symptomatic/asymptomatic and mild/severe). These analyses were performed using the traditional statistical method of logistic regression.

In the patient status column of the data file, we used asymptomatic and symptomatic rows. Deceased, hospitalized, and severe patients were categorized as severe, while outpatients, not hospitalized, and mild were categorized as mild. R version 4.1.1 and python 3.9.13 were used in our analysis.

Evaluation of SARS-CoV-2 mutation infectivity through evaluation modeling

We developed an evaluation model to comprehensively evaluate SARS-CoV-2 at the RNA, amino acid, and protein structure levels. SCPS, PCS, MR, and BPES were designed and calculated for each variant.

Sub-Clustering of Protein Structure (SCPS)

To better understand the molecular structure of SARS-CoV-2 lineages, we used the PAE, pLDDT, pdb, and hydrogen bonds that caused abnormal physical and chemical properties of the compounds to proceed with clustering.

Using clustering, the SARS-CoV-2 lineages were divided into six clusters and the main variants were mostly observed in the first and third clusters. For the cluster distribution for each lineage, “cluster 2” showed the largest difference in folding compared to the Wuhan-Hu-1 variant, and the lowest difference in folding compared to the Delta and Omicron (BA.4/BA.5) variants. Meanwhile, “cluster 4” contained the lineages with the most extreme folding showed a relaxation of folding in the main variants.

Based on this analysis, protein structure prediction-based clustering was defined as follows:

where T1 is cluster 1 (437–449), T2 is cluster 2 (450-460), T3 is cluster 3 (461-471), T4 is cluster 4 (472-479, 501-508), T5 is cluster 5 (480-489), and T6 is cluster 6 (490-500), which are positioned in the RBM of the RBD.

Omicron and Delta variants take up 88.3% of the 7,335,614 viral sequences from GISAID. Therefore, we applied a weight of 0.9 to T2 and T4 clusters and 0.1 to the remaining clusters.

Polarity Change Score (PCS)

To calculate the amino acid substitutions in the RBM sequence, we counted each lineage: Wuhan-Hu-1, Alpha, Beta, Delta, and Omicron (BA.1, BA.2, BA2.75, BA.4, BA.5, BQ.1, and XBB). For the counting method, we started from the beginning of the RBM, sliding three positions after each count until the final 72nd position was reached.

Three consecutive amino acids with two hydrophilic amino acids in Wuhan-Hu-1 changes due to mutation (0.9), when the inverse occurs (0.5), and for mutations unrelated to consecutive hydrophilic amino acids (0.1).

Mutation Rate (MR)

Based on the base trans rate and codon compositional bias, we investigated the ratio of A, T, G, and C in SARS-CoV-2 spike protein nucleotides and calculated the rate of change of bases depending on the mutations occurring. On average, the SARS-CoV-2 virus comprises 30.2% A, 19.9% G, 32.4% U, 17.6% C. For 17 lineages {Alpha (B.1.1.7), Beta (B.1.351), Gamma (P.1), Delta (B.1.617.2), Zeta (P.2), Epsilon (B.1.427), Epsilon (B.1.429), Eta (B.1.525), Lambda (C.37), Mu (B.1.621), Omicron (BA.1), Omicron (BA.2), Omicron (BA.4), Omicron (BA.5), Omicron (BA.2.75), Omicron (BQ.1), and Omicron (XBB)}, we checked the differences between the codon positions where mutations occurred. We set the 1st codon at 28%, 2nd codon at 49%, and the 3rd codon at 22% to create an amino acid mutation property equation based on biological properties.

For each codon position, MR_i (w,m) calculates the frequency of mutation (m) against the reference (w). mr(x) evaluates the frequency change for each nucleotide in a codon. N_j(x) is the frequency for j-the nucleotide. C_j(x) is the frequency for j-the codon position: C₁ (1st codon 28%), C₂ (2nd codon 49%), and C₃ (3rd codon 22%).

The total mutation rate interval is [0.0072, 0.0201] for 64 codons. 16 amino acids (G, E, A, V, R, K, N, T, I, Q, H, L, W, C, Y, and S) are amino acid substitutions that occur in SARS-CoV-2. The mutation rate interval for these amino acid substitutions is [0.0115, 0.0137].

Biochemical Properties Eigen Score (BPES)

We used amino acid properties such as residue, pH, and hydrophobic properties using principal component analysis (PCA). From VOCs, we analyzed the 72 amino acids included in the RBM region. Regarding amino acid properties, amino acids with residues C, D, E, K, R, and Y were in the order R:12.48, K:10.53, and Y:10.07. pH was in the order of R:10.76, K:9.74, and H:7.59. The hydrophobic amino acids were in the order of R: -4.5 and K: -3.9, where K and R showed similarly high scores. K and R received high scores for all the three amino acid properties (Supplementary Table 24). Based on the three biochemical properties of amino acids (residue, pH, and hydrophobicity) and the PCA clustering, we defined the eigen score as follows.

For each amino acid property, BPES_i (w,m) calculates the value of the mutation (m) against the reference (w). For bpes(x), R is the eigen vector of a residue score, H is the eigen vector of a hydrophobic score, P is the eigen vector of pH score, || || is norm of the vector, x is the amino acid, and ꞏ is the dot product.

Amino acid Property Eigen Selection Score (APESS)

APESS is the comprehensive value of SCPS, PCS, MR, and BPES. n is the number of positions, w means wild, and m means mutated.

The Kernel Density Estimator (kdeplot) from Python’s seaborn package was used to create a distribution function for APESS scores calculated across randomly sampled 30,000 SARS-CoV-2 sublineages sourced from GISAID. For all 7,335,614 lineages, we also calculated the scores and studied the distribution. This distribution function revealed that the SARS-CoV-2 lineages could be differentiated and characterized based on their APESS scores (Fig. 5A). Based on this observation, we applied Gaussian mixture models (GMM) to identify unique components of APESS score. We divided the range of the APESS score into two regions: the lower region for scores less than or equal to 3.2 and the upper region for scores greater than 3.2. To select the optimal number of components for GMM in each region, we used Bayesian Information Criterion (BIC), resulting in 4 and 5 as the optimal number of components in the lower and upper regions, respectively. We fit GMM for each region to construct prediction models with the corresponding optimal number of components.

Analysis of gene regulation from SARS-CoV-2

We analyzed the differences in gene expression for various SARS-CoV-2 lineages when they infect the host. We downloaded GSE235262 from the Gene Expression Omnibus (GEO) to obtain data on human gene expression. Comparisons were made between controls (uninfected people), Alpha (B.1.1.7), Delta (B.1.617.2), Omicron (B.1.529), Omicron (BA.2), and Omicron (BA.4/BA.5 and its sublineages). Only the genes that showed significant expression levels (log2 (TPM+1)) were selected. We used enrichR (https://currentprotocols.onlinelibrary.wiley.com/doi/10.1002/cpz1.90) for the enrichment analysis of these selected genes then visualized the comparison between Delta and Omicron (BA.4/BA.5).

SARS-CoV-2 variant prediction and machine learning model

The model predicted the probability of mutations in the RBM region of Alpha, Beta, Delta, and Omicron (BA.1, BA.2, BA.2.75, BA.4, BA.5, BQ.1, and XBB) variants and the 24 mutagenesis variants generated using machine learning. Eighteen variants were used for probability prediction. These were N440K, V445P, G446S, N460K, S477N, T478K, E484A, F486P, F490S, Q498R, N501Y, Y505H, L452R, E484K, Q493R, G496S, F486V, and K444T). Additionally, a total of 24 VOCs (S477R, S477K, N439R, Y501R, N437R, S438R, S459R, S469R, S494R, T470R, T500R, Q493K, Q493R, T470K, T500K, S469K, S494K, Y501K, N437K, N439K, N460R, N460K, S438K, and S459K) excluding overlapping variants and variants lacking expression in data were also used for probability prediction (Supplementary Fig. 26).

The prediction was carried out using data collected from 2019-12-23 to 2023-05-17 (1,242 days), and then with six different types of data, excluding specific clades:

• The first group consists of clades before 21M (Omicron B.1.1.529) including 19A, 19B, 20A, 20B, 20C, 20D, 20E, 20F, 20G, 20H, 20I, 20J, 21A, 21B, 21C, 21D, 21E, 21F, 21G, 21H, 21I, and 21J.

• The second group consisted of clades before 21K (Omicron BA.1) and 21L (Omicron ∼BA.2), including the prior group and the addition of 21M, 21K, and 21L.

• The third group consisted of clades before 22D (Omicron BA.2.75) and 22C (Omicron BA.2.12.1), including the prior group and the addition of 22C and 22D.

• The fourth group consists of clades before 22A (Omicron BA.4) and 22B (Omicron, BA.5), including the prior group and the addition of 22A and 22B.

• The fifth group consisted of clades before 22E (Omicron BQ.1), including the prior group and the addition of 22E.

• The sixth group consisted of clades before 22F (Omicron XBB), 23A (Omicron XBB 1.5), and 23 B (Omicron XBB), including the first group and the addition of 22F, 23A, and 23 B.

Machine learning determined the clade information and presence of mutations in 24 variants. The possibility of each variant is extracted through learning. Machine learning was performed using a LightGBM, which is a highly efficient gradient boosting decision tree, XGBoost, which is a scalable tree boosting system, Random Forest, which is a random forest guided tour, and an ensemble model that combined the three previous models. To calculate the probability of mutation occurrences, our custom script and following libraries were used in Python 3.9.13. Data processing employed scikit-learn for scaling, data splitting, and accuracy assessment. In the case of machine learning, LightGBM and XGBoost utilized the ‘lightgbm’ and ‘xgboost’ libraries, respectively, to build their models. Meanwhile, the Random Forest and ensemble models utilized scikit-learn’s subfunctions to create their models. The deep learning model used TensorFlow and its lower-level library, Keras, for model development.”

The Flow of Data Preprocessing and Machine Learning

We retrieved information from SARS-CoV-2 samples collected on GISAID from December 23, 2019, to May 17, 2023, spanning 1,242 days. We filtered the data using various columns, including date, Nextstrain clade, species, substitution, AA substitution, deletion, and insertion, among others.

First, for species data, we extracted the data where the host was either “homo” or “homo sapiens”. Then, we removed the samples with date values that did not follow the ‘YYYY-MM-DD’ format. Based on the insertion, deletion, and AA substitution columns, we removed the samples where deletions, insertions, or nonsense mutations occurred. Following this, we removed columns other than date, clade, and lineage, and then eliminated samples with missing values (NaN).

For the extracted samples, we organized the amino acid mutation information from the AA substitution column for the 437-508 region of the S protein. We created columns indicating the presence (1) or absence (0) of 40 target mutations. We removed mutations not found in the GISAID data, leaving the columns representing the presence or absence of 24 mutations (N440K, V445P, G446S, N460K, S477N, T478K, E484A, F486P, F490S, Q498R, N501Y, Y505H, L452R, E484K, Q493R, G496S, F486V, K444T, S494R, T500R, Q493K, T470K, N439K, and N460R).

To prepare the filtered data for training, we transformed the values for the dates into integers, representing the number of days that had passed since the initial collection date. We then normalized these values to fall within the range of [0-1]. Using one-hot encoding, we preprocessed the Nextstrain clade information. From the preprocessed data, we created six datasets, each based on different clade groups. Finally, we used the random state function to generate the training and testing datasets used in the learning process.

The training data is used for learning as the input for five different models created. The learning model then predicts the likelihood of the occurrence of 24 target mutations (Supplementary Fig. 27).

Cell culture

HEK293T and 293T-hACE2 cells were maintained in Dulbecco’s Modified Eagle Medium (DMEM, Welgene, South Korea) supplemented with 10% fetal bovine serum (FBS, Gibco) and 1% penicillin-streptomycin (Gibco) at 37 ℃ in a humidified 5% CO₂. A HEK293T cell line stably expressing human ACE2, 293T-hACE2, was generated using a lentivirus-mediated gene transduction system under the antibiotic pressure of hygromycin.

Molecular cloning and plasmid construction

The human codon-optimized full-length SARS-CoV-2 spike (Wuhan-Hu-1 strain) gene was obtained from Sino Biological Inc. (Beijing, China). The last 19 amino acids of the SARS-CoV-2 spike were truncated to increase protein expression, and a flag tag was inserted at the C-terminal end. Spike gene mutation constructs, including D614G, N437R, N460K, D467P, D467I, and Q493R, were generated using a PlatinumTM SuperFi ll PCR Master Mix (Invitrogen, Waltham, MA, USA).

Spike protein expression test

Viral spike genes were cloned into pcDNA3.1(+), and plasmid DNA was transfected into HEK293T cells using the PEI reagent (Sigma). After 48 h, the cells were lysed with a 1% NP-40 lysis buffer (150 mM NaCl, 1% NP-40, 50 mM Tris-HCl), and the extracted proteins were analyzed by SDS-PAGE and immunoblotting. A mouse anti-FLAG antibody (Sigma) was used to detect spike protein expression.

Pseudotyped virus and viral entry

The pseudotyped virus was generated from HEK293T cells (Invitrogen) by co-transfection with human immunodeficiency virus backbone plasmids expressing firefly luciferase as described previously. We used packaging plasmids (pLP1, pLP2, and pLP/VSV-G; all from Invitrogen) and pLVX-Luc-IRES-ZsGreen1_ Cat. 632187 (Luc stands for luciferase, and IRES stands for internal ribosomal entry site). For S protein pseudotyping, the full-length cDNA of the S gene (Sino Biological Inc.) was cloned into pcDNA3 and used for transfection instead of pLP/VSV-G. A plasmid carrying the gene encoding the indicated mutation in the S protein was generated using a PlatinumTM SuperFi ll PCR Master Mix (Invitrogen, Waltham, MA, USA) based on the wild-type construct, and the point mutation was confirmed by sequencing. Viral supernatants were harvested 48 h after transfection and normalized using the Lenti-X reverse transcription-quantitative PCR (qRT-PCR) titration kit _ Cat. 631235 according to the manufacturer’s protocol. Infected 293T-hACE2 cells were lysed 48 h after infection and the efficiency of viral entry was measured by comparing the luciferase activity of pseudotyped viruses bearing the wildtype or mutant S protein. A VSV-G-pseudotyped lentivirus was used as a positive control. Relative luciferase activity in the cell lysates was measured using a luciferase assay kit (Promega).

Cloning, expression, and purification of ACE2, SARS-CoV-2 RBD, and RBD variants

The SARS-CoV-2 spike protein RBD, and its mutations (N460K and Q493R), along with the N-terminal peptidase domain of human ACE2 was expressed using the Bac-to-Bac® baculovirus system (Invitrogen). The RBD (Residues Arg319-Phe541) and the N-terminal peptidase domain of human ACE2 (Residues Ser19-Asp615) with an N-terminal gp67 signal peptide were cloned upstream of the cleavable deca-histidine tag and FLAG tag. RBD variants were introduced into the wild-type RBD construct using PCR-based site-directed mutagenesis. The resulting constructs were expressed in Spodoptera frugiperda (Sf9) insect cells and secreted into the medium after being cultured at 27°C for 72 h.

Subsequently, ACE2, RBD, and its variants were isolated from the filtered supernatant using a HisTrap Excel column (Cytiva, Marlbourough, MA, USA). The C-terminal tags were cleaved using an in-house-generated HRV3C protease. The protease and cleaved C-terminal tags were removed using Ni-NTA resin (Takara, 635662), and the proteins were further purified by size-exclusion chromatography on a HiLoad 26/600 Superdex 200 pg (Cytiva). Purification was performed in a buffer containing 20 mM Tris-HCl (pH 7.5), 150 mM NaCl, and 1 mM DTT. To maintain stability during storage, each protein sample was frozen at -80°C and supplemented with 10% glycerol. For surface plasmon resonance (SPR), proteins were thawed from frozen storage immediately before the experiment.

Binding affinity measurement

The SPR experiments were performed at room temperature using BiaCore T200 with a series of S CM5 sensor chips (Cytiva). The surfaces of the sample and reference flow cells were activated using a 1:1 mixture of NHS (N-hydroxysuccinimide) and EDC (3-(N,N-dimethylamino) propyl-N-ethylcarbodiimide). The reference flow cell was left blank. All the surfaces were blocked with ethanolamine (pH 8.0). HBS-EP + (Cytiva) was used as running buffer.

For the binding affinity assay, the purified N-terminal domain of ACE2 was diluted in 10 mM sodium acetate buffer (pH 4.0) and immobilized on the chip at approximately 800 response units. The SARS-CoV-2 RBD and RBD mutations along the gradient were observed on the chip surface. After each cycle, the sensor surface was regenerated using 10 mM glycine-HCl pH 2.5). The data were fitted to a 1:1 interaction steady-state binding model using BIAevaluation 3.1 software. Curve fitting was performed via nonlinear regression using a one-site-specific binding equation in GraphPad Prism version 8.4.0 (GraphPad Software, Boston, MA, USA).

AIVE platform

The AIVE platform (https://ai-ve.org/) provides real-time protein structure predictions and APESS score calculations. AIVE reports the analysis results and mutation data from known sequences, randomly sampled SARS-CoV-2 sequences, or custom sequences provided by users. Based on the mutations in the SARS-CoV-2 lineages, Alphafold2 was utilized to analyze the protein structures. Using the Wuhan-Hu-1 sequence as a reference, we performed a comparative analysis of the main VOCs: Alpha, Beta, Delta, and Omicron. Reports are generated using visualization tools that provide Alphafold2 outputs such as pLDDT and PAE. Monomer folding and docking were visualized and scored. Differences caused by mutations are shown in the input amino acid sequences based on polarity and protein properties. APESS can predict the infectivity of sequences through the distribution of values obtained through the GMM prediction model and using previously known lineages. The results of the user’s APESS sequence scores can be visualized. AIVE runs on a high-performance computing (HPC) system with 3 RTX8700 GPUs, 96 CPUs, and 256 GB main memory.

Acknowledgements

We thank the Global Science experimental Data hub Center (GSDC) and the Korea Research Environment Open NETwork (KREONET) service for data computing and network provided by the Korea Institute of Science and Technology Information (KISTI).

Funding

This work was supported in part by grants from the National Research Foundation of Korea (NRF) grant funded by the Korea government (Ministry of Science and ICT) (NRF-2021M3H9A2097227, NRF: 2021M3A9I2080490, NRF-2022R1A2C3008162, and RS-2023-00220840), the Korea Health Technology R&D Project through the Korea Health Industry Development Institute (KHIDI), funded by the Ministry of Health & Welfare, Republic of Korea (RS-2023-00265923), and the Basic Medical Science Facilitation Program through the Catholic Medical Center of the Catholic University of Korea funded by the Catholic Education Foundation.

Authors’ contributions

D.H conceived and designed these studies. W.J.C, D.Y.S, J.Y.L, H.J.P, D.S.C, S.P.J, K.J.Y and J.P gathered SARS-CoV-2 data and performed bioinformatic analyses. S.P.J, J.P, W.J.C, D.Y.S, E.J.M, and D.H built mathematics models and performed statistical validation. J.P, W.J.C, D.S.C, D.Y.S, (relevant team) analyzed protein 3D structures and utilized machine learning for predictions. W.J.C, D.Y.S, J.P and D.H designed and built the AIVE system. U.K, G.Y.Y, H.K, T.K, S.G, H.S.C and N.H.C performed experiments in the experiment of luciferase and estimation of virus viability. D.Y.S, J.P and D.H wrote the manuscript with input from all other authors.

Competing interests

The authors declare no competing interests.

Data Availability

All data supporting the findings reported here are available in the paper and Supplementary Information. The raw and processed data, custom scripts, and codes used in this study are deposited in the GitHub repository (https://github.com/Honglab-Research/AIVE, https://github.com/Honglab-Research/AIVE-prediction) and AIVE (https://ai-ve.org).