# In Silico Analysis of High-Risk Missense Variants in Human ACE2 Gene and Susceptibility to SARS-CoV-2 Infection

Authors: Asmae Saih; Hana Baba; Meryem Bouqdayr; Hassan Ghazal; Salsabil Hamdi; Anass Kettani; Lahcen Wakrim
Journal: BioMed Research International
 (2021)

---

## Abstract

SARS-CoV-2 coronavirus uses for entry to human host cells a SARS-CoV receptor of the angiotensin-converting enzyme (ACE2) that catalyzes the conversion of angiotensin II into angiotensin (1-7). To understand the effect of ACE2 missense variants on protein structure, stability, and function, various bioinformatics tools were used including SIFT, PANTHER, PROVEAN, PolyPhen2.0, I. Mutant Suite, MUpro, SWISS-MODEL, Project HOPE, ModPred, QMEAN, ConSurf, and STRING. All twelve ACE2 nsSNPs were analyzed. Six ACE2 high-risk pathogenic nsSNPs (D427Y, R514G, R708W, R710C, R716C, and R768W) were found to be the most damaging by at least six software tools (cumulative score between 6 and 7) and exert deleterious effect on the ACE2 protein structure and likely function. Additionally, they revealed high conservation, less stability, and having a role in posttranslation modifications such a proteolytic cleavage or ADP-ribosylation. This in silico analysis provides information about functional nucleotide variants that have an impact on the ACE2 protein structure and function and therefore susceptibility to SARS-CoV-2.


---

## Body

## 1. Introduction

During December 2019 in Wuhan, China, a novel infectious disease called the coronavirus disease 2019 or COVID-19 was detected and found to be caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2). Later on, it was declared as a pandemic by the World Health Organization (WHO) on March 2020 [1, 2]. The most common symptoms were fever, cough, shortness of breath, sputum production, and fatigue [3]. SARS-CoV-2 infects the respiratory tract and engenders acute respiratory disease [3]. SARS-CoV-2 utilizes angiotensin-converting enzyme 2 (ACE2) as an entry receptor and the cellular serine protease TMPRSS2 for SARS-CoV-2 spike (S) protein priming [4].TheACE2 gene is located on the chromosome X, specifically at the position Xp22. It contains 18 exons that show a striking resemblance in exons size and organization to those of the ACE gene [5, 6]. The complete cDNA sequence of ACE2 encodes an 805 amino acid protein, which is composed mainly of an N-terminal signal peptide (17 amino acid residues), a peptidase domain (positioning from 19-615 amino acids), and a C-terminal collectrin, which plays a significant role as a regulator of renal amino acid transportation, insulin exocytosis, and β-cell proliferation [5, 7].ACE2 is expressed mostly in the vascular endothelial cells of the kidney and constitutes a vital site of blood pressure regulation in the renin-angiotensin system (RAS) and heart function. In the heart, ACE2 constitutes the first pathway for the metabolism of angiotensin II and is the most essential factor for progressive cardiac disease [8–10]. Interestingly, it seems that ACE2 has diverse biological functions, such as regulation of blood pressure by the renin-angiotensin aldosterone system (RAAS), and metabolizing angiotensin II (Ang II) to Ang (1-7), a biological active peptide that binds to the Mas receptor to exert many beneficial vasodilatory, antifibrotic, antithrombotic, and antiproliferative actions [11, 12]. Very recently, it remains possible that ACE2 gene polymorphism may influence both the susceptibility to SARS-CoV-2 infection and COVID-19 disease outcome [5]. There are a few reports related to computational analysis of missense variants on the ACE2 gene. Hence, the present study is aimed at identifying deleterious variants in the ACE2 gene by using various bioinformatics tools to extrapolate the possible associations of ACE2 polymorphisms with COVID-19 susceptibility.


## 2. Materials and Methods

### 2.1. Retrieval of Variant Datasets

