ORIGINAL RESEARCH article

Front. Immunol., 14 August 2025

Sec. Viral Immunology

Volume 16 - 2025 | https://doi.org/10.3389/fimmu.2025.1637955

From immune evasion to broad in silico binding: computational optimization of SARS-CoV-2 RBD-targeting nanobody

  • 1. Department of Physics, School of Science, Tianjin University, Tianjin, China

  • 2. State Key Laboratory of Synthetic Biology, Tianjin University, Tianjin, China

  • 3. Frontiers Science Center for Synthetic Biology and Key Laboratory of Systems Bioengineering (Ministry of Education), Tianjin University, Tianjin, China

Abstract

Introduction:

The rapid evolution of SARS-CoV-2 Omicron variants highlights the urgent need for therapeutic strategies that can target viral evolution and leverage host immune recognition mechanisms. This study uses molecular dynamics (MD) simulations to analyze the immune evasion mechanisms of class 1 nanobodies against emerging SARS-CoV-2 variants, and to develop an efficient in silico pipeline for rapid affinity optimization.

Methods:

We employed MD simulations and binding free energy calculations to investigate the immune evasion mechanisms of four class 1 nanobodies (R14, DL4, VH ab6, and Nanosota9) against wild-type (WT) and Omicron variants, including BA.2, JN.1, and KP.3/XEC. Building on these findings, we established a streamlined nanobody optimization pipeline integrating high-throughput mutagenesis of complementarity-determining regions (CDRs) and hotspot residues, protein-protein docking, and MD simulations.

Results:

MD analysis confirmed that the immune evasion mechanism of KP.3/XEC is significantly associated with the Q493E mutation, which weakens electrostatic interactions between the nanobodies and the receptor binding domain (RBD). Through our pipeline, we identified high-affinity mutants including 3 for R14, 3 for DL4, 11 for VH ab6, and 9 for Nanosota9. The optimized R14 variant L29W/S52C/A101V demonstrated exceptional performance, achieving a 62.6% binding energy improvement against JN.1 (-76.88 kcal/mol compared to -47.3 kcal/mol for original R14 nanobody) while maintaining < 15% affinity variation across variants (compared to > 40% for original R14 nanobody).

Discussion:

This study demonstrates that in silico affinity enhancement is a rapid and resource-efficient approach to repurpose nanobodies against SARS-CoV-2 variants, significantly accelerating affinity optimization while reducing experimental demands. This computational approach expedites the optimization of nanobody binding affinities while minimizing experimental resource requirements. By enhancing nanobody efficacy, our method provides a viable framework for developing targeted countermeasures against evolving SARS-CoV-2 variants and other pathogens.

1 Introduction

Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) continues to evolve with sustained global transmission, posing significant challenges to public health (1, 2). The virus relies on its spike glycoprotein (S protein) to bind with the angiotensin-converting enzyme 2 (ACE2) receptor on host cells, a critical step in viral entry (25). This interaction represents a critical target for therapeutic interventions, including monoclonal antibodies and antiviral drugs (6, 7). However, the rapid emergence of new variants carrying mutations in the receptor binding domain (RBD) of the S protein has raised concerns about the effectiveness of existing treatments and vaccines (1, 8, 9).

Notably, recent Omicron subvariants, including JN.1, KP.3, and XEC, exhibit enhanced transmissibility and immune evasion capabilities. The JN.1 variant, which first appeared in early 2024, has spread rapidly around the world, showcasing mutations that enhance its transmissibility. Specifically, JN.1 carries critical RBD mutations (R346T and F456L) that strengthen ACE2 receptor binding affinity and promote immune escape (1, 10, 11). Similarly, the KP.3 variant, a sublineage of JN.1, presents additional alterations in its RBD, particularly Q493E, which may further enhance its ability to evade neutralizing antibodies (1, 11). Meanwhile, the XEC variant, a recombinant of JN.1 subvariants KS.1 and KP.3.3, demonstrates a growth advantage in the population, suggesting enhanced virulence, greater immune evasion and caused breakthrough infectivity (12, 13). The World Health Organization (WHO) has classified JN.1 as a Variant of Interest (VOI) and KP.3 and XEC as Variants Under Monitoring (VUM), highlighting their potential global impact (14).

Nanobodies (Nbs), small single-domain antibodies derived from camelids, represent a promising therapeutic approach against SARS-CoV-2 (15). Structurally, nanobodies are composed of four highly conserved framework regions (FRs) interspersed with three hypervariable regions known as complementarity-determining regions (CDRs). These CDRs are responsible for antigen recognition and binding specificity (16, 17), making them critical for targeting viral epitopes. Compared to conventional monoclonal antibodies, nanobodies exhibited superior stability, simpler production processes, and enhanced capacity to recognize conformational epitopes properties (17) that make them particularly suitable for targeting the RBD of the S protein. Additionally, studies have shown that aerosol delivery achieves over 80% pulmonary deposition efficiency of nanobodies in murine models, ensuring localized antiviral activity at efficacious doses without systemic toxicity (18, 19). This delivery method enables direct lung targeting, maximizing local drug concentrations while minimizing systemic exposure. However, a critical limitation persists: most existing therapeutic antibodies, including nanobodies, were developed against early SARS-CoV-2 variants. The continuous accumulation of RBD mutations (17, 20) poses substantial challenges for both current treatment efficacy and future therapeutic development.

