Mutated amino acids of spike glycoprotein of Alpha and Delta variants were retrieved from the GISAID database (https://www.gisaid.org/). Amino acid sequence of native spike glycoprotein (Accession ID: QHD43416.1) of SARS-CoV-2 and therapeutic mAbs were retrieved from the NCBI (https://www.ncbi.nlm.nih.gov/) and CoV-AbDab database (http://opig.stats.ox.ac.uk/webapps/covabdab/), respectively (27).

Homology modeling is a template-dependent/independent method popularly used to model protein structure from its amino acid sequence. Based on the templates available in the database repository, an automated modeling server, SWISS-MODEL, was used to model all the mutant S proteins (28). Similarly, the web application ABodyBuilder, a tool for small-scale homology modeling (29), was applied for designing the Fv regions of the mAbs used in this study. The ABodyBuilder algorithm works through the following pathway, i.e., selection of template>prediction of orientation>prediction of the side chain>modeling of complementarity determining region loop to design computational model of the mAbs (30). Antibody i-Patch, a web tool that works on the antibody-specific statistics to determine the paratope and complementarity-determining regions (CDR) within the antibody, was also utilized to validate the efficacy of the paratopes towards an antigen (31). In addition, the webtool EpiPred was explored to determine the epitopes by using the homology model of the antibody as an input (32, 33). Each structure was verified using SAVES server that examines the stereochemical quality of a predicted protein structure through analyzing the residue-by-residue geometry and geometry of the overall structure (34, 35).

The protein–protein interactions between the S protein variants and human mAbs were determined through molecular docking study using High Ambiguity Driven protein–protein DOCKing (HADDOCK) v. 2.4. This flexible docking program executes the docking process by using the information from known and/or predicted protein interfaces in ambiguous interaction restraints (AIR) and the outputs generated by this server were found to support experimentally validated (NMR and cryo-EM) structures of protein complexes (36). Herein, the PDB files of each SARS-CoV-2 variant and mAb were subjected to binding interactions in HADDOCK platform and the protein–protein interaction was predicted through a binding score provided by the webserver as well as analyzing the interacting residues within the output file. All the docked complexes were analyzed for biophysical interactions using Discovery Studio-2020 and PyMOL.

Stability of the molecular topology, conformational topology, and dynamic behavior of the mAb–spike protein complex was studied through normal mode analyses (NMA) (37). Different protein complexes were analyzed for molecular flexibility and modal trajectories of structural dynamics using iMODS server from ChaconLab (38). NMA is a quantitative expression of motion of protein in a complex. In our study, NMA analyses were conducted to assess the biophysical attributes of S protein–mAb complexes. The iMODS server is a popular and customizable server that discloses coarse-grain (CG) levels and also determines the dihedral coordinates of C-alpha atoms in the proteins. iMods also delivers an NMA mobility that illustrates the collective motions and affine-model arrows that direct domain dynamics. Structural deformability usually depends on the helical content of a protein that is usually studied by analyzing the B-factor (39). Here, we studied the amplitude of atomic displacement in each spike protein–mAb complex by determining the B-factor. The relative impact of each deformation movement in the motion of interacting proteins within the S protein–mAb complex was studied by measuring the eigenvalues. Lower eigenvalues are indicative of greater stability while higher values reveal the reverse (4). The molecular flexibility, amplitude of fluctuation of each mode, and overall confirmational change in each S protein–mAb complex were studied by variance plot and covariance matrix. The docked complexes comprising S protein variants and mAbs were prepared in pdb format and submitted to the iMODS server to determine the normal modes within internal coordinates and details of the mobility (B-factor), structural deformability, covariance map, and linking matrix with eigenvalues. The data extracted from iMODS were analyzed and plotted using the statistical software package R. Changes in the conformation in S protein after binding of mAb were assessed by superimposing the unbound S protein structures on the S protein–mAb complex following our earlier reports (9, 39).

Selected variant S proteins of SARS-CoV-2 and the therapeutic mAbs were examined utilizing PRODIGY web tool to determine its binding affinity (ΔG) and dissociation constant (K _d) (40).

The rationale behind using in silico docking analysis was to check the binding efficacy of the mAbs against the S protein variants. The clear demonstration of the biomolecular interaction among S protein variants and mAbs was possible only through molecular docking. Hence, it was necessary to screen the possible efficacious mAbs and to obtain a precise picture on their respective CDRs to predict the possible interacting domains present in those antibodies. Such CDRs are important components of the chimeric antibody generated in silico.

After identifying the potential CDR patches in the screened mAbs, a series of chimeric antibodies were hypothesized for better efficacy. The binding domain of the mAbs were subjected to multiple sequence alignment and the peptide fragments showing binding affinity against SARS-CoV-2 variants were combined in silico. The resultant CDRs were fitted in a single amino acid chain to model the antibody structure by employing Therapeutic Antibody Profiler (TAP) (41, 42). The designed chimeric antibody was characterized for its immune-biological properties and efficacy against the S protein in silico.

For all the eight mAbs, 10 spike protein variants from Alpha and 10 from Delta lineages (shown in Figure S1 ) were tested for all possible interactions in silico ( Table S1 ). A total of 160 mAb–spike protein complexes were obtained using Haddock 2.4 and were further screened by setting the docking score >100 as strong interaction while <50 as weak interaction (43). With these criteria, we studied 3 strong and 3 weak spike protein–mAb complexes for each mAb. The rationale behind such a selection was to filter the mAbs that are effective against a greater number of SARS-CoV-2 variants and also to understand the possible variants that could escape binding of the mAb. Regdanvimab, bamlanivimab, sotrovimab, etesevimab, and cilgavimab showed high docking scores against most of the Alpha and Delta variants ( Figure S2 ). Particularly, regdanvimab displayed a high binding score against both Alpha and Delta variants. In contrast, tixagevimab, casirvimab, imdevimab, and cilgavimab possessed less docking scores. Docking score is a preliminary measurement of binding efficacy of the mAbs against the SARS-CoV-2 variants, and therefore, to validate the initial postulation as well as to obtain molecular insight on these interactions, we further refined our screening by including a total of 16 protein complexes comprising S protein variants from the Alpha lineage and 8 mAbs ( Table S1 ). This was repeated for the Delta lineage.

The biophysical basis of mAb-SARS-CoV-2 spike protein interactions was investigated through determining the binding free energy and binding affinity. Eventually, three strongly bound and three weakly bound protein complexes were selected from each lineage ( Table 1 ). P681H–tixagevimab, S982A–regdanvimab, and V70–cilgavimab complexes displayed significantly strong binding interactions with respective binding free energy (ΔG) of −15.1, −14.7, and −14.7 kcal/mol. Analysis of binding affinity in terms of the dissociation constant (K _d) supported the inference of binding energy. The data revealed that P681H–tixagevimab, S982A–regdanvimab, and V70–cilgavimab possess K _d values of 7.8E−12, 1.6E−11, and 1.6E−11 respectively. These data suggested that tixagevimab, regdanvimab, and cilgavimab could efficiently bind to the SARS-CoV-2 variants of B.1.1.7 lineage. In contrast, weak interactions were noted for S982A–tixagevimab, T716I–cilgavimab, and S982A-etesevimab, indicating that the variants with spike protein mutations like S982A could escape the inhibitory effect of certain mAbs.

For B.1.617.2 lineage, P681R–bamlanivimab, P681R–tixagevimab, and R158–sotrovimab complexes revealed strong interactions with respective binding free energies (ΔG) of −14.1, −13.9, and −13.3 kcal/mol. Furthermore, SARS-CoV-2 that harbors S protein mutations T19R and G142D were predicted as the variants that could escape from the inhibitory effects of aforesaid mAbs.

Taking clues from earlier studies, the top three strongly bound and three weakly bound mAb–spike protein variant complexes from each lineage were studied for protein–protein interactions at the molecular level ( Table 1 ). The topological differences in binding pattern of different mAbs to S protein variants as well as involvement of various non-covalent bonds/forces are summarized in Figure S3 and Table 2 . While studying the most strongly bound complex for the B.1.1.7 lineage, the P681H–tixagevimab complex was found to be stabilized by the most number (a total of 26 H-bonds) of hydrogen bonds. In addition, hydrophobic interaction (π–Alkyl: 1) was found to provide additional stability in holding the two proteins ( Figure S3A and Table 2 ). The other strongly bound complexes, viz., S982A–regdanvimab and V70–cilgavimab, were also found to possess similar numbers of H-bonds. However, additional interactive forces (4 hydrophobic interactions; π–Alkyl: 3 and alkyl: 1; and electrostatic bonds: 2) were noted in V70–cilgavimab ( Figure S3B and Table 2 ). However, the topology of binding was dissimilar (as evident from the interacting amino acids) to that of tixagevimab and that could explain the reduced binding affinity. In comparison to the high-affinity antibodies, we found lesser degree of interactions in terms of hydrogen bonding and/or other noncovalent interactions leading to weak binding of mAb to S variants as evident for T716I–cilgavimab, S982A–tixagevimab, and S982A-etesevimab ( Figures S3A, B and Table S3 ).

Biomolecular interactions between mAbs and S protein variants were also studied for the B.1.617.2 lineage. P681R–bamlanivimab, P681R–tixagevimab, and R158–sotrovimab displayed a strong binding as evidenced from their respective binding free energy ( Table 1 , Figures S3C, D and Tables 2 and S3 ). Considering the topology of binding ( Figure S3C ), bamlanivimab exhibited a very strong binding with S protein variant P681R, and this binding was stabilized by 16 H bonding, 1 electrostatic interaction, and 6 hydrophobic interactions (π–σ: 2; π–π: 1; π–Alkyl: 2 and alkyl: 1). Interestingly, tixagevimab was also found to interact with the same variant with a similar binding topology but with less affinity than that of bamlanivimab ( Figure S3C ). The abundance of H-bonds was also similar, and hence, lesser number of hydrophobic interactions could be responsible for the differences in the avidity of tixagevimab and bamlanivimab, while sotrovimab possessed strong affinity to the R158 variant and this interaction was found to be stabilized by non-covalent forces including H-bonding ( Table 2 ). Weak interactions were observed for T19R–bamlanivimab, G142D–tixagevimab, and T19R–regdanvimab ( Figure S3D and Table S3 ).

We scrutinized the conformational changes in S1 protein following its binding to mAb. Changes in the conformation of heterotrimeric S1 protein after binding of an antibody is key for the neutralization of SARS-CoV-2 (44). We observed clear changes in the conformation of S protein variants following interaction with the strongly binding mAbs ( Figure S4A ). This observation supports the earlier data that demonstrated that tixagevimab is the strongest neutralizing mAb against SARS-CoV-2 variants. However, binding of cilgavimab to T716I revealed no sign of conformational changes and indicated that T716I could evade cilgavimab ( Figure S4B ). The mode of binding of mAbs and subsequent changes in the conformation were also verified by comparing tixagevimab–P681H and cilgavimab–T716I protein complexes as well as a superimposed form of them ( Figure S4C ).

We have also studied the molecular dynamics of S protein–mAb complexes. The comparative analyses of molecular dynamics of the most stable S protein–mAb complex, i.e., P681H–tixagevimab, and the weakest bound complex, i.e., T716I–cilgavimab, revealed significant differences in various intramolecular and intermolecular parameters ( Figure S5 ). Clear difference in the direction of molecular motion was observed between the two protein complexes ( Figures S5A–F ). A relatively lower level of deformability described the compactness of P681H–tixagevimab while the reverse was observed for T716I–cilgavimab. Furthermore, protein components in the T716I–cilgavimab complex had higher mobility, indicating weak association, while the P681H–tixagevimab complex had lesser mobility, thus denoting the strong association between components of the complex ( Figures S5E, F ). Eigenvalue is a measure of deformation due to fluctuation in protein motion (38). A high eigenvalue is indicative of a localized displacement while low eigenvalue indicates cumulative conformational changes in the protein structure (38, 39). A low eigenvalue of 3.407305e^-05 for P681H–tixagevimab revealed lower energy deformation of the structure resulting in greater stability of the complex and conformational changes in the S protein after forming complex with mAb ( Figure S5G ). In contrast, T716I–cilgavimab displayed the reverse characteristics and interpreted as relatively unstable complex.

Structural integrity of the S protein–mAb was also investigated by molecular flexibility by measuring the theoretical fluctuation (variance) of each mode in a protein complex (37). A higher degree of flexibility was observed for P681H–tixagevimab, suggesting a higher affinity of two proteins towards each other while the reverse was observed for T716I–cilgavimab ( Figures S5I, J ). Covariance maps demonstrated abundance of correlated motions of Cα atoms in P681H–tixagevimab, thus indicating motion stiffness and rigidity in the protein complex compared to T716I–cilgavimab ( Figures S5K, L ).

Our in silico analyses predicted that a number of S protein variants depict very weak binding to mAbs, and therefore, these mAbs might not be effective against newly emerged strains especially Delta plus strain (B.1.617.2.1). This has prompted us to design a chimeric mAb that could be effective in encountering most of the variants ( Table S3 ). We have combined the CDRs of the strongly interactive mAb to prepare a series of chimeric mAbs and tested against the spike protein variants ( Tables S3 and S4A, B ). Eventually, when we incorporated CDRH3 of regdanvimab (ARIPGFLRYRNRYYYYGMDV) within the framework of sotrovimab without hampering its CDR (ARDYTRGAWFGESLIGGFDN), the resultant paratope was found to display very strong binding with most of the escape variants from both lineages ( Figures 1A–C and Table S5 ). Even the physico-biochemical properties, like solubility and hydropathicity of the selected chimeric mAb, supported its durability in neutralizing the variants ( Table S6 ). When we compared the efficacy of this mAb against the newly found Delta plus strain (B.1.617.2.1) and Beta strain (B.1.351), a very satisfactory binding interaction was observed as compared to the eight therapeutic mAbs ( Figure 1C and Table S7 ). However, binding efficacy of the chimeric antibody against Gamma strain (P.1) was comparatively lower than that of the Delta, Delta plus, and Beta strain ( Figure 1C , Tables S5 and S7 ).