Structural analysis of ACE2 variant N720D demonstrates a higher binding affinity to TMPRSS2.

AIMS: Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), a novel member of the betacoronaviruses family affecting the lower respiratory tract mainly through binding to angiotensin converting enzyme 2 (ACE2) via its S-protein. Genetic analysis of (ACE2) gene revealed several variants that have been suggested to regulate the interaction with S protein. This study investigates the N720D variant, positioned in the collectrin-like domain (CLD) at proximity to type II transmembrane serine protease (TMPRSS2) cleavage site. MAIN
METHODS: The effect of N720D variant on ACE2 structure and thermodynamic stability was studied by DynaMut. HDOCK was utilised to model TMPRSS2 protease binding to ACE2 WT and D720 variant cleavage site. PRODIGY was used to calculate binding affinities and MD simulation tools calculated the at 100 ns for ACE2 apo structure and the ACE2-TMPRSS2 complex. KEY
FINDINGS: The N720D variant is a more dynamic structure with a free energy change (ΔΔG): -0.470 kcal/mol. As such, introducing a tighter binding affinity of Kd = 3.2 × 10-10 M between TMPRSS2 and N720D variant. RMSD, RMSF calculations showed the N720D variant is less stable, however, RMSF values of the D720-TMPRSS2 complex reflected a slower dynamic motion. SIGNIFICANCE: The hotspot N720D variant in the CLD of ACE2 affected the stability and flexibility of ACE2 by increasing the level of motion in the loop region, resulting in a more favourable site for TMPRSS2 binding and cleavage. Consequently, this would facilitate S-protein binding and can potentially increase viral entry highlighting the importance of variants affecting the ACE2-TMPRSS2 complex.

Introduction

Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), similar to the earlier SARS-CoV and Middle East respiratory syndrome coronavirus (MERS-CoV), is associated with severe lower respiratory tract disease symptoms such as severe pneumonia and bronchiolitis [36,50]. SARS-CoV-2 is part of the betacoronaviruses family and has a long ORF1ab polyprotein at the 5′ end, followed by structural proteins, glycosylated spike (S), an envelope protein (E), membrane (M), and nucleocapsid (N) proteins [39]. SARS-CoV-2 S-protein binds to angiotensin-converting enzyme 2 (ACE2), facilitating the invasion of the host cells [24,38,43]. Primarily ACE2 is found in the lower respiratory tract [17], which may explain the severe respiratory syndrome associated with SARS-CoV-2 infection [12,48].

ACE2 is an 805 amino acid protease that exists in a monomer, dimer equilibrium [19,42], and is well-known for its role in hypertension. ACE2 exerts its functions through cleaving either Angiotensin I or Angiotensin II into the inactive peptides Angiotensin (1–9) and Angiotensin (1–7), respectively. Angiotensin (1–9) gets further metabolized into Angiotensin (1–7). Angiotensin (1–7) is a vasodilator; hence, ACE2 counteracts the vasoconstrictor effects of ACE-Angiotensin II axis [7,11,33,34,40]. ACE2 consists of an N-terminal peptidase domain (PD), a C-terminal collectrin-like domain (CLD), and a transmembrane helix [9,42,49]. S-protein binding to ACE2 dimer occurs at the receptor-binding domain (RBD) when the S1 subunit interacts with the PD of ACE2 [13,42]. The S-protein binding affinity of SARS-CoV-2 to ACE2 is 10–20-fold higher than SARS-CoV [24], which may explain the high virulence and the increased rate of human to human transmission. Active ACE2 is cleaved on the N-terminus by either metalloprotease ADAM17 (residues 652 to 659) or type II transmembrane serine protease TMPRSS2 (residues 697 to 716) [14]. TMPRSS2 cleavage is required for ACE2 to interact with the S-protein [14]. This was further corroborated on wild type (WT), and knockout (KO) TMPRSS2 in mice infected with SARS-CoV since TMPRSS2 cleavage resulted in an enhanced viral entry [16].