Conventional nanobody discovery methods, such as phage display and immunization of camelids, present significant limitations in terms of time investment and scalability (21). These techniques typically require months to years of experimental screening and optimization, compounded by challenges in reproducibility and large-scale production (2123). The 12-month half-life of SARS-CoV-2 antibody efficacy reported in 2024 virological surveys underscores the imperative to shorten discovery timelines (24). In contrast, computational design and optimization strategies offer a cost-effective and rapid alternative. By leveraging methods such as ab initio modeling, in silico mutagenesis, and machine learning, researchers can rapidly generate high-affinity binders with improved stability and solubility (2528). This approach not only accelerates the discovery of novel nanobodies but also enables the repurposing of existing ones to target emerging variants.

This study evaluates four class 1 (29) nanobodies [R14 (30), DL4 (31), VH ab6 (32), and Nanosota9 (33)], which were sequentially recognized as therapeutic candidates against evolving SARS-CoV-2 variants. Among them, VH ab6 and Nanosota9 have been recognized as broad-spectrum nanobodies because they demonstrated cross-variant neutralization activity. We systematically investigate their binding mechanisms across SARS-CoV-2 wild-type (WT) and variants, including BA.2, JN.1, and KP.3/XEC. Through detailed energetic and conformational analyses to elucidate the impact of viral mutations on binding affinity.

Furthermore, the development of nanobodies typically requires experimental screening of large antibody libraries, which is not only costly and time-consuming but also faces challenges such as low expression, poor solubility, and multispecificity (2123). These persistent limitations underscore the critical need for both innovative development approaches and efficient repurposing strategies for existing nanobodies (34, 35). To address these challenges, we present a streamlined computational pipeline (Figure 1) that generates a library of high-affinity and stable mutants by targeting substitutions in the CDRs and hotspot residues. The crystal structure of the nanobody in complex with the SARS-CoV-2 RBD served as template for high-throughput computational mutagenesis. High-affinity mutants were selected through protein-protein docking and molecular dynamics (MD) simulations, enabling the identification of favorable mutations with enhanced binding properties. This approach provides a robust framework for the development of effective therapeutic agents against SARS-CoV-2 and its evolving variants.

Figure 1

2 Methods

2.1 System preparation

The crystal structures of four SARS-CoV-2 RBD targeting nanobodies (3033) R14 (PDB ID: 7WD1) (30), DL4 (PDB ID: 7F5G) (31), VH ab6 (PDB ID: 8DLX) (32), and Nanosota9 (PDB ID:9CO9) (33) and SARS-CoV-2 WT (PDB ID: 7WD1) (30), BA.2 (PDB ID: 7ZF8) (36), JN.1 (PDB ID: 8Y5J) (37) RBDs were obtained from the RCSB Protein Data Bank (38). The KP.3/XEC RBD structure was generated by computational mutagenesis using RosettaDesign (39) with JN.1 RBD as template. These nanobodies were subsequently docked with the RBDs using the HDOCK server (40), and favorable binding conformations were selected based on docking score and confidence score (>0.99).

2.2 Molecular dynamics simulation

MD simulation was conducted using GROMACS 2022.5 with the CHARMM36 force field (41). The complex solvated in a cubic box using the TIP3P water model, and counter ions were added to neutralize the protein in the aqueous system. Energy minimization was executed through the steepest descent algorithm over 50,000 iterations to alleviate any high-energy contacts. Then, a two-step equilibration process was performed for 100 ps. First, canonical ensemble (NVT, constant number of particles, volume, and temperature) equilibration was conducted, maintaining the temperature at 310.15 K using the v-rescale thermostat, followed by isothermal-isobaric ensemble (NPT, constant number of particles, pressure, and temperature) equilibration with the c-rescale method for isobaric coupling at a reference pressure of 1.0 bar. Throughout the simulation, periodic boundary conditions (PBC) were applied, and long-range electrostatic interactions were calculated using the Particle Mesh Ewald (PME) method with a cutoff radius of 12 Å. Finally, a 200 ns production MD run was performed on the optimized system, with a time step of 2 fs, recording energy and log data every 1000 steps. The MD trajectory was analyzed to calculate the root mean square deviation (RMSD) and the number of hydrogen bonds, providing insights into the stability and interactions of the complexes.

2.3 Binding free energy calculations

Binding free energy calculations were performed using the Molecular Mechanics/Poisson-Boltzmann (Generalized Born) Surface Area (MM/PB(GB)SA) method. To ensure robust sampling, equilibrium trajectory data were extracted from the production phase of the MD simulations. Specifically, for each system, the last 50 ns of equilibrium trajectory data were used, with snapshots collected at 50 ps intervals. This methodology provided a total of 1,000 frames for binding free energy calculations.

The binding free energy () of the Nb-RBD complex can be regarded as the sum of enthalpy term (ΔH) and entropy term (-TΔS):

In this context, the enthalpy change is composed of ensemble average interaction energy and solvation energy (), in which can be separated into the electrostatic interaction energy (ΔEele) and van der Waals interaction energy (ΔEvdW).

The solvation energy () is further divided into polar solvation energy () and nonpolar solvation energy ():

In this study, a single trajectory approach was adopted due to its lower noise level and its suitability for systems with minimal structural rearrangement during binding. The entropy contribution was generally excluded from the calculations (42), as entropy differences associated with relative binding affinities were expected to be small, with only minor variations arising from mutations. MM/PB(GB)SA calculations were performed using gmx_MMPBSA version 1.6.3 (43), a novel tool developed for endpoint free energy calculations from MD trajectories.

To further elucidate the binding interactions, an energy decomposition analysis was conducted to assess the contribution of each residue to the binding free energy. This analysis highlights key residues that stabilize or destabilize the Nb-RBD complex, identifying potential targets for therapeutic intervention.