The data related toACE2 SNPs (rsIDs) was retrieved from genome aggregation database (gnomAD A3) (https://gnomad.broadinstitute.org/gene/ENSG00000130234?dataset=gnomad_r3) and Ensembl database (https://www.ensembl.org/Homo_sapiens/Gene/Variation_Gene/Table?db=core;g=ENSG00000130234;r=X:15561033-15602148). The amino acid (AA) sequence was retrieved from UniProt (https://www.uniprot.org/uniprot/Q9BYF1.fasta). Information on human ACE2 protein and gene was collected from the Online Mendelian Inheritance in Man (OMIM) database (https://www.omim.org/entry/300335).


### 2.2. Prediction of Deleterious nsSNPs ofACE2 Gene

To predict the deleterious effect of nonsynonymous SNPs (nsSNPs) of the humanACE2 gene on protein function, four tools were used, namely, SIFT, PolyPhen 2.0, PANTHER, and PROVEAN. SIFT (Sorting Intolerant from Tolerant) (https://sift.bii.a-star.edu.sg/) is a program that employs sequence homology to predict the impact of an AA substitution on protein function. It classifies AA substitutions with a score less than 0.05 as deleterious [13]. PolyPhen 2.0 (Polymorphism Phenotyping v2), (http://genetics.bwh.harvard.edu/pph2/) evaluates the impact of an AA substitution on the protein function and structure. It classifies the substitution as probably damaging (score=1.0), and possibly damaging or benign (score=0.0) [14]. PANTHER cSNP (Protein ANalysis THrough Evolutionary Relationship_coding_SNP) (http://pantherdb.org/tools/csnpScoreForm.jsp) calculates the likelihood of a single AA change on protein function, and it is based on the PANTHER-PSEP (Position_Specific Evolutionary Preservation) method [15]. PROVEAN (Protein Variation Effect Analyzer) (http://provean.jcvi.org/index.php) is a web server that was used to predict the impact of AA substitutions on the biological function of a protein [16].


### 2.3. Effect on the Stability of Protein

The stability of the protein was checked using the MUpro bioinformatics tool and I. Mutant Suite. MUpro (http://mupro.proteomics.ics.uci.edu/) predicts if a mutation increases or decreases the stability of protein structure [17]. I. Mutant Suite, a support vector machine-based algorithm, is available at (http://gpcr2.biocomp.unibo.it/cgi/predictors/I-Mutant3.0/I-Mutant3.0.cgi). It predicts the mutant protein stability starting from protein sequence alone [18].


### 2.4. Phylogenetic Conservational Analysis ofACE2

Conservation prediction of theACE2 protein sequence was analyzed with ConSurf (https://consurf.tau.ac.il/), a web server used for identifying functional regions in proteins by analysing the evolutionary dynamics of AA substitutions among homologous sequences. The conservation score of 1-3 is variable, 5-6 as an intermediate scale, and 7-9 as a highly conserved AA positions [19].


### 2.5. Prediction of Posttranslational Modification Sites forACE2

ACE2 protein posttranslational modification sites were predicted using ModPred (http://www.modpred.org/), a sequence-based predictor of 23 types of posttranslational modification (PTM) sites on proteins [20].


### 2.6. Statistical Analysis and Cumulative Score Calculation of Pathogenic nsSNPs

In order to predict the high-risk pathogenic nsSNPs of humanACE2, the cumulative score for all software tools (SIFT, PolyPhen2.0, PROVEAN, PANTHER, MUpro, I. Mutant, and ModPred) was calculated by using the Sum function in Excel [21]. Then, we set a restricted cumulative score value; when the result of the seven software tools were combined, the amino acid substitution that was evaluated to be deleterious by at least 6 tools would be identified as ACE2 high-risk pathogenic nsSNPs. For the final, correlation analysis was performed using SPSS v19 software. ANOVA and Student t-tests were applied to compare the predictions of the different tools. A p value less than 0.05 was significant [22].


### 2.7. Modelling

The 3D structure of full-lengthACE2 protein was retrieved from The Research Collaboratory for Structural Bioinformatics Protein Data Bank (RCSB-PDB, http://rcsb.org), (ID 6M17) with a resolution of 2.9 Å. The 3D structure of the wild-type and mutant protein was generated using SWISS-MODEL (https://swissmodel.expasy.org/), and the quality of the model was checked using the Qualitative Model Energy Analysis (QMEAN) server (https://swissmodel.expasy.org/qmean/). SWISS-MODEL is a fully automated program that was used to predict the 3D structure of proteins. It generates 3D models by using homology modelling techniques [23]. The FASTA AA sequence of ACE2 protein was an input for the SWISS-MODEL. The predicted model of ACE2 from SWISS-MODEL was an input for the QMEAN analysis. QMEAN is a server that provides access to three scoring functions (QMEAN [24], QMEANBrane [25], and QMEANDisco) [26]. It estimates the quality of protein structure models in protein structure prediction [27]. PyMol is an open source program used for the three-dimensional visualization of macromolecules including proteins, nucleic acids, and small molecules [28]. The align command in PyMol is used to superpose two or more protein structures, and the superposition was evaluated based on RMSD (root mean square deviation) calculation [29].


### 2.8. Prediction of Structural Effect of Point Mutation onACE2

Project HOPE (Have (y) Our Protein Explained) is a next-generation web server available at “https://www3.cmbi.umcn.nl/hope/” that is used to analyze the single point mutations. HOPE collects information from data sources and produces a mutation report enriched with figures that illustrates the effects of the mutation [30].


### 2.9. Prediction of Protein Interactions

STRING (Search Tool for Recurring Instances of Neighbouring Genes) is a web server available at “https://string-db.org/” that provides a platform for searching functional associations between proteins [31].


## 2.1. Retrieval of Variant Datasets

The data related toACE2 SNPs (rsIDs) was retrieved from genome aggregation database (gnomAD A3) (https://gnomad.broadinstitute.org/gene/ENSG00000130234?dataset=gnomad_r3) and Ensembl database (https://www.ensembl.org/Homo_sapiens/Gene/Variation_Gene/Table?db=core;g=ENSG00000130234;r=X:15561033-15602148). The amino acid (AA) sequence was retrieved from UniProt (https://www.uniprot.org/uniprot/Q9BYF1.fasta). Information on human ACE2 protein and gene was collected from the Online Mendelian Inheritance in Man (OMIM) database (https://www.omim.org/entry/300335).


## 2.2. Prediction of Deleterious nsSNPs ofACE2 Gene

To predict the deleterious effect of nonsynonymous SNPs (nsSNPs) of the humanACE2 gene on protein function, four tools were used, namely, SIFT, PolyPhen 2.0, PANTHER, and PROVEAN. SIFT (Sorting Intolerant from Tolerant) (https://sift.bii.a-star.edu.sg/) is a program that employs sequence homology to predict the impact of an AA substitution on protein function. It classifies AA substitutions with a score less than 0.05 as deleterious [13]. PolyPhen 2.0 (Polymorphism Phenotyping v2), (http://genetics.bwh.harvard.edu/pph2/) evaluates the impact of an AA substitution on the protein function and structure. It classifies the substitution as probably damaging (score=1.0), and possibly damaging or benign (score=0.0) [14]. PANTHER cSNP (Protein ANalysis THrough Evolutionary Relationship_coding_SNP) (http://pantherdb.org/tools/csnpScoreForm.jsp) calculates the likelihood of a single AA change on protein function, and it is based on the PANTHER-PSEP (Position_Specific Evolutionary Preservation) method [15]. PROVEAN (Protein Variation Effect Analyzer) (http://provean.jcvi.org/index.php) is a web server that was used to predict the impact of AA substitutions on the biological function of a protein [16].


## 2.3. Effect on the Stability of Protein

The stability of the protein was checked using the MUpro bioinformatics tool and I. Mutant Suite. MUpro (http://mupro.proteomics.ics.uci.edu/) predicts if a mutation increases or decreases the stability of protein structure [17]. I. Mutant Suite, a support vector machine-based algorithm, is available at (http://gpcr2.biocomp.unibo.it/cgi/predictors/I-Mutant3.0/I-Mutant3.0.cgi). It predicts the mutant protein stability starting from protein sequence alone [18].


## 2.4. Phylogenetic Conservational Analysis ofACE2

Conservation prediction of theACE2 protein sequence was analyzed with ConSurf (https://consurf.tau.ac.il/), a web server used for identifying functional regions in proteins by analysing the evolutionary dynamics of AA substitutions among homologous sequences. The conservation score of 1-3 is variable, 5-6 as an intermediate scale, and 7-9 as a highly conserved AA positions [19].


## 2.5. Prediction of Posttranslational Modification Sites forACE2

ACE2 protein posttranslational modification sites were predicted using ModPred (http://www.modpred.org/), a sequence-based predictor of 23 types of posttranslational modification (PTM) sites on proteins [20].


## 2.6. Statistical Analysis and Cumulative Score Calculation of Pathogenic nsSNPs

In order to predict the high-risk pathogenic nsSNPs of humanACE2, the cumulative score for all software tools (SIFT, PolyPhen2.0, PROVEAN, PANTHER, MUpro, I. Mutant, and ModPred) was calculated by using the Sum function in Excel [21]. Then, we set a restricted cumulative score value; when the result of the seven software tools were combined, the amino acid substitution that was evaluated to be deleterious by at least 6 tools would be identified as ACE2 high-risk pathogenic nsSNPs. For the final, correlation analysis was performed using SPSS v19 software. ANOVA and Student t-tests were applied to compare the predictions of the different tools. A p value less than 0.05 was significant [22].


## 2.7. Modelling

The 3D structure of full-lengthACE2 protein was retrieved from The Research Collaboratory for Structural Bioinformatics Protein Data Bank (RCSB-PDB, http://rcsb.org), (ID 6M17) with a resolution of 2.9 Å. The 3D structure of the wild-type and mutant protein was generated using SWISS-MODEL (https://swissmodel.expasy.org/), and the quality of the model was checked using the Qualitative Model Energy Analysis (QMEAN) server (https://swissmodel.expasy.org/qmean/). SWISS-MODEL is a fully automated program that was used to predict the 3D structure of proteins. It generates 3D models by using homology modelling techniques [23]. The FASTA AA sequence of ACE2 protein was an input for the SWISS-MODEL. The predicted model of ACE2 from SWISS-MODEL was an input for the QMEAN analysis. QMEAN is a server that provides access to three scoring functions (QMEAN [24], QMEANBrane [25], and QMEANDisco) [26]. It estimates the quality of protein structure models in protein structure prediction [27]. PyMol is an open source program used for the three-dimensional visualization of macromolecules including proteins, nucleic acids, and small molecules [28]. The align command in PyMol is used to superpose two or more protein structures, and the superposition was evaluated based on RMSD (root mean square deviation) calculation [29].


## 2.8. Prediction of Structural Effect of Point Mutation onACE2

Project HOPE (Have (y) Our Protein Explained) is a next-generation web server available at “https://www3.cmbi.umcn.nl/hope/” that is used to analyze the single point mutations. HOPE collects information from data sources and produces a mutation report enriched with figures that illustrates the effects of the mutation [30].


## 2.9. Prediction of Protein Interactions

STRING (Search Tool for Recurring Instances of Neighbouring Genes) is a web server available at “https://string-db.org/” that provides a platform for searching functional associations between proteins [31].


## 3. Results

### 3.1. SNP Dataset

In this study, 242 missense variants were collected fromEnsembl and gnomAD A3 databases. SIFT, PolyPhen 2.0, PROVEAN, and PANTHER algorithms were used to predict the functional effects of mutation on the protein. MUpro and I. Mutant Suite tools were used to identify the mutation effects on protein stability. HOPE, ConSurf, ModPred, SWISS-MODEL, and STRING were used to predict the mutation effects on protein structure, function, and protein-protein interactions. In total, 13 different bioinformatics programs and web servers were used to assess the mutation effects on the variants in this investigation, because the prediction of the effect of a mutation using one algorithm is not sufficient for assessing that mutation effect.


### 3.2. Prediction of Functional nsSNPs inACE2

The twelve missense SNPs were predicted to be deleterious using SIFT (where tolerance index score was ranged from 0 to 0.02) (Table1). The nsSNPs predicted by SIFT were validated by PolyPhen 2.0, PANTHER, and PROVEAN. In particular, PolyPhen 2.0 results showed that twelve AA switches (R219C, R219H, M383T, P389H, D427Y, R514G, R708W, R710H, R710C, R716C, L731F, and R768W) were predicted probably damaging (score>0.96) (Table 1). The PROVEAN analysis identified seven AA substitutions (M383T, P389H, D427Y, R514G, R708W, R710C, and R768W) were scored as deleterious (score below -2.5), and the rest were noted as neutral (score above -2.5) (Table 1). PANTHER predicted that nine nsSNPs were probably damaging (R219C, R219H, M383T, P389H, R514G, R710H, R710C, L731F, and R768W), two were noted as possibly damaging (R708W, and R716C), and one was scored as probably benign (D427Y) (Table 1).Table 1
Summary of prediction results for nsSNPs in the humanACE2 gene using various bioinformatics tools.
SIFTPolyPhen 2.0PROVEANPANTHERSNPs IDAA variantScorePredictionScorePredictionScorePredictionPSPEPredictionrs372272603R219C0Del1.000Pro Dam-2.483Not Del750Pro Damrs759590772R219H0Del0.988Pro Dam-1.520Not Del750Pro Damrs1396769231M383T0Del1.000Pro Dam-4.797Del911Pro Damrs762890235P389H0Del0.993Pro Dam-7.862Del1037Pro Damrs1316056737D427Y0Del0.970Pro Dam-3.390Del176Not Delrs1352194082R514G0.02Del0.989Pro Dam-6.483Del1037Pro Damrs776995986R708W0.01Del1.000Pro Dam-3.105Del361Pos Delrs901495523R710C0Del1.000Pro Dam-2.936Del750Pro Damrs370187012R710H0Del1.000Pro Dam-1.788Not Del750Pro Damrs144869363R716C0.01Del0.975Pro Dam-1.638Not Del220Pos Delrs759590772L731F0Del0.995Pro Dam-1.124Not Del750Pro Damrs372272603R768W0Del1.000Pro Dam-2.822Del750Pro DamAA variant: amino acid variant; Del: deleterious; Not Del: not deleterious; Pro Dam: probably damaging; Pos Dam: possibly damaging; PSPE: position-specific evolutionary preservation.


### 3.3. Effect on the Stability of Protein

MUpro and I. Mutant Suite were used to predict change in protein stability. The result of I. Mutant Suite showed that eight amino acid substitutions (R219H, D427Y, R514G, R708W, R710C, R710H, R716C, and R768W) were recorded as decreasing the stability of theACE2 protein, with a reliability index value ranged between 0 and 9, while the mutations (R219C, M383T, P389H, and L731F) were predicted to increase the stability of the ACE2 protein (Table 2). The result of stability predicted by MUpro showed that all the twelve nsSNPs (R219C, R219H, M383T, P389H, D427Y, R514G, R708W, R710H, R710C, R716C, L731F, and R768W) have decreased the stability of the ACE2 protein (Table 3).Table 2
Prediction of effect of nsSNPs on protein stability using I. Mutant Suite.
SNPs IDAmino acid substitutionStabilityRIDDG value prediction (kcal/mol)rs372272603R219CIncrease0-0.50rs759590772R219HDecrease3-0.58rs1396769231M383TIncrease2-0.32rs762890235P389HIncrease20.20rs1316056737D427YDecrease3-0.56rs1352194082R514GDecrease8-1.84rs776995986R708WDecrease4-0.68rs901495523R710CDecrease8-1.72rs370187012R710HDecrease9-1.91rs144869363R716CDecrease7-2.20rs759590772L731FIncrease40.01rs372272603R768WDecrease0-0.55RI: reliability index; DDG: the free energy change value.Table 3
Prediction of effect of nsSNPs onACE2 protein stability using MUpro.
Method 1: SVM (support vector machine)Method 2: neural networkMutationDelta GPredictionConfidence scoreEffectConfidence scoreEffectR219C-0.89Decrease-0.082Decrease-0.57DecreaseR219H-1.16Decrease-0.089Decrease-0.51DecreaseM383T-1.56Decrease-0.94Decrease-0.70DecreaseP389H-0.45Decrease-0.55Decrease-0.85DecreaseD427Y-0.79Decrease0.16Increase-0.67DecreaseR514G-1.15Decrease-0.76Decrease-0.98DecreaseR708W-0.99Decrease-0.88Decrease-0.98DecreaseR710C-1.17Decrease-0.83Decrease-0.87DecreaseR710H-1.66Decrease-0.07Decrease-0.98DecreaseR716C-1.34Decrease-0.10Decrease-0.82DecreaseL731F-0.87Decrease-0.04Decrease-0.66DecreaseR768W-1.89Decrease-1Decrease-0.99DecreaseDelta G: energy change value; Confidence score: (>0: increase the stability, <0: decrease the stability).


### 3.4. Phylogenetic Conservation

According to ConSurf analysis, R514, R708, R710, and R768 are highly conserved residues with conservation score equal to 9 and also predicted to be exposed and functional. M389 and L731 are slightly conserved residues and buried. R219, P389, and D427 are variable residues and exposed, and finally one residue (R716) is average (conservation score equal to 5). Results of ConSurf prediction ofACE2 SNPs are summarized in Figure 1.Figure 1
ConSurf analysis ofACE2_HUMAN angiotensin-converting enzyme 2 (UniProt ID: Q9BYF1).


### 3.5. Prediction of Posttranslational Modification (PTM) Sites

ModPred was used to predict the effect of nsSNPs on posttranslational modification (PTM) process of the humanACE2 protein. ModPred identified sites for methylation (R219), ADP-ribosylation (R219, R708, and R768), and Proteolytic cleavage (P389, D427, R708, R710, R716, and R768). The results of PTM sites prediction are shown in Table 4.Table 4
ModPred analysis for posttranslational modification sites (PTMs) prediction.
ResidueModificationScoreConfidenceR219Proteolytic cleavage0.70MediumADP-ribosylation0.67MediumMethylation0.76MediumP389Proteolytic cleavage0.52LowD427Proteolytic cleavage0.73MediumR708Proteolytic cleavage0.76MediumADP-ribosylation0.67MediumR710Proteolytic cleavage0.54LowR716Proteolytic cleavage0.77MediumR768Proteolytic cleavage0.84MediumADP-ribosylation0.66Low


### 3.6. Statistical Analysis and Cumulative Score Calculation of Pathogenic nsSNPs

As a result of combining seven algorithms, the amino acid substitution, namely, R708W, was scored as the most deleterious nonsynonymous SNP with a cumulative score of 7 by all the seven tools, while 6 (R219H, D427Y, R514G, R710C, R716C, and R768W), 4 (R219C, M383T, P389H, and R710H), and 1 (L731F) variant got a cumulative score of 6, 5, and 4, respectively (Table5). Among these variants, six nsSNPs (D427Y, R514G, R708W, R710C, and R768W) are newly evaluated as ACE2 high-risk pathogenic nsSNPs and were selected for further investigation. The prediction from SIFT, PROVEAN, PolyPhen 2.0, PANTHER, ModPred, I. Mutant, and MUpro were shown to be significant with a p value equal to 6.1308E-29 of ANOVA test and correlated (Supplementary figure 1). Student t-test results between the software tools were significant with a p value less than 0.0001, suggesting that the selected algorithms are accurate enough to evaluate the pathogenicity of these nsSNPs.Table 5
Cumulative score calculation of tested nsSNPs.
MutationSIFTPolyPhenPROVEANPANTHERI.MutantMUproModPredCumulative scoreR219C11010115R219H11011116M383T11110105P389H11110105D427Y∗∗11101116R514G∗∗11111106R708W∗11111117R710C∗∗11111106R710H11011105R716C∗∗11011116L731F11010104R768W∗∗11111106∗The highly pathogenic nsSNPs by all the seven tools, ∗∗the highly pathogenic nsSNPs by at least 6 tools, 0 indicates neutral prediction, 1 indicates deleterious prediction.


### 3.7. 3D Modelling and Biophysical Validation ofACE2

The crystal structure of theACE2 protein was retrieved from the Protein Data Bank (PDB ID: 6 M17, resolution at 2.9 Å) in complex with SARS-CoV-2 receptor-binding domain and sodium-dependent neutral amino acid transporter B(0)AT1. The structure was used as a template for the comparative modelling of mutant structures through SWISS-MODEL by submitting the FASTA AA sequence of ACE2 (percent identity equal to 100). The QMEAN server was used to predict the quality of the models, where the global QMEAN scores were 0.83±0.05. This indicated that the predicted models were of good quality (Table 6). The six tested mutants, namely, D427Y, R514G, R708W, R710C, R716C, and R768W, were subjected to PyMol using “Align command” and compared to the native structure with regard to conformational variations by calculating the RMSD (root mean square deviation) (Table 7 and Figure 2). Project HOPE was used to identify the structural effects of mutations of interest. Results of project-HOPE of ACE2 SNPs were displayed in Table 8.Table 6
Data generated from SWISS-MODEL.
SNPsQMEANCBeta interaction energyAll-atom pairwise energySolvation energyTorsion energyD427Y-1.091.081.050.76-1.52R514G-1.270.390.780.52-1.50R708W-1.280.330.700.52-1.50R710C-1.370.340.740.42-1.55R716C-1.330.420.760.45-1.54R768W-1.091.011.050.74-1.50CBeta interaction energy: distance-dependent potential using CBeta atoms as interaction center; All-atom pairwise energy: assessment of long-range interactions; Solvation energy: description of the burial status of the residues; Torsion energy: analysis of the local backbone geometry [32].Table 7
Structural alignment comparing mutant and wild-typeACE2 models.
PositionVariantAlignRMSD (Å)427D427Y6 m17.1.E0.010514R514G6 m17.1.E0.112708R708W6 m17.1.E0.112710R710C6 m17.1.E0.112716R716C6 m17.1.E0.112768R768W6 m17.1.E0.010Figure 2
Superimposed structure of native and mutant models of theACE2 protein. The superimposed structure of native amino acid aspartic acid (blue color) with mutant amino acid tyrosine (light blue color) at position 427. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid glycine (red color) at position 514. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid tryptophan (pale green color) at position 708. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid cysteine (green cyan color) at position 710. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid cysteine (violet color) at position 716. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid tryptophan (light orange color) at position 768.Table 8
ACE2 protein phenotype feature prediction by HOPE analysis.
ResidueStructurePropertiesD427Y(i) The WT is predicted to be located in its preferred secondary structure, a turn the mutant prefers to be in another secondary structure, therefore the local conformation will be slightly destabilized(ii) Mutation of the WT into none has the following effect SPD-PSN: slightly inhibits interaction with SARS-CoV spike gp(iii) The mutation is possibly damaging to the proteinR514G(i) The mutation introduces a glycine at this position. Glycine is very flexible and can disturb the required rigidity of the protein(ii) Residues in the vicinity of the mutated residue are annotated in the UniProt as being a binding site(iii) The mutation could affect the local structure and as a consequence affect the binding site(iv) The mutation is possibly not damaging to the protein(v) The mutant is smaller than the WT; this might lead to loss of interactionR708W(i) The mutant is more hydrophobic than the WT(ii) The WT charge was positive, and the mutant charge is neutral(iii) The mutant is bigger; this might lead to bumpsR710C(i) The mutant is more hydrophobic than the WT(ii) The WT is very conserved(iii) The mutant is smaller; this might lead to loss of interactionR716C(i) The mutant is smaller than the WT(ii) The mutation is located within a stretch of residues annotated in UniProt as a special region: essential of cleavage by TMPRSS11D and TMPRSS2(iii) The difference in amino acid properties can disturb this region and its functionR768W(i) The WT is very conserved. In some rare cases, the mutation might occur without damaging the protein(ii) The mutant is bigger than the WT; this might lead to bumps(iii) The mutation introduces a more hydrophobic residue at this position; this can result in loss of hydrogen bonds and or disturb correct folding


### 3.8. Analysis of Protein-Protein Interaction

STRING prediction showed that ACE2 interacts with angiotensin II receptor type 1 (AGTR1), angiotensin II receptor type 2 (AGTR2), prolyl carboxypeptidase (PRCP), renin (REN), angiotensinogen (AGT), membrane metalloendopeptidae (MME), dipeptidyl peptidase 4 (DPP4), meprin A subunit alpha (MEP1A), meprin A subunit beta (MEP1B), and X-prolyl aminopeptidase 2 (XPNPEP2) (Figure 3).Figure 3
STRING network analysis of the ACE2 gene.


## 3.1. SNP Dataset

In this study, 242 missense variants were collected fromEnsembl and gnomAD A3 databases. SIFT, PolyPhen 2.0, PROVEAN, and PANTHER algorithms were used to predict the functional effects of mutation on the protein. MUpro and I. Mutant Suite tools were used to identify the mutation effects on protein stability. HOPE, ConSurf, ModPred, SWISS-MODEL, and STRING were used to predict the mutation effects on protein structure, function, and protein-protein interactions. In total, 13 different bioinformatics programs and web servers were used to assess the mutation effects on the variants in this investigation, because the prediction of the effect of a mutation using one algorithm is not sufficient for assessing that mutation effect.


## 3.2. Prediction of Functional nsSNPs inACE2

The twelve missense SNPs were predicted to be deleterious using SIFT (where tolerance index score was ranged from 0 to 0.02) (Table1). The nsSNPs predicted by SIFT were validated by PolyPhen 2.0, PANTHER, and PROVEAN. In particular, PolyPhen 2.0 results showed that twelve AA switches (R219C, R219H, M383T, P389H, D427Y, R514G, R708W, R710H, R710C, R716C, L731F, and R768W) were predicted probably damaging (score>0.96) (Table 1). The PROVEAN analysis identified seven AA substitutions (M383T, P389H, D427Y, R514G, R708W, R710C, and R768W) were scored as deleterious (score below -2.5), and the rest were noted as neutral (score above -2.5) (Table 1). PANTHER predicted that nine nsSNPs were probably damaging (R219C, R219H, M383T, P389H, R514G, R710H, R710C, L731F, and R768W), two were noted as possibly damaging (R708W, and R716C), and one was scored as probably benign (D427Y) (Table 1).Table 1
Summary of prediction results for nsSNPs in the humanACE2 gene using various bioinformatics tools.
SIFTPolyPhen 2.0PROVEANPANTHERSNPs IDAA variantScorePredictionScorePredictionScorePredictionPSPEPredictionrs372272603R219C0Del1.000Pro Dam-2.483Not Del750Pro Damrs759590772R219H0Del0.988Pro Dam-1.520Not Del750Pro Damrs1396769231M383T0Del1.000Pro Dam-4.797Del911Pro Damrs762890235P389H0Del0.993Pro Dam-7.862Del1037Pro Damrs1316056737D427Y0Del0.970Pro Dam-3.390Del176Not Delrs1352194082R514G0.02Del0.989Pro Dam-6.483Del1037Pro Damrs776995986R708W0.01Del1.000Pro Dam-3.105Del361Pos Delrs901495523R710C0Del1.000Pro Dam-2.936Del750Pro Damrs370187012R710H0Del1.000Pro Dam-1.788Not Del750Pro Damrs144869363R716C0.01Del0.975Pro Dam-1.638Not Del220Pos Delrs759590772L731F0Del0.995Pro Dam-1.124Not Del750Pro Damrs372272603R768W0Del1.000Pro Dam-2.822Del750Pro DamAA variant: amino acid variant; Del: deleterious; Not Del: not deleterious; Pro Dam: probably damaging; Pos Dam: possibly damaging; PSPE: position-specific evolutionary preservation.


## 3.3. Effect on the Stability of Protein

MUpro and I. Mutant Suite were used to predict change in protein stability. The result of I. Mutant Suite showed that eight amino acid substitutions (R219H, D427Y, R514G, R708W, R710C, R710H, R716C, and R768W) were recorded as decreasing the stability of theACE2 protein, with a reliability index value ranged between 0 and 9, while the mutations (R219C, M383T, P389H, and L731F) were predicted to increase the stability of the ACE2 protein (Table 2). The result of stability predicted by MUpro showed that all the twelve nsSNPs (R219C, R219H, M383T, P389H, D427Y, R514G, R708W, R710H, R710C, R716C, L731F, and R768W) have decreased the stability of the ACE2 protein (Table 3).Table 2
Prediction of effect of nsSNPs on protein stability using I. Mutant Suite.
SNPs IDAmino acid substitutionStabilityRIDDG value prediction (kcal/mol)rs372272603R219CIncrease0-0.50rs759590772R219HDecrease3-0.58rs1396769231M383TIncrease2-0.32rs762890235P389HIncrease20.20rs1316056737D427YDecrease3-0.56rs1352194082R514GDecrease8-1.84rs776995986R708WDecrease4-0.68rs901495523R710CDecrease8-1.72rs370187012R710HDecrease9-1.91rs144869363R716CDecrease7-2.20rs759590772L731FIncrease40.01rs372272603R768WDecrease0-0.55RI: reliability index; DDG: the free energy change value.Table 3
Prediction of effect of nsSNPs onACE2 protein stability using MUpro.
Method 1: SVM (support vector machine)Method 2: neural networkMutationDelta GPredictionConfidence scoreEffectConfidence scoreEffectR219C-0.89Decrease-0.082Decrease-0.57DecreaseR219H-1.16Decrease-0.089Decrease-0.51DecreaseM383T-1.56Decrease-0.94Decrease-0.70DecreaseP389H-0.45Decrease-0.55Decrease-0.85DecreaseD427Y-0.79Decrease0.16Increase-0.67DecreaseR514G-1.15Decrease-0.76Decrease-0.98DecreaseR708W-0.99Decrease-0.88Decrease-0.98DecreaseR710C-1.17Decrease-0.83Decrease-0.87DecreaseR710H-1.66Decrease-0.07Decrease-0.98DecreaseR716C-1.34Decrease-0.10Decrease-0.82DecreaseL731F-0.87Decrease-0.04Decrease-0.66DecreaseR768W-1.89Decrease-1Decrease-0.99DecreaseDelta G: energy change value; Confidence score: (>0: increase the stability, <0: decrease the stability).


## 3.4. Phylogenetic Conservation

According to ConSurf analysis, R514, R708, R710, and R768 are highly conserved residues with conservation score equal to 9 and also predicted to be exposed and functional. M389 and L731 are slightly conserved residues and buried. R219, P389, and D427 are variable residues and exposed, and finally one residue (R716) is average (conservation score equal to 5). Results of ConSurf prediction ofACE2 SNPs are summarized in Figure 1.Figure 1
ConSurf analysis ofACE2_HUMAN angiotensin-converting enzyme 2 (UniProt ID: Q9BYF1).


## 3.5. Prediction of Posttranslational Modification (PTM) Sites

ModPred was used to predict the effect of nsSNPs on posttranslational modification (PTM) process of the humanACE2 protein. ModPred identified sites for methylation (R219), ADP-ribosylation (R219, R708, and R768), and Proteolytic cleavage (P389, D427, R708, R710, R716, and R768). The results of PTM sites prediction are shown in Table 4.Table 4
ModPred analysis for posttranslational modification sites (PTMs) prediction.
ResidueModificationScoreConfidenceR219Proteolytic cleavage0.70MediumADP-ribosylation0.67MediumMethylation0.76MediumP389Proteolytic cleavage0.52LowD427Proteolytic cleavage0.73MediumR708Proteolytic cleavage0.76MediumADP-ribosylation0.67MediumR710Proteolytic cleavage0.54LowR716Proteolytic cleavage0.77MediumR768Proteolytic cleavage0.84MediumADP-ribosylation0.66Low


## 3.6. Statistical Analysis and Cumulative Score Calculation of Pathogenic nsSNPs

As a result of combining seven algorithms, the amino acid substitution, namely, R708W, was scored as the most deleterious nonsynonymous SNP with a cumulative score of 7 by all the seven tools, while 6 (R219H, D427Y, R514G, R710C, R716C, and R768W), 4 (R219C, M383T, P389H, and R710H), and 1 (L731F) variant got a cumulative score of 6, 5, and 4, respectively (Table5). Among these variants, six nsSNPs (D427Y, R514G, R708W, R710C, and R768W) are newly evaluated as ACE2 high-risk pathogenic nsSNPs and were selected for further investigation. The prediction from SIFT, PROVEAN, PolyPhen 2.0, PANTHER, ModPred, I. Mutant, and MUpro were shown to be significant with a p value equal to 6.1308E-29 of ANOVA test and correlated (Supplementary figure 1). Student t-test results between the software tools were significant with a p value less than 0.0001, suggesting that the selected algorithms are accurate enough to evaluate the pathogenicity of these nsSNPs.Table 5
Cumulative score calculation of tested nsSNPs.
MutationSIFTPolyPhenPROVEANPANTHERI.MutantMUproModPredCumulative scoreR219C11010115R219H11011116M383T11110105P389H11110105D427Y∗∗11101116R514G∗∗11111106R708W∗11111117R710C∗∗11111106R710H11011105R716C∗∗11011116L731F11010104R768W∗∗11111106∗The highly pathogenic nsSNPs by all the seven tools, ∗∗the highly pathogenic nsSNPs by at least 6 tools, 0 indicates neutral prediction, 1 indicates deleterious prediction.


## 3.7. 3D Modelling and Biophysical Validation ofACE2

The crystal structure of theACE2 protein was retrieved from the Protein Data Bank (PDB ID: 6 M17, resolution at 2.9 Å) in complex with SARS-CoV-2 receptor-binding domain and sodium-dependent neutral amino acid transporter B(0)AT1. The structure was used as a template for the comparative modelling of mutant structures through SWISS-MODEL by submitting the FASTA AA sequence of ACE2 (percent identity equal to 100). The QMEAN server was used to predict the quality of the models, where the global QMEAN scores were 0.83±0.05. This indicated that the predicted models were of good quality (Table 6). The six tested mutants, namely, D427Y, R514G, R708W, R710C, R716C, and R768W, were subjected to PyMol using “Align command” and compared to the native structure with regard to conformational variations by calculating the RMSD (root mean square deviation) (Table 7 and Figure 2). Project HOPE was used to identify the structural effects of mutations of interest. Results of project-HOPE of ACE2 SNPs were displayed in Table 8.Table 6
Data generated from SWISS-MODEL.
SNPsQMEANCBeta interaction energyAll-atom pairwise energySolvation energyTorsion energyD427Y-1.091.081.050.76-1.52R514G-1.270.390.780.52-1.50R708W-1.280.330.700.52-1.50R710C-1.370.340.740.42-1.55R716C-1.330.420.760.45-1.54R768W-1.091.011.050.74-1.50CBeta interaction energy: distance-dependent potential using CBeta atoms as interaction center; All-atom pairwise energy: assessment of long-range interactions; Solvation energy: description of the burial status of the residues; Torsion energy: analysis of the local backbone geometry [32].Table 7
Structural alignment comparing mutant and wild-typeACE2 models.
PositionVariantAlignRMSD (Å)427D427Y6 m17.1.E0.010514R514G6 m17.1.E0.112708R708W6 m17.1.E0.112710R710C6 m17.1.E0.112716R716C6 m17.1.E0.112768R768W6 m17.1.E0.010Figure 2
Superimposed structure of native and mutant models of theACE2 protein. The superimposed structure of native amino acid aspartic acid (blue color) with mutant amino acid tyrosine (light blue color) at position 427. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid glycine (red color) at position 514. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid tryptophan (pale green color) at position 708. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid cysteine (green cyan color) at position 710. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid cysteine (violet color) at position 716. The superimposed structure of native amino acid arginine (blue color) with mutant amino acid tryptophan (light orange color) at position 768.Table 8
ACE2 protein phenotype feature prediction by HOPE analysis.
ResidueStructurePropertiesD427Y(i) The WT is predicted to be located in its preferred secondary structure, a turn the mutant prefers to be in another secondary structure, therefore the local conformation will be slightly destabilized(ii) Mutation of the WT into none has the following effect SPD-PSN: slightly inhibits interaction with SARS-CoV spike gp(iii) The mutation is possibly damaging to the proteinR514G(i) The mutation introduces a glycine at this position. Glycine is very flexible and can disturb the required rigidity of the protein(ii) Residues in the vicinity of the mutated residue are annotated in the UniProt as being a binding site(iii) The mutation could affect the local structure and as a consequence affect the binding site(iv) The mutation is possibly not damaging to the protein(v) The mutant is smaller than the WT; this might lead to loss of interactionR708W(i) The mutant is more hydrophobic than the WT(ii) The WT charge was positive, and the mutant charge is neutral(iii) The mutant is bigger; this might lead to bumpsR710C(i) The mutant is more hydrophobic than the WT(ii) The WT is very conserved(iii) The mutant is smaller; this might lead to loss of interactionR716C(i) The mutant is smaller than the WT(ii) The mutation is located within a stretch of residues annotated in UniProt as a special region: essential of cleavage by TMPRSS11D and TMPRSS2(iii) The difference in amino acid properties can disturb this region and its functionR768W(i) The WT is very conserved. In some rare cases, the mutation might occur without damaging the protein(ii) The mutant is bigger than the WT; this might lead to bumps(iii) The mutation introduces a more hydrophobic residue at this position; this can result in loss of hydrogen bonds and or disturb correct folding


## 3.8. Analysis of Protein-Protein Interaction

STRING prediction showed that ACE2 interacts with angiotensin II receptor type 1 (AGTR1), angiotensin II receptor type 2 (AGTR2), prolyl carboxypeptidase (PRCP), renin (REN), angiotensinogen (AGT), membrane metalloendopeptidae (MME), dipeptidyl peptidase 4 (DPP4), meprin A subunit alpha (MEP1A), meprin A subunit beta (MEP1B), and X-prolyl aminopeptidase 2 (XPNPEP2) (Figure 3).Figure 3
STRING network analysis of the ACE2 gene.


## 4. Discussion

Angiotensin-converting enzyme 2 (ACE2) plays an important role in the renin-angiotensin aldosterone system by metabolizing angiotensin II to angiotensin (1-7) [11]. This important enzyme was identified as a functional receptor for the severe acute respiratory syndrome coronavirus (SARS-CoV) and the novel coronavirus (SARS-CoV-2) [11]. Some studies suggest that the receptor-binding domains of SARS-CoV-2-S and SARS-CoV-S bind with identical affinities to ACE2 [25]. It has been shown that ACE2 genomic variants may play a key role in susceptibilities to COVID-19 [33].The SARS-CoV-2 starts its infection by binding toACE2 via its receptor-binding domain (RDB) [34]. Recently, it has been reported that the SARS-CoV-2 uses two mechanisms of host cell entry: the first mechanism on the ACE2 mediated virus endocytosis by using the clathrin- and caveolae-dependent pathways. The second one is dependent on the transmembrane serine protease 2- (TMPRSS2-) mediated membrane fusion [34].The study by Shang et al. 2020 [35], Walls et al. [36], Wrapp et al. 2020, and Yan et al. 2020 have demonstrated that there are important ACE2 residues, namely, S19, Q24, T27, F28, D30, K31, H34, E35, E37, R357, E329, N330, K353, G354, D355, R357, P389, and R393, that play crucial roles in mediating the SARS-CoV-2 spike protein-ACE2 interaction. Other results were broadly in line with ACE2 variants that are identified to alter the virus host interaction and consequently alter susceptibility to SARS-CoV-2 [37]. Some of them, like S19P, I21V, E23K, K26R, T27A, N64K, T92I, Q102P, and H378R, are predicted to increase susceptibility and therefore render individuals more susceptible to the SARS-CoV-2 [37]. Other ACE2 variants such K31R, N35I, H34R, E35K, E37K, D38V, Y50F, N51S, M62V, K68E, F72V, Y83H, G326H, G352V, D355N, Q388L, and D509Y are predicted to decrease susceptibility and render individuals more resistant to the SARS-CoV-2 [37]. A number of 61 deleterious variants in the ACE2 gene including R219C, R219H, M383T, P389H, D427Y, R514G, R708W, R710H, R710C, R716C, L731F, and R768W, have been reported to influence susceptibility to COVID-19 (Hou et al. 2020). To confirm these results, the twelve variants, namely, R219C, R219H, M383T, P389H, D427Y, R514G, R708W, R710H, R710C, R716C, L731F, and R768W, were tested. This in silico analysis might be helpful in understanding the effect of missense variants on protein structure, function, and stability of ACE2 in relation with COVID-19.A 2020 study by Fahd Al-Mulla et al. has shown that the arginine residues at 708 and 716 positions play an important role inACE2 cleavage by TMPRSS2 and TMPRSS11D, whereas mutations on R708 and R716 seem to reduce directly ACE2 cleavage by TMPRSS2 [32]. Another line of research (Hossein Lanjanian et al. 2021) suggests that the arginine residue at 710 position of the ACE2 receptor plays a crucial role in mediating ACE2-TMPRSS2 interaction by involving hydrogens bonds [38]. A second study achieved by Behrooz Darbani (2020) has identified the rs1316056737 (R514G) as an interaction inhibitor variant that might impact the interaction between the ACE2 receptor and the viral spike protein [39].In the current study, various algorithms, namely, SIFT, PolyPhen 2.0, PROVEAN, PANTHER, I. Mutant Suite, MUpro, and ConSurf, were used to identify the most deleterious nsSNPs of theACE2 gene. In our in silico analysis, we identified six nsSNPs (D427Y, R514G, R708W, R710C, R716C, and R768W) from twelve nsSNPs. The six nsSNPs were predicted to be deleterious by at least six algorithms, decreased the stability of ACE2 by both MUpro and I. Mutant, and located in conserved regions with a score conservation of 9, expected D427Y, which may affect the structure and function of ACE2 protein.ModPred was used to determine the posttranslational modification (PTM) sites ofACE2 identified R710 and R716 as PTM sites for proteolytic cleavage and R708 for ADP-ribosylation and proteolytic cleavage. Consequently, mutations at R708, R710, and R716 might affect PTM of the ACE2 gene.The structural deviations between the native and modelled mutants were analyzed using the RMSD by measuring the average distance between the atoms of the superimposed proteins [40]. It has also been shown that RMSD values greater than 0.15 were evaluated to be significant and would have an impact on protein function and structure [40]. Therefore, the six models R427Y, R514, R708W, R710C, R716C, and R768W showed lower RMSD values, i.e., 0.010, 0.112, 0.112, 0.112, 0.112, and 0.010, respectively, which indicate minimal structural dissimilarity between the native and mutant models of ACE2.According to HOPE, the AA substitution (R514G) introduces a glycine (flexible) and could disturb the required rigidity of theACE2 protein and could affect the binding site where the mutation being located, while the mutation R716C is located in a region essential for cleavage by TMPRSS11D and TMPRSS2 proteases, and the difference in AA properties can disturb this region and its function.Furthermore, by using theSTRING server, ACE2 protein interacts with ten various proteins, namely, AGTR1, AGTR2, PRCP, REN, AGT, MME, DPP4, MEP1A, MEP1B, and XPNPEP2. The majority of those proteins play a crucial role in the regulation of blood pressure and kidney pathways, suggesting the association of ACE2 in kidney diseases. All the ten various proteins had a score ranging from 0.991 to 0.858, suggesting that proteins are partially biologically connected with an PPI enrichment p value equal to 1.29E-09. The proteins AGT, REN, AGTR1, and AGTR2 were found to have an important role in RAS signalling pathways such as regulation of blood pressure, body fluid and electrolyte homeostasis, and sodium retention by the kidney. The interaction analysis using the STRING server suggested that the ACE2 protein not only participate in hypertension and cardiovascular diseases but also a candidate in COVID-19 infection.However, these findings showed that D427Y, R514G, R708W, R710C, R716C, and R768W, are the structurally and functionally most significant variants in the humanACE2 gene. (Table 7). It has been found that the residues ARG 708, 710, and 716, are located in the dimeric interface of human ACE2, and they are found to be important for its cleavage by the protease TMPRSS2, and this processing is indispensable for augmentation of SARS-S-driven entry into host cells [9]. The limitation of this in silico analysis is that the six deleterious variants should be confirmed with future extensive experiments and clinical wet lab approaches to figure out the mechanism of these mutations in susceptibility to COVID-19.


## 5. Conclusion

Our presentin silico analysis of nsSNPs of the human ACE2 gene concluded that the mutations D427Y, R514G, R708W, R710C, R710H, and R716C are the most deleterious nsSNPs among the reported ACE2 gene variants. All six nsSNPs were predicted to be damaging and also affecting conserved positions, protein stability, and the posttranslational modification sites. Therefore, these mutations would likely affect the function on ACE2 with regard to host susceptibility to SARS-CoV-2 infection. Hence, the results of this study confirm previous findings and may be helpful for further understanding the role of these six ACE2 nsSNPs in susceptibility to COVID-19.



---

*Source: 6685840-2021-04-09.xml*