A recent large-scale investigation of ACE2 gene coding sequences variants (1700 variants) and allele frequencies (AF) between different populations reported 32 different variants, including seven hotspot variants [5]. These missense allelic variants included K26R, I468V, A627V, N638S, S692P, N720D, and L731I/F [5,22]. The most prominent mutations on the ACE2 gene were K26R and N720D, where both mutations presented a higher frequency in the European population in comparison to the Chinese and Middle Eastern populations [2,5,22,27]. Cao et al. suggested that these genetic variants were not significantly different between Caucasian and Asian populations but can potentially alter viral binding and entry to host cells [5]. Whereas, Al-Mulla et al. showed that the prevalence of the most ACE2 activating variant N720D was much higher among Europeans (2.5%), when compared to Iranians (0.6%), Kuwaitis (0.3%), Qataris (0.2%) and other global populations (0.4%) and MAF of these variants associated with higher death rate from SARS-CoV-2 infection in the European population compared to the Middle Eastern population [2].

In addition to its high allelic frequency, the significance of the N720D variant emerges from its location in the CLD domain on a dynamic loop region 4 amino acids away from the TMPRSS2 cleavage site [2,4,27]. Since TMPRSS2 cleavage of ACE2 is an essential step in SARS-CoV-2 S-protein binding, any structural changes caused by the N720D variant might affect this process. The aim of this study is to elucidate the impact of structural changes caused by ACE2 N720D variant on TMPRSS2 cleavage and, in turn, S-protein binding. To achieve these objectives, we used computational structural biology tools to predict the effect of the N720D variant on the stability and flexibility of ACE2 structure. In addition, we modelled the binding affinities of TMPRSS2 to ACE2 WT and variant.

Results and discussion

ACE2 N720D mutation

Both K26R on the PD and N720D on the CLD domains have the highest frequency of mutations on ACE2 gene in the European population in comparison to the Chinese and Middle Eastern populations. The PD of ACE2 binds to S-protein with a high affinity of ~15 nM Kd [38]; therefore, a properly functioning ACE2 is essential for viral cell entry. As such, studies have observed the effect K26R and other variants on the structure of the PD that might affect the S-protein binding, and subsequently, the SARS-CoV-2 entry into cells [15,22,27]. Whereas, the N720D variant is located on the same interface of the loop region in close vicinity to the TMPRSS2 cleavage site. Having the mutation in the loop region close to the cleavage site may result in conformational dynamics that affect the binding affinity of TMPRSS2 to ACE2. Here we analysed the structure of the recently published Cryo-EM ACE2 dimer with a 2.9 Å resolution bound to neutral amino acid transporter B0AT1 (PDB ID:6M18) (Fig. 1A) [42]. In our thermodynamic and binding analysis, we modelled out the B0AT1 since, Yan et al. stated that B0AT1 is only found in the kidneys and intestine and not the lungs, in addition, it might hinder access of TMPRSS2 to the cleavage site [42].

Fig. 1

The effect of N720D mutation on the ACE2 protein structure. A) ACE2 dimer complex promoters 1 and 2 (blue and white, respectively) showing the PD site where the SARS-CoV-2 S protein RBD binds. The CLD region is also shown where ACE2 is cleaved by TMRSS2. B) The single domain of the ACE2 mutation N720D where the red region of the protein depicts the more flexible region of the protein due to the N720D mutation with a (ΔΔG): −0.470 kcal/mol and ΔΔS: 0.070 kcal·mol−1·K−1. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Conformational dynamics of N720D