2.4 Generation of mutant library

The HDOCK generated Nb-RBD complexes were subjected to single-point mutations using DDMut-PPI (44). Single-point mutations were introduced at CDRs and hotspot residues, substituting each residue with the 19 alternative amino acids, excluding itself. The resulting changes in binding affinities were calculated as ΔΔG (ΔΔGGWTGmutant, kcal/mol), with ΔΔG > 0 indicating increased affinity and ΔΔG < 0 indicating decreased affinity compared to the wild type.

Empirical evidence and previous studies indicate that combinations of individually beneficial mutations are more likely to result in additive or synergistic improvements in binding affinity (45, 46). Therefore, we prioritized the combination of single-point mutations with ΔΔG > 0 for subsequent analysis. Stabilizing mutations (ΔΔG > 0) identified across variant complexes were subsequently combined for multiple mutations analysis using DDMut-PPI’s combinatorial module. Utilizing this approach, an in silico library comprising high-binding and stable single-point and multiple-points mutations in the nanobodies R14, DL4, VH ab6, and Nanosota9 were generated.

2.5 Interaction assessment and selection of best mutants

To validate our streamlined computational pipeline, mutations identified through computational affinity maturation were introduced into R14. The resulting mutant structures were docked with the RBD of the KP.3/XEC variant using HDOCK, employing its scoring function to filter for the best-performing mutants. The evaluation focused on predicted binding affinities and structural compatibility with the RBD.

Subsequently, the optimized combinations for R14 were subjected to a 200 ns production MD simulation with WT and three variants (BA.2, JN.1, and KP.3/XEC) complexes. The dynamic stability and interaction profiles of the mutants were analyzed using RMSD and MM/PB(GB)SA method. RMSD was utilized to assess conformational changes in the protein backbone, while MM/PB(GB)SA was employed to calculate binding free energies. These analyses enabled the identification of high-affinity and stable mutants, further refining the selection based on dynamic stability and interaction profiles. This approach ensured a robust assessment of the best-performing mutants.

3 Results

3.1 Anti-SARS-CoV-2 Omicron variants spectrum of nanobodies

In this study, we evaluated the antiviral spectrum of four RBD-targeting nanobodies (R14, DL4, VH ab6, and Nanosota9) against SARS-CoV-2 Omicron variants. Despite the RBD accumulating mutations at 3–5 times the rate of other spike protein domains (47), RBD-targeting antibodies remain clinically dominant due to their direct disruption of ACE2 binding (48). To elucidate their binding characteristics, we conducted a sequence alignment of RBDs from BA.2, JN.1, and KP.3/XEC. The analysis identified three critical mutations (L455S, F456L, and Q493E) within the RBD regions that directly interact with the nanobodies. Epitope mapping (Figure 2A) revealed that VH ab6 and Nanosota9 target relatively conserved epitopes on the Omicron RBD, whereas R14 and DL4 bind to more variable regions. These findings indicate that the binding affinities of VH ab6 and Nanosota9 exhibit minimal variability across different Omicron subvariants, enabling effective neutralization of the most predominant variants (Figures 2F-I).

Figure 2

Furthermore, structural analysis revealed that R14 primarily interacts with KP.3/XEC through main-chain functional groups, suggesting that side-chain variations exert minimal influence on binding (Figure 2B). This finding accounts for the partial neutralization efficacy of R14 against KP.3/XEC. In contrast, VH ab6 and Nanosota9 involve changes in both main-chain and side-chain interactions, leading to a more pronounced influence of side-chain groups on overall antibody affinity (Figures 2D, E). These structural insights align with biochemical and virological findings regarding the antiviral efficacy of each nanobody against earlier and recent Omicron subvariants.

3.2 Immune escape mechanisms of SARS-CoV-2 Omicron variants against nanobodies

To further investigate the interactions between nanobodies and SARS-CoV-2 Omicron variants, we performed MD simulations on sixteen systems, including four nanobodies (R14, DL4, VH ab6, and Nanosota9) and four RBDs (WT, BA.2, JN.1, and KP.3/XEC). We assigned the complexes formed by R14 with these RBDs as R14-RBD systems (R14-WT, R14-BA.2, R14-JN.1 and R14-KP.3/XEC). The naming convention for the complexes formed by the other nanobodies followed the same pattern. Throughout the simulations, all systems maintained stable binding states, with RMSD fluctuations consistently below 6 Å (Supplementary Figures S1A-D). Based on the equilibrium portion of the MD trajectory (150–200 ns), the free energy landscape (FEL) analysis (Supplementary Figure S2) was performed to explore the distribution of conformations and stability of Nb-RBD systems. In general, lower energy states indicate the corresponding conformations possess greater stability. Notably, the complexes formed between JN.1 and KP.3/XEC RBDs and the nanobodies exhibited increased flexibility at the binding interface compared to the complexes formed between the WT and BA.2 RBDs and the nanobodies. This enhanced flexibility suggests more dynamic interactions, potentially affecting binding characteristics and antibody efficacy.

To assess the effect of mutations on binding modes and energy changes, we calculated the binding free energies using the MM/PBSA method. The binding affinities of R14 and DL4 were significantly reduced, with the R14-KP.3/XEC having a binding affinity of -45.41 kcal/mol, representing increases of 36.94 kcal/mol and 39.85 kcal/mol, respectively, compared to their WT systems (Figures 3B, C). In contrast, the reduction in binding affinity for VH ab6 and Nanosota9 was more gradual, with increases of 13.96 kcal/mol and 3.99 kcal/mol, respectively, compared to their WT complexes, consistent with the analysis in section 3.1 regarding the antiviral spectrum of nanobodies (Figures 3D, E). These findings were further validated through MM/GBSA cross-verification (Supplementary Table S1), confirming the robustness of our predictions. Although the absolute binding energies varied between methods due to differences in the solvation model, both approaches showed consistent trends in changes to binding affinity across variants.

