Structural models of human ACE2 variants with SARS-CoV-2 Spike protein for structure-based drug design

Emergence of coronaviruses poses a threat to global health and economy. The current outbreak of SARS-CoV-2 has infected more than 28,000,000 people and killed more than 915,000. To date, there is no treatment for coronavirus infections, making the development of therapies to prevent future epidemics of paramount importance. To this end, we collected information regarding naturally-occurring variants of the Angiotensin-converting enzyme 2 (ACE2), an epithelial receptor that both SARS-CoV and SARS-CoV-2 use to enter the host cells. We built 242 structural models of variants of human ACE2 bound to the receptor binding domain (RBD) of the SARS-CoV-2 surface spike glycoprotein (S protein) and refined their interfaces with HADDOCK. Our dataset includes 140 variants of human ACE2 representing missense mutations found in genome-wide studies, 39 mutants with reported effects on the recognition of the RBD, and 63 predictions after computational alanine scanning mutagenesis of ACE2-RBD interface residues. This dataset will help accelerate the design of therapeutics against SARS-CoV-2, as well as contribute to prevention of possible future coronaviruses outbreaks.


Database search.
To identify all relevant variants of ACE2, we performed a search in multiple databases and created workflow for assembling the variants in complex with the SARS-CoV-2 RBD (Fig. 1). For variants naturally occurring in the human population, we searched gnomAD 16 and identified 155 unique missense ACE2 variants, 140 of which are mapped on the structural model (Online-only Table 1). For variants of ACE2 with

CC
GG C A C U 5 ' 5 ' 3 ' T 3 ' 3 ' 5 ' CG GC G C Fig. 1 Schematic overview of the structure-based benchmark of ACE2 variants-S protein complexes. All available variants are collected from (a) missense mutations identified in the human genome; (b) overexpressed constructs of ACE2 variants reported in the literature; and (c) designed alanine scanning mutagenesis variants of ACE2, targeting the interface residues with the S protein (upper panels). In the bottom panel, a structurebased benchmark including all variants is assembled for use in drug development, and optimization of the interface of the variant by including the Zn +2 ion is performed using the HADDOCK software. Zn 2+ is represented magnified, because it was considered for calculations. www.nature.com/scientificdata www.nature.com/scientificdata/ In silico alanine scanning mutagenesis. The initial model 11 was used to calculate interface residues by considering all residue-residue pairs of the wild-type ACE2 and the SARS-CoV-2 RBD within 10 Å distance of each other. These positions were then individually mutated to alanine residues (Ala). In total, we selected 63 residues, plus 6 which were alanine residues in the wild-type sequence and served as positive controls ( Table 2). interface refinement. Heterodimers of ACE2 variants and SARS-CoV-2 RBD were extracted from all three datasets and submitted to water refinement with the HADDOCK webserver v2.2 22 as previously described 23,24 , with the goal to optimize interface geometry and energetics. Briefly, ACE2/RBD heterodimers without glycans but in the presence of Zn 2+ were uploaded to the HADDOCK refinement interface and submitted with default parameters. Weighting for the sorting of structures (scoring) after water explicit refinement [25][26][27] were set for van

