Bexultan Kazybay1, Ashfaq Ahmad2, Chenglin Mu3, Diana Mengdesh1, Yingqiu Xie4. 1. Biology Department, Nazarbayev University, Nur-Sultan, 010000, Kazakhstan. 2. Department of Bioinformatics, Hazara University, Mansehra, 21300, Pakistan. 3. Zhejiang University, Hangzhou Global Scientific and Technological Innovation Center, Hangzhou, China. 4. Biology Department, Nazarbayev University, Nur-Sultan, 010000, Kazakhstan. Electronic address: xieautumnus@yahoo.com.
Abstract
Variants of SARS-CoV-2 lineages including the most recently circulated Omicron, and previous pandemic B.1.351, B.1.1.7, which have been public concerns, contain a N501Y mutation located in the spike receptor binding domain. However, the potential interactions with host cells linking N501Y mutation to pathogenic relevance remain elusive. Recently, we and others report that kinases such as PI3K/AKT signaling are essential in SARS-CoV-2 entry. Here we analyzed the predicted potential kinases interacting with the mutation. Bioinformatics tools including structure-prediction based molecular docking analysis were applied. We found kinases such as EGFR might potentially act as new factors involving the N501Y mutation binding through possible phosphorylation at Y501 and enhanced affinity in certain variants. To our surprise, the Omicron receptor binding domain harboring N501Y mutation did not enhance binding to EGFR which might be due to the mutations of charged polar to uncharged polar side chains located on the interaction interfaces. Similarly, potent gains of phosphorylation in B.1.351 and B.1.1.7 by mutations were predicted and interaction networks were analyzed with enrichment of pathways. Given kinases might be elevated in cancer patients, the N501Y mutation containing lineages may be possibly much more infectious and additional care for cancer management might be taken into consideration by precision prevention, therapy or recovery.
Variants of SARS-CoV-2 lineages including the most recently circulated Omicron, and previous pandemic B.1.351, B.1.1.7, which have been public concerns, contain a N501Y mutation located in the spike receptor binding domain. However, the potential interactions with host cells linking N501Y mutation to pathogenic relevance remain elusive. Recently, we and others report that kinases such as PI3K/AKT signaling are essential in SARS-CoV-2 entry. Here we analyzed the predicted potential kinases interacting with the mutation. Bioinformatics tools including structure-prediction based molecular docking analysis were applied. We found kinases such as EGFR might potentially act as new factors involving the N501Y mutation binding through possible phosphorylation at Y501 and enhanced affinity in certain variants. To our surprise, the Omicron receptor binding domain harboring N501Y mutation did not enhance binding to EGFR which might be due to the mutations of charged polar to uncharged polar side chains located on the interaction interfaces. Similarly, potent gains of phosphorylation in B.1.351 and B.1.1.7 by mutations were predicted and interaction networks were analyzed with enrichment of pathways. Given kinases might be elevated in cancer patients, the N501Y mutation containing lineages may be possibly much more infectious and additional care for cancer management might be taken into consideration by precision prevention, therapy or recovery.
Recently, Omicron (B.1.1.529) was reported or predicted to be associated with rapid increase of COVID-19 cases (https://www.who.int/news/item/28-11-2021-update-on-omicron). The new, widely reported and WHO concerned SARS-CoV-2 variant has caused many countries travel restrictions or border closures even more studies are still ongoing or urgently needed for the clinical investigations and border closure benefit [1]. The number of reported Omicron cases from early discovery dates to nowadays showed that in the short period about less than 2 weeks, the circulating of Omicron has been found in at least 65 countries and territories with more than thousands of confirmed cases and more than a thousand suspected cases (Fig. 1
) (https://www.gisaid.org/hcov19-variants/; https://en.wikipedia.org/wiki/SARS-CoV-2_Omicron_variant).
Fig. 1
Superimposition of the spike RBD and ACE2 complexes. Docked complexes of the Spike RBD wild type (pink) and N501Y (light green) are shown in cartoon representation with hACE2 receptor. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)
Superimposition of the spike RBD and ACE2 complexes. Docked complexes of the Spike RBD wild type (pink) and N501Y (light green) are shown in cartoon representation with hACE2 receptor. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)Omicron contains many more mutations compared with other variants such as SARS-CoV-2501.V2 (B.1.351 lineage) and VOC 202012/01 (B.1.1.7 lineage) [2]. All the three variants retain N501Y mutation on spike protein receptor binding domain (RBD) which is required for binding to human receptor ACE2 (hACE2) for viral entry to host cells. The B.1.1.7 is reported to be more infectious with higher spreading rate and increased binding affinity to ACE2 [2] in particular due to N501Y, compared with other mutations such as E484K [3,4]. The N501Y mutation localizes at RBD that might help in attaining higher binding to host cells and result elevated transmission and infectious. However, not much are reported as an evidence of kinase or other signaling pathways resulting the severe complications or increased death rate by variant containing N501Y mutation has been concluded. Nevertheless, further in-depth studies on the newly emerged variant mutations induced host cell signaling is highly needed to aid the drug development and travel infection prevention.It has been reported that tyrosine (Y) residue is easily phosphorylated [2]. Therefore, we hypothesized that the N501Y mutation could be more prone to the host-kinases to be modified or phosphorylated. Here we aimed to investigate the potential of the phosphorylation and predict kinases with molecular docking analysis for binding to kinases using first reported, wild type SARS-CoV-2 RBD as control and repurpose the potential therapeutic methods.
Methods
To test the hypothesis and predict the tyrosine kinases that might phosphorylate the mutated Spike RBD at N501Y, several bioinformatics tools are applied. Group-based Prediction System of version 5.0 web server (GPS 5.0) [5] is based on the highest probability of phosphorylation indicated by the “score”, and the “cut-off” value, which should not exceed 15% to give significant results. Another is NetPhos 3.1 web server [6,7] - DTU (www.cbs.dtu.dk › services › NetPhos) with the threshold score which is equal to 0.5.For modeling and docking assignments, wild type hACE2, coronavirus spike RBD and kinase domain of human EGFR were retrieved from the Protein data bank (7EKG and 4I23). AlphaFold and Modeller were employed to generate mutant structures for Omicron and N501Y RBD [8,9]. In order to generate protein-protein complexes, the HDock program was utilized that samples the potent binding modes by Fast Fourier transform (FFT)–based global search optimization method and considering the available homologue information from the PDB in generating binding modes [10,11].For the generated complexes, contact based prediction methods were used to predict and evaluate the binding-affinities(ΔG) and dissociation constant Kd (M), explained in Ref. [12]. For protein-protein interaction networks analysis, we applied a STRING database [13] and retrieved it by Cytoscape7.3 software [14]. ClusterONE plug-in of Cytoscape 7.3 was used to identify the most significant cluster (nodes as 45; density as 0.498; p-value <0.001) [15,16]. The mutation data were obtained from resources [17,18]. The most important KEGG pathways were identified by the data taken from the Kobas database. Five KEGG pathways were chosen according to the corrected p-value ≤3.40e-9 and kappa score 0.4. Analysis was also performed by an online tool [19].
Results and discussion
Docking analysis of a kinase with N501Y mutation in SARS-CoV-2 lineages
First we applied NetPhos 3.1 server for predicting kinases that potent potentially phosphorylate Y501, and we found the phosphorylation site was predicted whereas the gain of phosphorylation at N501 (wild type) was not found (Supplementary Fig. S1). SRC and epidermal growth factor receptor (EGFR) were potent kinases to phosphorylate Y501. Through another server GPS 5.0, EGFR and Fibroblast growth factor receptor (FGFR) were found. Interestingly, MET was predicted too even with a reduced score. It has been shown that both EGFR and MET kinases play an essential role in lung cancer [20]. On the contrary, we did not find kinases using other tools like Quokka, KinasePhos and ELM which might be explained by the limited database or different algorithm used. Thus, through the two programs prediction, EGFR potentially interacts with and phosphorylates SARS-CoV-2 spike RBD at N501Y.Next, we employed docking tools to generate desired complexes of the EGFR with wild and N501Y spike RBD. First, we calculated the binding affinities of the wild type PDB complex (7EKG) and dock complex of the mutant (N501Y) spike RBD against the hACE2 receptor as a control. Our results depicted that in both cases the complex retrieved same binding affinity of −12.6 kcal/mol, suggesting a negligible impact of the point mutation (Fig. 1, Table 1
). Subsequently, dock complexes of kinase domain of EGFR-RBD wild and mutant were generated. The N501Y EGFR-RBD complex showed slightly higher affinity of −12.3 kcal/mol compared to −11.6 kcal/mol with wild type RBD (Fig. 2, Fig. 3
, Table 1). Interestingly, the Omicron RBD showed the lowest binding affinity of −10.4 kcal/mol among all the tested complexes (Fig. 4
). Omicron RBD contains 11 mutations, where 4 of them K417 N, E484A, G496S and N501Y were found on the interface region (Fig. 4). Besides, the disappearance of the acidic and basic residues, with the interface of the Omicron RBD presented an interface region with higher uncharged polar side chains retaining amino acids.
Table 1
Statistics of docking analysis for all the generated complexes.
Complex Name
Docking Score
Template Used
LG-score/Quality
ACE2-RBD (wild)
−253.22
7EKH
5.051
ACE2-RBD (N501Y)
−271.33
7EKH
5.051
EGFR-RBD (wild)
−267.34
4C02
5.284
EGFR-RBD (N501Y)
−257.35
4C02
5.284
EGFR-Omicron RBD
−266.98
4C02
5.037
LG-score indicates the reliability of the docking results and was calculated by ProQ program. Based on the LG-score, all the docked complexes were ranked in top hierarchy (LG-score > 5.0) [11].
Fig. 2
Complex structure of the spike RBD against the kinase domain of EGFR. The predicted complex of spike RBD (pink) with kinase domain of EGFR (blue) are shown in cartoon representation. The stick representation are the interface residues located in the region of 5 Å. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)
Fig. 3
Complex structure of the mutant spike RBD against the kinase domain of EGFR. The predicted complex of mutant spike RBD (pink) with kinase domain of EGFR (blue) are shown in cartoon representation. The residues in stick representation are of the interface located in the region of 5 Å. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)
Fig. 4
Complex structure of the Omicron spike RBD against the kinase domain of EGFR. The predicted complex with green and red spheres indicate the mutations observed in Omicron RBD and on the interface region, respectively. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)
Statistics of docking analysis for all the generated complexes.LG-score indicates the reliability of the docking results and was calculated by ProQ program. Based on the LG-score, all the docked complexes were ranked in top hierarchy (LG-score > 5.0) [11].Complex structure of the spike RBD against the kinase domain of EGFR. The predicted complex of spike RBD (pink) with kinase domain of EGFR (blue) are shown in cartoon representation. The stick representation are the interface residues located in the region of 5 Å. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)Complex structure of the mutant spike RBD against the kinase domain of EGFR. The predicted complex of mutant spike RBD (pink) with kinase domain of EGFR (blue) are shown in cartoon representation. The residues in stick representation are of the interface located in the region of 5 Å. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)Complex structure of the Omicron spike RBD against the kinase domain of EGFR. The predicted complex with green and red spheres indicate the mutations observed in Omicron RBD and on the interface region, respectively. (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.)Overall, these results indicate that RBD retains the capability to interact EGFR, as the wild type and mutant complex showed binding affinity almost equivalent to the control hACE2-RBD complex.
Network prediction of the mutations in SARS-CoV-2
In addition to N501Y, we also performed prediction of other mutations of B.1.351 and B.1.1.7 lineages by potent phosphorylation using servers search as mentioned above. These include gain of possible phosphorylation by mutations of the two lineages with mutations of I2230T, N501Y; or H2799Y, D4527Y, N501Y to gain phosphorylation, compared with wild type sequences. All potent kinases were listed and further analyzed for gain of potent phosphorylation with protein-protein interaction maps (Supplementary Figs. S2 and S3). We found EGFR, ErbB and MAPK pathways are top ranked in the variants host cell signaling based on KEGG analysis (Fig. 5
A). Gene Ontology analysis on biological processes suggest the phosphorylation, cell proliferation, apoptosis and DNA replication pathways are involved in the two variants (Fig. 5B). Thus, the gain of potent phosphorylation by kinases might induce the possibly signaling transduction events resembling cancer induced cell death or phosphorylation events.
Fig. 5
Signaling analysis of potent network kinases and suggested precision therapy or prevention. KEGG pathways (A) and Gene Ontology (GO) biological process (B) were shown. The top five KEGG pathways were identified by Kobas database and conducted by OriginLab 2020b software and bioinformatics tools online. C. Suggested precision therapy, and prevention medicine based on the role of EGFR/FGFR as an example of targeting potent virus-host cell interactions to disrupt potent membrane fusion for precision prevention infection especially in kinase elevated cancer patients.
Signaling analysis of potent network kinases and suggested precision therapy or prevention. KEGG pathways (A) and Gene Ontology (GO) biological process (B) were shown. The top five KEGG pathways were identified by Kobas database and conducted by OriginLab 2020b software and bioinformatics tools online. C. Suggested precision therapy, and prevention medicine based on the role of EGFR/FGFR as an example of targeting potent virus-host cell interactions to disrupt potent membrane fusion for precision prevention infection especially in kinase elevated cancer patients.
Potent kinase inhibitors in treatment and prevention of the N501Y mutations induced gain of phosphorylation
The EGFR, a receptor tyrosine kinase, is known to be frequently mutated and upregulated in some cancers in response to their internal and external stresses to survive and resist against different types of therapies. EGFR has also been reported acting as a receptor in viral entry for influenza, and hepatitis C viruses [21]. According to Ref. [22] SARS-CoV-2 infection causes activation of aberrant STAT pathway which results in an acute lung injury upregulating EGFR that leads to constitutive STAT3 activation. Moreover, EGFR signaling can disrupt antiviral events induced by INF-α, thus promoting severity of viral infection [23]. A recent study by proteomics found that the SARS-CoV-2 proteins are mostly phosphorylated and network analysis identified EGFR as a central hint of the proteomic interactome [21]. The study also reported that sorafenib, a RAF inhibitor and RO5126766, a RAF/MEK dual inhibitor can disrupt viral entry and infection which is promising in anti-COVID-19 [21]. For mechanistic analysis, EGFR likely plays essential roles in cell fusion during infection of respiratory syncytial virus predicted by the other virus studies.Given the specificity of over-expression of EGFR or SRC kinases in some cancers and its predicted role in COVID-19, it suggests an additional consideration with the treatment of cancer patients, as cancer patients with elevated EGFR or SRC might be susceptible to the variants of SARS-CoV-2 as predicted above. Our data also suggest that the mutated N501Y and other sites of the RBD has the potential to be phosphorylated by different kinases including EGFR. Thus kinase inhibitors could possibly be considered and studied as an agent for precision treatment of cancer patients with SARS-CoV-2 infection after experimental validation. For example, erlotinib, as EGFR inhibitor, could be recommended because its combination with INF-α showed antiviral activity [24].Moreover, we may apply some herbal food which contains EGFR inhibitors compound of natural products for preventing the new variants induced COVID-19 (Fig. 5C). For example, we may select Genistein rich food to prevent the infection, including different species of soybean, black soybean, and other herbs or herbal food analyzed by a database of TCMSP Version 2.3 [[25], [26], [27]]. While processing of natural products of herbal food or traditional Chinese Medicine, possibly generates nanozyme of phosphatase which can block the kinases activity by combination, as we reported, the herbal food resources of kinase inhibitors may benefit the prevention of the variants. For example, Huangjing, which is already recommended as a component of anti-COVID-19 herbal food by clinical agencies [28,29], could be recommended for prevention of infection of variant. Moreover, integration with the traditional Chinese Medicine Huangjing, and EGFR inhibitor for treatment by applying the phosphatase activity of nanocale Huangjing of nanozyme, as we recently proposed in cancer cells [27], would also be promising. Thus, cancer patients may benefit the integration method if undergoing cancer precision therapy, with the target inhibitors as repurposed or rehabilitation medicine during Omicron circulation [30].Recently, we reported that PI3K/AKT signal pathway is essential in COVID-19 and capivasertib as an inhibitor restricts SARS-CoV-2 entry into host cells [31,32]. In addition, based on our hypothesis mentioned above, to further investigate the potential of combination therapy by EGFR and AKT inhibitors would be promising in anti-variants. Thus the crosstalk pathways would also be deserved to clarify.Given kinases are usually elevated in cancer patients, the potent kinase associated new variant mutation binding may induce more cancer signaling and sever infection, even there is no clear evidence that the SARS-CoV-2 can promote tumorigenesis. In addition, cancer patients should be noticed that the potent kinase associated infection by spike RBD binding to cell membrane receptor kinases. Thus the potent risk may exist that cancer patient might be more vulnerable for infection due to elevated kinase expression levels.Most recently, while we are revising the manuscript, a study was published on which N501Y mutation's superiority among different mutations tested for transmission [33]. The tested mutations include spike mutants of T716I, N501Y, A570D, P681H, D1118H, 982A, deletion 69–70, deletion 145, with lineage B.1.1.7 variant spike [33]. The experimental data by infected hamsters suggest that in addition to deletion 69–70, N501Y showed much more advantages of fitness compared with wild type in replication efficiency in the upper airway, and shedding in nasal cavity or secretions which are related to air-borne transmission [33]. It has been shown that N501Y alone can exert the function as a powerful driver for elevated transmission but not inducing immune disruption suggesting the vaccine escape is minimum [[33], [34]].
Conclusion
In summary, though our prediction and hypothesis are limited which needs experimental investigation, we provided potential new mutations mediated interactions and its predicted phosphorylation landscape with host cell kinases based on current evidence from published reports in database. Our results depict that EGFR receptor is one of the potential binary partner of the spike RBD that binds almost with equal affinity as RBD-hACE2 complex. Keeping the fact of EGFR kinase activity, this report suggests the possibility of phosphorylating N501Y mutation that is located on the interface region to further potentiate its importance. We also suggest a prevention approach, and applicable solutions to fight against the new variants in the current pandemic era. Cancer patients might consider the kinase elevation induced binding to spike protein for increased risk of infection and applying kinase inhibitors or natural products for both precision cancer therapy and preventive or rehabilitative medicine. Potent host cell kinases elevation as biomarker for prevention of viral infection might be considered for cancer patients during Omicron circulation.
Data and code availability
Not applicable.
Ethics declarations
Not applicable.
Ethics approval and consent to participate
Not applicable.
Consent for publication
Not applicable.
Funding
This research is not supported specifically by any funding related to bioinformatics and COVID-19. However, the authors received a cancer research funding during the work period: Faculty Development Competitive Research Grants of (16797152).
Author information
All authors agreed to the current form.
Contributions
Bexultan Kazybay: results of supplementary Figures of network and old version draft of structure prediction, pathways, and draft writing; Ashfaq Ahmad: bioinformatics analysis, new docking data including omicron docking and methods, results writing with editing; Chenglin Mu: bioinformatics analysis, and part of methods writing; Diana Mengdesh: data gathering of kinase prediction and results with part of draft writing; Yingqiu Xie: designing as supervisor, data analysis and draft writing.
Authors: Joachim Lupberger; François H T Duong; Isabel Fofana; Laetitia Zona; Fei Xiao; Christine Thumann; Sarah C Durand; Patrick Pessaux; Mirjam B Zeisel; Markus H Heim; Thomas F Baumert Journal: Hepatology Date: 2013-09-05 Impact factor: 17.425
Authors: Minkyung Baek; Frank DiMaio; Ivan Anishchenko; Justas Dauparas; Sergey Ovchinnikov; Gyu Rie Lee; Jue Wang; Qian Cong; Lisa N Kinch; R Dustin Schaeffer; Claudia Millán; Hahnbeom Park; Carson Adams; Caleb R Glassman; Andy DeGiovanni; Jose H Pereira; Andria V Rodrigues; Alberdina A van Dijk; Ana C Ebrecht; Diederik J Opperman; Theo Sagmeister; Christoph Buhlheller; Tea Pavkov-Keller; Manoj K Rathinaswamy; Udit Dalwadi; Calvin K Yip; John E Burke; K Christopher Garcia; Nick V Grishin; Paul D Adams; Randy J Read; David Baker Journal: Science Date: 2021-07-15 Impact factor: 47.728
Authors: John Jumper; Richard Evans; Alexander Pritzel; Tim Green; Michael Figurnov; Olaf Ronneberger; Kathryn Tunyasuvunakool; Russ Bates; Augustin Žídek; Anna Potapenko; Alex Bridgland; Clemens Meyer; Simon A A Kohl; Andrew J Ballard; Andrew Cowie; Bernardino Romera-Paredes; Stanislav Nikolov; Rishub Jain; Demis Hassabis; Jonas Adler; Trevor Back; Stig Petersen; David Reiman; Ellen Clancy; Michal Zielinski; Martin Steinegger; Michalina Pacholska; Tamas Berghammer; Sebastian Bodenstein; David Silver; Oriol Vinyals; Andrew W Senior; Koray Kavukcuoglu; Pushmeet Kohli Journal: Nature Date: 2021-07-15 Impact factor: 49.962