Figure 3

Overall, the binding affinities of JN.1 and KP.3/XEC variants to nanobodies were markedly lower than those of WT and BA.2 variants. This discrepancy primarily arises from a significant reduction in electrostatic and van der Waals (vdW) energies (Figure 3A). The former is attributed to the introduction of positively charged residues by mutations, while the latter results from unstable binding modes between the mutated RBD and nanobodies, leading to fewer atomic contacts. In addition, hydrogen bond network analysis (Supplementary Figure S3) revealed that the JN.1 and KP.3/XEC variants had more low-occupancy hydrogen bonds (occupancy < 70%) and fewer high-occupancy hydrogen bonds (occupancy ≥ 70%), which further impacted the stability of the complex.

Subsequently, we calculated the binding free energy contributions of residues near the binding interface using residue decomposition methods. Residues with an absolute binding free energy difference ≥ 1 kcal/mol compared to the WT system were defined as hotspot residues. In the R14-KP.3/XEC complex, a greater number of hotspot residues were identified than in the R14-JN.1 complex (Figure 4A, Supplementary Figures S4A, B), likely due to the unique F456L and Q493E mutations in KP.3/XEC. Notably, the Q493E mutation in the R14-KP.3/XEC complex exhibited a substantial energy change, attributed to the transformation from neutral glutamine to negatively charged glutamate at physiological pH, enhancing electrostatic repulsion between nanobody and RBD. This mutation also disrupted existing hydrogen bonds, further destabilizing the complex (Supplementary Figure S4). Similarly, the charged mutation at residue 493 significantly decreased electrostatic interactions in the DL4, VH ab6, and Nanosota9 systems (Figure 3A). In the DL4-KP.3/XEC complex, residues F489, G495, and R497 exhibited notable energy reductions, clustering in a similar region, suggesting substantial conformational changes. The VH ab6 and Nanosota9 complexes displayed comparable binding modes, primarily influenced by the Q493E mutation in KP.3/XEC (Figures 4C, D, Supplementary Figures S4E-H). However, the charged residues introduced by the mutations did not approach the binding interface closely enough to significantly impact binding affinity compared to R14 and DL4.

Figure 4

3.3 Generation of mutant library and assessment of binding affinity

The sixteen nanobody-RBD systems were analyzed using DDMut-PPI to predict high-affinity single-point mutations targeting CDRs and hotspot residues (Supplementary Tables S2-3). This computational affinity maturation approach generated a virtual library containing 2,508 mutations for R14, 1,824 for DL4, 2,432 for VH ab6, and 2,508 for Nanosota9. Screening for mutations with predicted ΔΔG > 0 (indicating enhanced affinity and stability, Supplementary Figure S6) identified 4 beneficial mutations for R14, 5 for DL4, 15 for VH ab6, and 9 for Nanosota9 (Figure 5, Supplementary Table S4).

Figure 5

Moreover, mutations within specific CDRs and hotspot residues not only alter the binding characteristics at the mutation sites but also affect the binding modes of other CDRs and framework regions. Analysis of combinatorial mutations revealed substantial variations in binding energy, with the multi-points mutant library containing 28 high-affinity variants for R14, 40 for DL4, 13,748 for VH ab6, and 176 for Nanosota9. Notably, VH ab6 displayed markedly greater mutational plasticity, accommodating significantly more viable multi-points mutation than other nanobodies while maintaining stable binding conformations.

3.4 Selection of high-affinity mutants using docking and MD simulation

High-affinity substitution mutations identified through preliminary screening using DDMut-PPI, including both single and multiple mutations, were introduced into the nanobodies (Supplementary Table S5). For R14 variants, docking against four SARS-CoV-2 variant RBDs using HDOCK identified three top-ranking mutation combinations exhibiting enhanced interaction (Supplementary Table S6).

To validate the feasibility of this nanobody repurposing pipeline, we conducted molecular dynamics simulations on the complexes of three optimized variants of R14 with four different SARS-CoV-2 variants to assess structural stability and energy changes following mutation. The RMSD was utilized to determine the conformational changes occurring in the protein backbone throughout the simulation, indicating the dynamic stability of the complex. Although the RMSD of the optimized structure was higher than that of the initial structure, reflecting internal structural changes that occurred during the optimization, it remained stable throughout the simulation. The RMSD fluctuations of RBD and R14 mutated nanobody complex remained consistently below 6 Å (Supplementary Figure S6), induced by minimal fluctuations and suggesting stable binding. The increased RMSD post-optimization, while stable, indicates a reorganization of the protein structure that may enhance the interaction with the target, thus aligning with the improved binding energy observed in the analysis.