Stability analysis of ACE2 showed that the variant N720D is more dynamic in structure with a free energy change (ΔΔG): −0.470 kcal/mol and entropy of (ΔΔSVib: 0.070 kcal·mol−1·K−1) (Fig. 1B). Hussain et al. measured the free energy change for a series of notable mutations on the PD of ACE2. Most mutations studied on the PD, predicted a less stable ACE2 structure, which in turn may increase the susceptibility of the individual to SARS-CoV-2 [4,22]. The variant K26R, which presented the highest allelic frequency, resulted in a less stable ACE2 structure at 25 °C (ΔΔG −0.34 kcal/mol) and 35 °C (ΔΔG −0.30 kcal/mol) [15]. However, the N720D mutation is positioned in the loop region preceded by the TMPRSS2 cleavage site between residues 697 and 716, which are the third and fourth helices in the CLD of ACE2 (Fig. 1A). Protein loop regions have higher motional dynamics (femtoseconds to seconds) due to the limited covalent interactions in the region [18]. As a result of the mutation, a decrease in ΔΔG and an increase in entropy in the system triggers local unfolding, hence exposing cleavage sites to the surrounding moiety [30] to be readily attacked by proteases such as TMPRSS2.

ACE2 TMPRSS2 binding

The loop region close to the cleavage site of ACE2 D720 variant displayed a less stable and more dynamic structure. Therefore, we modelled the binding of TMPRSS2 to the cleavage site of ACE2 N720 and D720 variant. The arginine and lysine residues within ACE2 TMPRSS2 cleavage site are essential for cleavage and required for augmentation of SARS-S-driven entry [14]. The ACE2 residues R697 (70.025 Å), K702 (108.684 Å) and R705 (106.398 Å) on the cleavage site to TMPRSS2 HDOCK analysis (Fig. 2A) [44,45]. Whereas, the binding affinities were predicted by PRODIGY using the ACE2-TMPRSS2 binding model (Fig. 2B).

Fig. 2

A) Modelled TMPRSS2 bound to ACE2 TMPRSS2 cleavage site. B) The dissociation constant Kd for TMPRSS2 bound to ACE2 N720 (purple) and D720 (green). A tighter binding affinity of Kd = 3.2 × 10−10 M between TMPRSS2 and D720 variant than N720 WT. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

The binding energy scores (Kd) of ACE2-TMPRSS2 complex (WT and D720 variant) were predicted as a function of temperature [41]. As the temperature increase from 25 °C to 40 °C, there was a 2-fold increase in Kd score, indicating a decrease in binding affinity for both ACE2 WT and D720 variant TMPRSS2 complexes (Fig. 2B). At 25 °C the D720-TMPRSS2 complex showed a tighter binding affinity (3.2 × 10−10 M) in comparison to the WT complex (4.9 × 10−10 M), which is indicative of a more favourable binding between TMPRSS2 and ACE2 D720 variant. Furthermore, the predicted Kd for D720 was less affected as the temperature increased from 25 °C to 40 °C with a ΔKd 9.0 × 10−10 M, in contrast to WT with a ΔKd 1.4 × 10−9 M. Such an increase in thermal energy disrupts the noncovalent interactions between ACE2 and TMPRSS2 resulting in a lower binding affinity [10]. Despite that, the decrease in binding affinity across the temperatures was less apparent for the ACE2 D720 variant complex with TMPRSS2. As such, this tight binding supports the increased susceptibility [21] of ACE2 N720D variant to TMPRSS2 cleavage.

Molecular dynamic simulations of ACE2 bound to TMPRSS2

To understand the conformational and dynamic features of N720D variant on the stability of ACE2 and the effect on binding to TMPRSS2, molecular dynamics (MD) simulation is an imperative method to explore the behavior of each system in real-time. In this study, 100 ns simulations were run to compare the dynamic behavior of all the systems, including the WT N720-Apo, D720-Apo, N720-TMPRSS2, and D720-TMPRSS2 where the RMSD, RMSF using AMBER package.