Figshare and SBGrid. Structure files and associated data of human ACE2 variants in complex with SARS-
CoV-2 RBD generated in this work have been deposited in Figshare 28 . The same data have also been deposited in SBGrid 29 .
Two folders are shared, (a) 6M0J for the models derived from the crystal structure 30 and (b) 6M17 for the models derived from the cryo-EM structure 11 . In each folder the following subdirectories are placed: variants, ALA_scan, and UniProt, and specifically for 6M17, an additional subdirectory is included, PyMOL_models_6M17. This directory includes .pdb and .cif files of variants which were created by considering the complete cryo-EM model, with cofactors (ions, sugars) and all interfaces. In addition, the initial .pdb files that were used to produce all reported variants are placed in each folder (6M0J_chains_AE.pdb or 6M17_chains_BE.pdb).
For the common subdirectories (variants, ALA_scan, UniProt), structure is as follows: The subdirectory variants contains data for ACE2 residue variants naturally occurring in the human population 16 , UniProt contains data for variants with in vitro mutations reported in the literature [17][18][19][20] , and ALA_scan contains data for variants resulted from the performance of computational alanine scanning mutagenesis at the interface of SARS-CoV-2 RBD and the human ACE2 receptor.
In detail, each subdirectory (variants, ALA_scan, UniProt) includes three files: the results file after the HADDOCK refinement 22 (.html file), the parameter file that was used for the structure calculation (.web), and the top scoring refined structure file (.pdb file).The user can reproduce any run by uploading the.web file using the online server (https://haddock.science.uu.nl/services/HADDOCK2.2/haddockserver-file.html).
The nomenclature of each file in subdirectories variants, ALA_scan and UniProt corresponds to XXXX_ R1NUMR2. XXXX stands for the PDB ID from which the model was extracted, R1 is the one-letter residue code of the native residue of the ACE2 receptor, NUM is the residue number according to the Uniprot sequence of human ACE2 receptor and R2 is the one-letter residue code of the variant to which the residue R1 was changed. Results of the energetic calculations with HADDOCK for each generated variant of the complex are summarized in Online-only Table 1, and Tables 1 and 2. Github. An online structure viewer of the resulting models from all refinement runs and their energetics is available at: https://kastritislab.github.io/human-ace2-variants/. The structure viewer allows the user to visualize interface contacts, compare structural information, and be informed about the corresponding energetics for any model reported in this work.

technical Validation
Data redundancy and structural mapping. Variants in the 3 datasets are distinct, showing minor overlap in terms of amino acid substitution (Fig. 2a). The computational alanine scanning shows a minor overlap with reported mutagenesis studies, where only 13% of the total mutations can be found in both datasets. In addition, only 1 out of the 39 in vitro designed ACE2 variants can be found in the human population (Fig. 2a). Mutations from missense variants are distributed across the entire ACE2 surface (Fig. 2b), including the interfaces with the SARS-CoV-2 RBD and B 0 AT1 partners 11 (Fig. 2c). This structural mapping highlights the usefulness of ACE2 variants for structure-based design, as different residues affect the physical-chemical parameters of the receptor, and consequently, its underlying affinity towards different protein-protein interactions. www.nature.com/scientificdata www.nature.com/scientificdata/ Stereochemical quality. The stereochemical quality of derived models of ACE2 variants is of equivalent quality as their template structures, since we performed mostly single amino acid substitutions and refined them using restrained molecular dynamics simulations in explicit water 26 . This protocol is well-known to improve the quality of experimental structures and docking models 26,27 . Modeling from different templates. To assess the consistency of the HADDOCK water refinement protocol, we additionally constructed homology models using the crystal structure of the ACE2 in complex with the SARS-CoV-2 RBD 30 . Although the root-mean-square deviation (RMSD) between the C α atoms of the residues from the two calculated structures is low (RMSD = 1.054 Å), we observed high variability in rotamer states, in particular for interface residues. The buried surface area (BSA, Å 2 ) of both structures is within the distribution of BSAs for transient protein-protein interactions with known affinities 24 (Fig. 3a). Interestingly, the crystallographic structure and designed variants have larger BSA as compared to the cryo-EM counterparts (Fig. 3a). This is expected since structures determined by X-ray diffraction are more tightly packed due to the crystal state of the protein. In contrast, the cryo-EM interface is smaller, likely because the specimen was captured in vitreous ice and was free in solution. In addition, model building during cryo-EM map interpretation is performed within an averaged Coulomb electrostatic potential map, which may lead to low resolution or absent densities in flexible regions and, therefore, less tight interface packing. Consistency in energy calculations. Usage of these two templates for generating variants and performing energy calculations constitutes an independent test for the robustness of the refinement protocol. Overall, for all generated models, high values for the corresponding Pearson-product momentum correlation coefficients are observed for HADDOCK score and underlying desolvation energies (Fig. 3b-c). This shows that energetic components for the HADDOCK score in both structures have similar contributions, desolvation energy being the most dominant. Only favourable energies are calculated for the variants when using the crystallographic model as an initial structure (Fig. 3c), whereas both favourable and unfavourable energies are calculated for the variants using the cryo-EM model (Fig. 3b). This is due to the presence of both transmembrane and soluble domains of the ACE2 in the cryo-EM model, whereas the crystallographic model includes only soluble domains. Desolvation energies, therefore, reflect contributions of solvation in the structures, in the presence or absence of the transmembrane regions.  www.nature.com/scientificdata www.nature.com/scientificdata/ Overlap with external datasets. To identify systematically present variations in our datasets, we overlapped the reported variations for which we communicate the respective structural models with 3 additional datasets described below: • The experimental Procko dataset (PROCKO). A recent preprint tested affinity of 2,223 ACE2 missense mutants with the RBD of the S protein of SARS-CoV-2 after one round of selection 31 . Interestingly, overlap of those data with the 3 datasets described above is minor (78 common out of 242 mutations) (Fig. 4). In particular, overlap with genome variants is even lower (20 out of 141). This highlights the complexity underlying genome variation in the human population and the distinct evolutionary pressure of the ACE2 gene as compared to in vitro deep mutagenesis experiments. Still, our structural models for the 20 overlapping mutations which have available affinity values (S19P , E23K, K26R, E37K, F40L, N64K, M82I, G326E, E329G, G352V,  H378R, M383T, Q388L, P389H, T445M, I446M, F504I, F504L, S511P, R514G) can act as a starting point for further characterization. • ACE2 mutations from cancer patients derived from COSMIC v91 32 . Due to the higher risk of severe COVID-19 symptoms manifesting in cancer patients 33 , we have specifically focused on retrieving genetic variants of ACE2 available in COSMIC v91 32 (Fig. 4). Interestingly, 15 genetic variants reported in gnomAD (R115Q, R115W, H195Y, R219H, D368N, E375D, F400L, D609N, R671Q, R708W, R710C, R716H, R716C, N720D, R768W) are also identified in cancer patients (Fig. 4, shown in bold). This result provides a hypothesis on the role of these mutations in SARS-CoV-2 infection to be further investigated. • ACE2 mutations from COVID-19 patients included in LOVD 3.0 34 (Fig. 4). LOVD 3.0 reports additional variants for the ACE2 receptor and includes the N720D mutation which has been identified as a variant in COVID-19 patients in the Italian population 35 . N720D is found in genomic data (gnomAD), cancer (COS-MIC v91) and COVID-19 patients (LOVD 3.0). Another ACE2 protein variation identified in COVID-19 patients is the K26R, which is also included in the gnomAD data, but not in cancer patients. This mutation has been successfully expressed by Procko 31 and appears to increase binding affinity for the RBD of the S protein (Fig. 4). Interestingly, our respective 3D interaction model shows one of the lowest HADDOCK scores (−108.9 ± 5.1 a.u.), strongest van der Waals interactions (−57.8 ± 5.6 kcal.mol −1 ) and most favourable desolvation energy (−11.4 ± 7.9 kcal.mol −1 ) compared to all other analyzed mutations (Online-only Table 1).
Considering the communicated correlation of HADDOCK score components with binding affinities for 144 protein-protein interactions 24 , the above-mentioned calculated energetic values corroborate the Procko results on the increased affinity for K26R, and therefore, possible higher infectivity of SARS-CoV-2. This is also corroborated by our distance calculations showing that K26R is only ~10 Å away from the interaction interface (Fig. 2).