Binding affinity serves as a key metric for evaluating the strength and specificity of nanobody-antigen interactions. We designate the complexes formed by the original R14 nanobody and its three optimized R14 as R14-ORI, R14-OPT1 (L29W/S52C/A101V), R14-OPT2 (L29W/A101V), and R14-OPT3 (L29W/S52C/A101L), respectively. Comprehensive binding free energy analysis (Figure 6) revealed substantial improvements in binding affinity against all SARS-CoV-2 variants mentioned above in this study, following nanobody optimization. The R14-ORI exhibited progressively weaker binding from WT (-82.35 kcal/mol) to KP.3/XEC (-45.41 kcal/mol), consistent with evolving immune evasion. Optimized Nb-RBD complexes demonstrated markedly enhanced and more uniform affinities. R14-OPT1 showed particularly strong JN.1 binding affinity, with a free energy of -76.88 kcal/mol compared to -47.3 kcal/mol for R14-ORI, while maintaining stable WT interactions at -67.04 kcal/mol. R14-OPT2 achieved exceptional KP.3/XEC binding affinity at -78.09 kcal/mol with consistent performance across variants. R14-OPT3 displayed optimal BA.2 affinity at -79.47 kcal/mol despite a moderate reduction for JN.1 to -55.69 kcal/mol. Notably, the binding free energy for the R14-OPT2-KP.3/XEC system exhibited slight discrepancies compared to the DDMut-PPI and HDOCK predictions, which can be attributed to the stronger stability of this system, resulting in lower overall perturbations in its structure compared to other similar systems.

Figure 6

Furthermore, R14-OPT1 exhibited the most balanced binding characteristics, with affinities ranging from -64.49 kcal/mol (KP.3/XEC) to -76.88 kcal/mol (JN.1), a variation of 16.1% compared to the 44.9% variation observed for R14-ORI (-45.41 to -82.35 kcal/mol). These energy landscapes suggest R14-OPT1 represent the optimal compromise between variant coverage and binding potency, particularly given its superior JN.1 recognition (-76.88 kcal/mol versus -47.3 kcal/mol for R14-ORI). Overall, the binding affinities of the optimized R14 nanobody with these four variants were comparatively favorable. These findings highlight the potential of integrating computational methods for affinity maturation in the development of effective therapeutics against SARS-CoV-2.

4 Discussion

When facing Omicron variants, countermeasures, such as vaccines and therapeutic drugs, display weaker or even lost effectiveness (1, 8, 9). Our integrated computational analyses reveal that the Q493E mutation in emerging SARS-CoV-2 variant KP.3/XEC drives immune evasion through electrostatic disruption at nanobody-RBD interfaces. This immune evasion effect is particularly pronounced for nanobodies like R14 and DL4, which engage variable epitopes, whereas VH ab6 and Nanosota9 maintain broader efficacy by targeting evolutionarily constrained regions. Critically, the structural and computational research of this part is consistent with the existing experimental studies (3033), which also confirms the feasibility of our next workflow. Building on these mechanistic insights, we developed an efficient computational pipeline that synergistically combines integrates high-throughput mutagenesis, protein-protein docking, and MD simulations to engineer optimized nanobody variants. Unlike conventional strategies focusing exclusively on single-point mutations in either CDRs or hotspot residues, our approach comprehensively targets both regions while incorporating structural dynamics from molecular simulations (7, 28, 49, 50). This pipeline specifically targets substitutions in the CDRs and hotspot residues, generating a mutant library of high-affinity and stable mutants. This integrated methodology generates mutant libraries enriched for high-affinity, stable nanobodies, ultimately providing atomistic insights to accelerate structure-guided nanobody design and therapeutic development.

Notably, R14 emerged as the primary candidate for optimization due to its unique therapeutic and structural properties. Aerosolized R14 maintained neutralizing activity and prevented infection (30), and it exhibited exceptional conformational stability across variants (RMSD <6 Å). This resilience stems from its predominant reliance on main-chain interactions with RBD, minimizing vulnerability to side-chain mutations. Residue decomp energetics revealed that Q493E-induced repulsion is a quantifiable liability (ΔΔG = +8 kcal/mol), offering a clear strategy for compensatory engineering through mutations.

The substantial binding affinity improvements achieved through computational maturation (e.g., R14-OPT1’s 62.6% ΔΔG enhancement against JN.1) hold significant implications for neutralization potency. Empirical calibrations indicate that, under typical Cheng–Prusoff conditions, every ~1.4 kcal/mol reduction in ΔGbind corresponds to an order-of-magnitude drop in half maximal inhibitory concentration (IC50) (51, 52). Extrapolating this relationship, the >15 kcal/mol ΔΔG enhancement observed for our top variants against KP.3/XEC points to multi-log neutralization gains, although the exact IC50 shift must ultimately be confirmed in functional assays. When benchmarked against literature-reported broad-spectrum nanobodies VH ab6 (32) and Nanosota9 (33), R14-OPT1 achieved higher binding affinity despite targeting a more plastic epitope. This demonstrates that computational repurposing can confer breadth even on epitopes traditionally considered mutationally vulnerable.

Despite providing comprehensive analysis, this study has several limitations that need to be considered. While our pipeline predicted high-affinity nanobodies, direct comparison to clinical-stage anti-coronavirus nanobodies was constrained by their absence in late development pipelines (33). Furthermore, experimental validation of the identified mutations and the exploration of additional mutation combinations remains essential. And our strategy holds promises for extension beyond SARS-CoV-2. Residues that are critical for the development of broad-spectrum nanobodies are conserved across betacoronaviruses, such as SARS-CoV and MERS-CoV. Future work should leverage our pipeline to engineer pan-sarbecovirus therapeutics, potentially integrated with machine learning approaches to identify key residues and optimize nanobody interactions against a broader range of viral variants.

5 Conclusion