The stability of all the four systems were calculated by Root Mean Square Deviation (RMSD). Fig. 3A depicts the RMSD of the apo structures, with the average RMSD for both WT N720 (blue) and D720 (red) variant was 3.0 Å during the 100 ns simulation. From 0 to 38 ns the WT N720 showed a more rigid structure with a higher RMSD than D720. However, after 40 ns the RMSD value decreased to be 2.5 Å till the end of simulation with no significant convergence. Thus, stability remains uniform for the rest of the time observed during the 100 ns simulation indicating the protein is relatively stable [1,3]. On the other hand, the D720 variant system (Red) appeared relatively lower than N720 WT for the first 38 ns. However, onwards the systems gained equilibrium [25], and higher fluctuations were observed at 55 ns and 90-97 ns indicating significant changes in the conformation of D720 variant. Such major convergence at these two-time points is the results of the D720 being introduced in the system as such, causing instability. Therefore, the variation in conformation may be one of the factors contributing to the destabilization of the protein, making them close to the TMPRSS2 site more accessible cleavage.

Fig. 3

Molecular dynamic simulations RMSD plots at 100 ns simulations of ACE2 N720 and D720 (A) and ACE2 N720-TMPRSSS2 and D720-TMPRSS2 complex (B). RMSF plot of ACE2 N720 and D720 (C) and ACE2 N720-TMPRSSS2 and D720-TMPRSS2 complex (D).

As for the effect of N720D mutation on stability and dynamics of ACE2-TMPRSS2 complex, RMSD values are shown in Fig. 3B. The WT N720-TMPRSS2 complex remained stable during the simulation with an average RMSD of 1 Å, with the system showing slight convergence at 22 ns and from 35 to 38 ns. Whereas, the D720-TMPRSS2 complex showed a very dynamic behavior, with average RMSD being higher than the WT. The system starts to converge at 30–35 ns; then a major convergence is observed at 60–80 ns with the RMSD continually increasing. Thus, these results signify that the mutation has produced a direct effect on the stability and binding of ACE2. As the RMSD of the mutant complex is still increasing with major convergence at intervals. This is probably due to the integrated mutation at position 720 which causes changes in the conformational dynamics of the proteins and ultimately affecting the interaction with the TMPRSS2. Longer simulation is required for the system to enter the production phase and to better understand the impact of the substitution.

To demonstrate residual flexibility, the root mean square fluctuations (RMSF) were calculated at 100 ns simulations for both apo and complex structures Fig. 3C and D. The N720 WT (blue) shows relatively lower residual flexibility in comparison to D720 variant. The D720 variant (red) shows higher fluctuation, specifically in the region between 320 and 350, 410–420, 605–618, and 695–710. This indicates the mutation has affected the internal dynamics of the protein and causes the residues to be more flexible compared to the N720 WT [1]. However, in the complex system, the N720-TMPRSS2 (purple) showed higher molecular mobility than the variant D720-TMRSS2 complex. The lower mobility observed close to the loop region between residues 695–720 could be the result of the higher binding affinity of TMPRSS2 to the D720 variant.

Thus, the RMSF displayed a higher degree of motion in loop regions and the TMPRSS2 binding site for mutant D720 (Fig. 3C). Loop regions located opposite or close to the cleavage site regulate the dynamic and structural changes of enzyme-protein binding. Consequently, any conformational changes occurring opposite the binding site can in turn, influence the enzyme cleavage process [47]. Therefore, higher RMSF exhibited by the N720D mutation caused higher conformational mobility in the loop region close to the cleavage site, which can result in a more favourable TMPRSS2-ACE2 binding interaction through either conformational selection mechanism or induced-fit mechanism [32,47].

Conclusion

The recent variants observed on the ACE2 gene in different populations have gained attention as they might affect the SARS-CoV-2 binding. Since the ACE2 PD bound to S-protein of SARS-CoV-2 has been observed by Cryo-EM studies [38,42], most studies have concentrated on variants on the PD [4,15,22]. Whereas, the effect of N720D variant positioned on the CLD of ACE2 structure has been less examined.