This study confirmed that the Q493E mutation is a key driver of immune evasion in KP.3/XEC through disruption of electrostatic nanobody-RBD interactions. To overcome this challenge, we developed a streamlined computational pipeline integrating high-throughput mutagenesis of CDRs and hotspot residues, protein-protein docking, and MD simulations. Using this approach, we repurposed four nanobodies (R14, DL4, VH ab6, and Nanosota9) to target the Omicron RBD. Notably, MD simulations validated the best-performing R14 optimized combination, R14-OPT1 systems showed a 62.6% binding energy improvement against JN.1, with a consistent binding across variants (<15% affinity variation). These results collectively demonstrate our pipeline’s capacity to significantly improve nanobody binding affinity (achieving >60% enhancement for key variants) while maintaining broad neutralization capacity. Our pipeline provides critical insights into the interactions between nanobodies and evolving viral variants, supporting the potential use of existing nanobodies as therapeutic agents.

Statements

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.

Author contributions

SC: Visualization, Conceptualization, Validation, Formal Analysis, Methodology, Data curation, Investigation, Writing – original draft. BS: Visualization, Formal Analysis, Validation, Conceptualization, Data curation, Writing – review & editing, Investigation, Methodology. FG: Project administration, Funding acquisition, Resources, Validation, Supervision, Writing – review & editing.

Funding

The author(s) declare financial support was received for the research and/or publication of this article. This work was supported by the National Natural Science Foundation of China (Grant numbers 32270692 and 31571358).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The author(s) declared that they were an editorial board member of Frontiers, at the time of submission. This had no impact on the peer review process and the final decision.

Generative AI statement

The author(s) declare that no Generative AI was used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2025.1637955/full#supplementary-material