In our current study, we demonstrated the potential importance of the N720D mutation on the stability and flexibility of ACE2. Since ACE2-TMPRSS2 complex structure is not experimentally determined, we modelled the effect of N720D on TMPRSS2-ACE2 interaction. The D720 variant was characterized by an increase in free energy that affected protein stability and flexibility. This was reflected by higher RMSD and RMSF indicating an increased instability due to conformational changes, which elevates the level conformational mobility in the loop region and TMPRSS2 binding site. With the variant D720, ACE2-TMPRSS2 complex depicted a tighter binding affinity, reflected by lower mobility close to the loop region. The effect of the variant D720 on the structural dynamics of ACE2 was emphasized by the calculated Kd for the ACE2-TMPRSS2 complex. The D720 ACE2 variant showed a tighter binding affinity in comparison to ACE2 WT. The more favourable TMPRSS2-ACE2 binding interaction supports the increased susceptibility of ACE2 N720D variant to TMPRSS2 cleavage.

Our report is based on modelling analysis of ACE2 N720D variant to develop an understanding of its interaction with TMPRSS2 and the consequent of SARS-CoV-2 binding. Our finding gives an insight into the potential reason for which N720D was the most frequent variant linked to the higher death rates in Europe. Further computational and functional studies are required to elucidate the clinical importance of ACE2 N720D variant.

In conclusion, based on our modelling analysis of ACE2 variant N720D, we postulate that this variant impacted the stability and flexibility of ACE2, resulting in a more favourable site for TMPRSS2 binding and cleavage. As a result, carriers of this variant have an increased S-protein binding and can potentially have increased viral entry caused by enhanced interaction between ACE2 and TMPRSS2.

Methods

The structures of ACE2-B0AT1 (PDB ID:6M18) [38] was used for thermodynamic and structural analysis in this report. The structure of TMPRSS2 was modelled from SWISS-Model [37] and I-TASSER [46]. DynaMut [28] webserver was used to predict the effect of genetic variants on the stability and flexibility on N720D on the ACE2 protein.

Protein-protein docking of ACE2 (PDB ID:6M18) WT(N720) and variant (D720) and TMPRSS2 was done by HDOCK server [44,45], which is based on a hybrid algorithm of template-based modelling and ab initio free docking. The ACE2 residues R697 (70.025 Å), K702 (108.684 Å), and R705 (106.398 Å) on the cleave site to TMPRSS2 residues H296, D345, S441, and D435. Model 1 with the lowest docking energy score and the highest ligand root-mean-square deviation (RMSD) was selected to analyse the binding energy scores (Kd) using PRODIGY server [41]. PRODIGY is a robust predictive system that utilises structural properties of protein-protein interactions, the number of interfacial contacts and non-interacting surfaces to calculate proteins binding affinity [35].

To understand the dynamics and interacting behavior, both the apo N720 and D720 ACE2 and ACE2-TMPRSS2 complexes were subjected to molecular dynamics simulation. Amber 18 package with AMBER ff14SB force field was used to execute the simulations [6,26]. Systems were solvated with TIP3P water box with 18.0 Å distance each side and were neutralized by adding Na+ ions. We used 300 K temperature and 1.0 bar pressure of 1.0 bar using Langevin thermostat and Berendsen Barostat controllers [8,23]. For hydrogen long-range interaction SHAKE algorithm and particle mesh Ewald summation (PME) algorithm was used [20,31]. For the non-bonded interactions 10.0 Å cut-off was fixed. A two-steps gentle minimization followed by heating and equilibration was performed. A 100 ns simulation was carried out at the NPT ensemble. The Cartesian coordinates were stored at every 10 ps. We obtained 5000 frames from each simulation. Systems stability and residual flexibility was also calculated using CPPTRAJ and PTRAJ [29]. For stability RMSD while for flexibility RMSF was calculated.