References

  • 1

    GillotCDavidCDognéJMCaboJDouxfilsJFavresseJ. Neutralizing antibodies against KP.2 and KP.3: why the current vaccine needs an update. Clin Chem Lab Med. (2024) 63:e82–5. doi: 10.1515/cclm-2024-0919

  • 2

    ZhaoXGaoF. Novel Omicron variants enhance anchored recognition of TMEM106B: A new pathway for SARS-CoV-2 cellular invasion. J Phys Chem Lett. (2024) 15:671–80. doi: 10.1021/acs.jpclett.3c03296

  • 3

    HoffmannMKleine-WeberHSchroederSKrügerNHerrlerTErichsenSet al. SARS-CoV-2 cell entry depends on ACE2 and TMPRSS2 and is blocked by a clinically proven protease inhibitor. Cell. (2020) 181:271280.e8. doi: 10.1016/j.cell.2020.02.052

  • 4

    ZhouPYangXLWangXGHuBZhangLZhangWet al. A pneumonia outbreak associated with a new coronavirus of probable bat origin. Nature. (2020) 579:270–3. doi: 10.1038/s41586-020-2012-7

  • 5

    YanFFGaoF. Comparison of the binding characteristics of SARS-CoV and SARS-CoV-2 RBDs to ACE2 at different temperatures by MD simulations. Brief Bioinform. (2021) 22:1122–36. doi: 10.1093/bib/bbab044

  • 6

    CaoYSuBGuoXSunWDengYBaoLet al. Potent neutralizing antibodies against SARS-CoV-2 identified by high-throughput single-cell sequencing of convalescent patients’ B cells. Cell. (2020) 182:7384.e16. doi: 10.1016/j.cell.2020.05.025

  • 7

    SunBGaoF. Investigation of escape mechanisms of SARS-CoV-2 Omicron sub-lineages and exploration of potential antibodies for XBB.1. J Infect. (2023) 87:354–7. doi: 10.1016/j.jinf.2023.07.014

  • 8

    LiuJLuoSXuXZhangELiangHZhangJZHet al. Evaluating the synergistic effects of multi-epitope nanobodies on BA.2.86 variant immune escape. J Phys Chem Lett. (2025) 16:396405. doi: 10.1021/acs.jpclett.4c03028

  • 9

    XuKAnYLiuXXieHLiDYangTet al. Neutralization of SARS-CoV-2 KP.1, KP.1.1, KP.2 and KP.3 by human and murine sera. NPJ Vaccines. (2024) 9:215. doi: 10.1038/s41541-024-01016-6

  • 10

    KakuYUriuKKosugiYOkumuraKYamasobaDUwaminoYet al. Virological characteristics of the SARS-CoV-2 KP.2 variant. Lancet Infect Dis. (2024) 24:e416. doi: 10.1016/S1473-3099(24)00298-6

  • 11

    KakuYYoMSTolentinoJEUriuKOkumuraKSatoK. Virological characteristics of the SARS-CoV-2 KP.3, LB.1, and KP.2.3 variants. Lancet Infect Dis. (2024) 24:e482–3. doi: 10.1016/S1473-3099(24)00415-8

  • 12

    LiuJYuYJianFYangSSongWWangPet al. Enhanced immune evasion of SARS-CoV-2 variants KP.3.1.1 and XEC through N-terminal domain mutations. Lancet Infect Dis. (2025) 25:e6–7. doi: 10.1016/S1473-3099(24)00738-2

  • 13

    JianFWangJYisimayiASongWXuYChenXet al. Evolving antibody response to SARS-CoV-2 antigenic shift from XBB to JN.1. Nature. (2025) 637:921–9. doi: 10.1038/s41586-024-08315-x

  • 14

  • 15

    KoenigPADasHLiuHKümmererBMGohrFNJensterLMet al. Structure-guided multivalent nanobodies block SARS-CoV-2 infection and suppress mutational escape. Science. (2021) 371:eabe6230. doi: 10.1126/science.abe6230

  • 16

    MuyldermansS. Nanobodies: natural single-domain antibodies. Annu Annu Rev Biochem. (2013) 82:775–97. doi: 10.1146/annurev-biochem-063011-092449

  • 17

    YangYLiFDuL. Therapeutic nanobodies against SARS-CoV-2 and other pathogenic human coronaviruses. J Nanobiotechnol. (2024) 22:304. doi: 10.1186/s12951-024-02573-7

  • 18

    XuJXuKJungSConteALiebermanJMueckschFet al. Nanobodies from camelid mice and llamas neutralize SARS-CoV-2 variants. Nature. (2021) 595:278–82. doi: 10.1038/s41586-021-03676-z

  • 19

    Van HeekeGAlloseryKDe BrabandereVDe SmedtTDetalleLde FougerollesA. Nanobodies® as inhaled biotherapeutics for lung diseases. Pharmacol Ther. (2017) 169:4756. doi: 10.1016/j.pharmthera.2016.06.012

  • 20

    LiMRenYAwZQChenBYangZLeiYet al. Broadly neutralizing and protective nanobodies against SARS-CoV-2 Omicron subvariants BA.1, BA.2, and BA.4/5 and diverse sarbecoviruses. Nat Commun. (2022) 13:7957. doi: 10.1038/s41467-022-35642-2

  • 21

    LedsgaardLKilstrupMKaratt-VellattAMcCaffertyJLaustsenAH. Basics of antibody phage display technology. Toxins (Basel). (2018) 10:236. doi: 10.3390/toxins10060236

  • 22

    PardonELaeremansTTriestSRasmussenSGWohlkönigARufAet al. A general protocol for the generation of Nanobodies for structural biology. Nat Protoc. (2014) 9:674–93. doi: 10.1038/nprot.2014.039

  • 23

    SchoofMFaustBSaundersRASangwanSRezeljVHoppeNet al. An ultra-potent synthetic nanobody neutralizes SARS-CoV-2 by locking Spike into an inactive conformation. Science. (2020) 370:1473–9. doi: 10.1101/2020.08.08.238469

  • 24

    WalmsleySNabipoorMQiFLovblomLERavindranRColwillKet al. Declining levels of neutralizing antibodies to SARS-CoV-2 Omicron variants are enhanced by hybrid immunity and original/Omicron bivalent vaccination. Vaccines (Basel). (2024) 12:564. doi: 10.3390/vaccines12060564

  • 25

    KiyoshiMCaaveiroJMMiuraENagatoishiSNakakidoMSogaSet al. Affinity improvement of a therapeutic antibody by structure-based computational design: generation of electrostatic interactions in the transition state stabilizes the antibody-antigen complex. PloS One. (2014) 9:e87099. doi: 10.1371/journal.pone.0087099

  • 26

    YanFFGaoF. EK1 with dual Q1004E/N1006I mutation: a promising fusion inhibitor for the HR1 domain of SARS-CoV-2. J Infect. (2022) 84:579613. doi: 10.1016/j.jinf.2021.12.022

  • 27

    CannonDAShanLDuQShirinianLRickertKWRosenthalKLet al. Experimentally guided computational antibody affinity maturation with de novo docking, modelling and rational design. PloS Comput Biol. (2019) 15:e1006980. doi: 10.1371/journal.pcbi.1006980

  • 28

    SinghVChoudharySBhutkarMNehulSAliSSinglaJet al. Designing and bioengineering of CDRs with higher affinity against receptor-binding domain (RBD) of SARS-CoV-2 Omicron variant. Int J Biol Macromol. (2025) 290:138751. doi: 10.1016/j.ijbiomac.2024.138751

  • 29

    JiangJBoughterCTAhmadJNatarajanKBoydLFMeier-SchellersheimMet al. SARS-CoV-2 antibodies recognize 23 distinct epitopic sites on the receptor binding domain. Commun Biol. (2023) 6:953. doi: 10.1038/s42003-023-05332-w

  • 30

    LiuHWuLLiuBXuKLeiWDengJet al. Two pan-SARS-CoV-2 nanobodies and their multivalent derivatives effectively prevent Omicron infections in mice. Cell Rep Med. (2023) 4:100918. doi: 10.1016/j.xcrm.2023.100918

  • 31

    LiTZhouBLiYHuangSLuoZZhouYet al. Isolation, characterization, and structure-based engineering of a neutralizing nanobody against SARS-CoV-2. Int J Biol Macromol. (2022) 209:1379–88. doi: 10.1016/j.ijbiomac.2022.04.096

  • 32

    MannarDSavilleJWSunZZhuXMartiMMSrivastavaSSet al. SARS-CoV-2 variants of concern: spike protein mutational analysis and epitope for broad neutralization. Nat Commun. (2022) 13:4696. doi: 10.1038/s41467-022-32262-8

  • 33

    YeGBuFSaxenaDTurner-HubbardHHerbstMSpillerBet al. Discovery of Nanosota-9 as anti-Omicron nanobody therapeutic candidate. PloS Pathog. (2024) 20:e1012726. doi: 10.1371/journal.ppat.1012726

  • 34

    SormanniPAprileFAVendruscoloM. Rational design of antibodies targeting specific epitopes within intrinsically disordered proteins. Proc Natl Acad Sci U.S.A. (2015) 112:9902–7. doi: 10.1073/pnas.1422401112

  • 35

    Adolf-BryfogleJKalyuzhniyOKubitzMWeitznerBDHuXAdachiYet al. RosettaAntibodyDesign (RAbD): A general framework for computational antibody design. PloS Comput Biol. (2018) 14:e1006112. doi: 10.1371/journal.pcbi.1006112

  • 36

    NutalaiRZhouDTuekprakhonAGinnHMSupasaPLiuCet al. Potent cross-reactive antibodies following Omicron breakthrough in vaccinees. Cell. (2022) 185:2116–31. doi: 10.1016/j.cell.2022.05.014

  • 37

    LiLShiKGuYXuZShuCLiDet al. Spike structures, receptor binding, and immune escape of recently circulating SARS-CoV-2 Omicron BA.2.86, JN.1, EG.5, EG.5.1, and HV.1 sub-variants. Structure. (2024) 32:1055–67. doi: 10.1016/j.str.2024.06.012

  • 38

    BermanHMWestbrookJFengZGillilandGBhatTNWeissigHet al. The protein data bank. Nucleic Acids Res. (2000) 28:235–42. doi: 10.1093/nar/28.1.235

  • 39

    LemanJKWeitznerBDLewisSMAdolf-BryfogleJAlamNAlfordRFet al. Macromolecular modeling and design in Rosetta: recent methods and frameworks. Nat Methods. (2020) 17:665–80. doi: 10.1038/s41592-020-0848-2

  • 40

    YanYZhangDZhouPLiBHuangSY. HDOCK: a web server for protein-protein and protein-DNA/RNA docking based on a hybrid strategy. Nucleic Acids Res. (2017) 45:W365–73. doi: 10.1093/nar/gkx407

  • 41

    CroitoruAParkSJKumarALeeJImWMacKerellADJret al. Additive CHARMM36 force field for nonstandard amino acids. J Chem Theory Comput. (2021) 17:3554–70. doi: 10.1021/acs.jctc.1c00254

  • 42

    HomeyerNGohlkeH. Free energy calculations by the molecular mechanics poisson-boltzmann surface area method. Mol Inform. (2012) 31:114–22. doi: 10.1002/minf.201100135

  • 43

    Valdés-TresancoMSValdés-TresancoMEValientePAMorenoE. gmx_MMPBSA: A new tool to perform end-state free energy calculations with GROMACS. J Chem Theory Comput. (2021) 17:6281–91. doi: 10.1021/acs.jctc.1c00645

  • 44

    ZhouYMyungYRodriguesCHMAscherDB. DDMut-PPI: predicting effects of mutations on protein-protein interactions using graph-based deep learning. Nucleic Acids Res. (2024) 52:W207–14. doi: 10.1093/nar/gkae412

  • 45

    TongZTongJLeiWXieYCuiYJiaGet al. Deciphering a reliable synergistic bispecific strategy of rescuing antibodies for SARS-CoV-2 escape variants, including BA.2.86, EG.5.1, and JN.1. Cell Rep. (2024) 43:114338. doi: 10.1016/j.celrep.2024.114338

  • 46

    MyungYPiresDEVAscherDB. mmCSM-AB: guiding rational antibody engineering through multiple point mutations. Nucleic Acids Res. (2020) 48:W125–31. doi: 10.1093/nar/gkaa389

  • 47

    PiccoliLParkYJTortoriciMACzudnochowskiNWallsACBeltramelloMet al. Mapping neutralizing and immunodominant sites on the SARS-CoV-2 spike receptor-binding domain by structure-guided high-resolution serology. Cell. (2020) 183:102442.e21. doi: 10.1016/j.cell.2020.09.037

  • 48

    ShangJYeGShiKWanYLuoCAiharaHet al. Structural basis of receptor recognition by SARS-CoV-2. Nature. (2020) 581:221–4. doi: 10.1038/s41586-020-2179-y

  • 49

    YanFFGaoF. RBD-ACE2 binding properties in five SARS-CoV-2 variants of concern with new perspectives in the design of pan-coronavirus peptide inhibitors. J Infect. (2023) 86:e51–4. doi: 10.1016/j.jinf.2022.09.011

  • 50

    ShahMShinJYWooHG. Rational strategies for enhancing mAb binding to SARS-CoV-2 variants through CDR diversification and antibody-escape prediction. Front Immunol. (2023) 14:1113175. doi: 10.3389/fimmu.2023.1113175

  • 51

    ChengHC. The power issue: determination of KB or Ki from IC50. A closer look at the Cheng-Prusoff equation, the Schild plot and related power equations. J Pharmacol Toxicol Methods. (2001) 46:6171. doi: 10.1016/s1056-8719(02)00166-1

  • 52

    GiraldoJSerraJRocheDRoviraX. Assessing receptor affinity for inverse agonists: Schild and Cheng-Prusoff methods revisited. Curr Drug Targets. (2007) 8:197202. doi: 10.2174/138945007779315687

Summary

Keywords

SARS-CoV-2, immune evasion, nanobody, computational pipeline, affinity maturation

Citation

Cao S, Sun B and Gao F (2025) From immune evasion to broad in silico binding: computational optimization of SARS-CoV-2 RBD-targeting nanobody. Front. Immunol. 16:1637955. doi: 10.3389/fimmu.2025.1637955

Received

30 May 2025

Accepted

29 July 2025

Published

14 August 2025

Volume

16 - 2025

Edited by

Maolin Lu, University of Texas at Tyler Health Science Center, United States

Reviewed by

Chaofan Li, University of Virginia, United States

Chang Liu, The University of Texas at Austin, United States

Caihong Bi, Harvard Medical School, United States

Updates

Copyright

*Correspondence: Feng Gao,

†These authors have contributed equally to this work

Disclaimer

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.

Outline

Figures

Cite article

Copy to clipboard


Export citation file


Share article

Article metrics