key: cord-0787984-hfyffqsv authors: Madanagopal, Premnath; Muthukumar, Harshini; Thiruvengadam, Kothai title: Computational study and design of effective siRNAs to silence structural proteins associated genes of Indian SARS-CoV-2 strains date: 2022-04-29 journal: Comput Biol Chem DOI: 10.1016/j.compbiolchem.2022.107687 sha: 59a9de60bab8dc653fd69d628e5b1860a55a2413 doc_id: 787984 cord_uid: hfyffqsv SARS-CoV-2 is a highly transmissible and pathogenic coronavirus that first emerged in late 2019 and has since triggered a pandemic of acute respiratory disease named COVID-19 which poses a significant threat to all public health institutions in the absence of specific antiviral treatment. Since the outbreak began in March 2020, India has reported 4.77 lakh Coronavirus deaths, according to the World Health Organization (WHO). The innate RNA interference (RNAi) pathway, on the other hand, allows for the development of nucleic acid-based antiviral drugs in which complementary small interfering RNAs (siRNAs) mediate the post-transcriptional gene silencing (PTGS) of target mRNA. Therefore, in this current study, the potential of RNAi was harnessed to construct siRNA molecules that target the consensus regions of specific structural proteins associated genes of SARS-CoV-2, such as the envelope protein gene (E), membrane protein gene (M), nucleocapsid phosphoprotein gene (N), and surface glycoprotein gene (S) which are important for the viral pathogenesis. Conserved sequences of 811 SARS-CoV-2 strains from around India were collected to design 21 nucleotides long siRNA duplex based on various computational algorithms and parameters targeting E, M, N and S genes. The proposed siRNA molecules possessed sufficient nucleotide-based and other features for effective gene silencing and BLAST results revealed that siRNAs' targets have no significant matches across the whole human genome. Hence, siRNAs were found to have no off-target effects on the genome, ruling out the possibility of off-target silencing. Finally, out of 157 computationally identified siRNAs, only 4 effective siRNA molecules were selected for each target gene which is proposed to exert the best action based on GC content, free energy of folding, free energy of binding, melting temperature, heat capacity and molecular docking analysis with Human AGO2 protein. Our engineered siRNA candidates could be used as a genome-level therapeutic treatment against various sequenced SARS-CoV-2 strains in India. However, future COVID-19, A viral disease, caused by the new strain of severe acute respiratory syndrome coronavirus (SARS-CoV), known as SARS-CoV-2, has challenged humanity at the beginning of 2020, impacting the lives of billions of people worldwide. Since its outbreak in late December 2019 in Wuhan, China, after a sudden epidemic of atypical pneumonia with unclear illness aetiology, it has caused substantial morbidity and mortality all over the world [1] . Fever, cough, fatigue, dyspnea, and headache [2] are the most common signs of this disease, but it can also be asymptomatic [3] . According to the World Health Organization (WHO), India has recorded 4.77 lakh Coronavirus deaths since the outbreak began in March 2020. Furthermore, 3.47 crore Coronavirus cases were reported in India [4] . Especially the second wave of COVID-19 has resulted in a rise in cases, a decline in crucial treatment supplies, and an increase in deaths, particularly among the young [5] . In severe cases, the patient may die due to massive alveolar damage and progressive respiratory failure [2] . Also, several occurrences of mucormycosis, popularly known as the black fungus, have been reported in patients with diabetes and patients with COVID-19, as well as patients who were recovering from infection [6] . The SARS-CoV-2 is enveloped by single-stranded positive-sense RNA and has 50% and 80% homology with Middle East Respiratory Syndrome virus and SARS-CoV, respectively and it comprises four structural proteins: envelope (E), membrane (M), nucleocapsid (N) and spike (S) [7, 8] . The E protein is the smallest and most cryptic of the major structural proteins. It is extensively expressed inside the infected cell during the replication cycle, but only a small proportion of it is incorporated into the virion envelope [9] . The M glycoprotein is the most common structural protein in coronaviruses [10] ; All other structural proteins can bind to the M protein. Binding with M protein aids in the stabilisation of N proteins and enhances viral assembly completion by stabilising the N protein-RNA complex within the internal virion [11] . In the coronavirus life cycle, the N protein plays two important roles. Its principal function is to form the viral RNA-protein (vRNP) complex with genomic RNA and to mediate vRNP packaging into virions. Second, the N protein is hypothesised to recruit host factors and promote RNA template switching in RTCs during the early stages of infection, facilitating viral RNA synthesis and translation [12] . The spike glycoprotein extends from the viral surface and is made up of two functional subunits (S1 and S2). The S1 subunit aids in the interaction of SARS-CoV-2 with the host cell's angiotensin-converting enzyme 2 (ACE2) receptor, whereas the S2 subunit aids in the fusion of viral and host cell membranes [13] . Since these structural proteins (E, M, N and S) are crucial for viral pathogenesis, these could be the main targets for designing therapeutics. RNA interference or RNAi, the biological mechanism by which double-stranded RNA (dsRNA) induces gene silencing by targeting complementary mRNA for degradation ( Figure 1 ) [14] , is a huge breakthrough in disease therapy and is changing the way scientists analyse transcriptional activity [15] and it is a prospective tool for the control of human viral infections. Small interfering RNA (siRNA) is a double-stranded RNA of around 21 base pairs long RNA duplex bearing overhangs of two nucleotides on the 3′ end [16] . The strand of siRNA which has complementarity with the target gene is the guide strand and the other strand is the passenger strand. Combinations of chemically synthesised siRNA duplexes targeting SARS-CoV genomic RNA result in therapeutic activity of up to 80% inhibition, according to studies [17] . Due to the complementarity of seven nucleotides in the seed region of siRNA with the off-target gene, there is a risk of off-target gene silence or unexpected gene downregulation. The off-target binding impact is mediated by the melting temperature or thermodynamic stability of the duplex of seed siRNA sequence (2-8 nucleotides of siRNA guide strand from 5′ end) and the target gene, according to studies. Thus, a siRNA with perfect complementarity solely for the target gene and low seed-target duplex thermostability (Tm less than 21.5°C) can efficiently remove siRNA off-target binding. Furthermore, selecting a siRNA with at least two mismatches with any other off-target region can lower the likelihood of siRNA binding to the unwanted off-target sequence [18, 19] . Along with the seed region, the non-seed region of the guide strand has been shown to be effective in mediating off-target impact, yet there was a negative correlation between Tm value and GC content and off-target effect [20] . Because structural proteins expressed by the E, M, N and S genes are involved in virus survival and infectivity, they can be used as a target for suppressing SARS-CoV-2 infection. In silico approaches to discovering innovative therapeutic approaches are becoming increasingly popular, and viruses are no exception. Therefore, in this present study, we have designed siRNAs specific to the various conserved region of the envelope protein, membrane glycoprotein, nucleocapsid and surface glycoprotein genes by analyzing 811 Indian SARS-CoV-2 strains and finally shortlisted 4 effective siRNA molecules for each gene which will inhibit the translations of these target genes and allow the host to discard this infection ( Figure 2 ). Table 1 ) [21] . Out of 812 gene sequences, 811 were Indian strains and the remaining one was the Wuhan-Hu-1 strain (China, Genbank accession: NC_045512.2) that was used as a reference genome sequence for multiple sequence alignment and site numbering for amino acids. Multiple sequence alignment of all 812 gene sequences of E, M, N and S genes were performed using the ClustalW [22] algorithm in the MEGA-X [23] program to find the conserved regions. siDirect version 2.0 [19] is used to design effective and target-specific siRNA molecules against SARS-CoV-2 E, M, N and S gene sequences. This tool utilizes some rules such as Ui-Tei [24] , Amarzguioui [25] and Reynolds [26] algorithms (Table 1) for designing siRNAs and the melting temperature (Tm) of the seed-target duplex was kept below 21.5°C as a default parameter. To avoid off-target silencing, it chooses siRNA sequences with at least two mismatches to any other non-targeted transcripts. The BLAST [27] tool was used to identify the possible off-target matches in a human genomic transcript against the whole GenBank [28] database by applying the expected threshold value 10 and BLOSUM 62 matrix [29] as a parameter. OligoCalc [30] was used to analyze the GC content of predicted siRNA molecules. The secondary structure and the free energy (ΔG) of folding for predicted siRNA molecules were computed using the MaxExpect [31] program of the RNA structure web server [32] . Higher energy values imply that those molecules are less likely to fold and it indicates that they are better candidates. The thermodynamic interaction between the target strand and the siRNA guide strand was predicted using the DuplexFold [33] program of the RNA structure web server [32] . It folds two sequences of RNA into their lowest hybrid free energy conformation. Higher interaction between the target and the guide strand will aid in a better predictor of siRNA effectiveness. The heat capacity plot and concentration plot were calculated for the predicted siRNAs using the DINA Melt web server [34] . The ensemble heat capacity (Cp) is plotted as a function of temperature, with the melting temperature Tm (Cp) indicated in Table 2 . The detailed heat capacity plot also shows the contributions of each species to the ensemble heat capacity, whereas the concentration plot-Tm (Conc), the point at which the concentration of doublestranded molecules of one-half of its maximum value defines the melting temperature Tm (Conc). The proper interaction between siRNA duplex and RISC complex, as well as guide siRNA strand and target mRNA inside RISC complex, are required to trigger a sufficient antiviral response via RNAi-mediated viral gene silencing [14] . Hence for interaction pattern analysis, molecular docking was carried over between siRNA (guide strand) and argonaute protein. The 3D modelled structure of guide siRNA was generated using a 3dRNA v2.0 web server [35] . 3dRNA is an automatic, fast, and high-accuracy RNA tertiary structure prediction method. It uses sequence and secondary structure information to build three-dimensional structures of RNA from template segments. We used an optimization procedure which is an integrated 3dRNA pipeline, to model our predicted siRNAs. The human argonaute 2 (Ago2) protein is involved in directing siRNA to the target RNA cleavage site for mRNA degradation, translational suppression, or a combination of the both. Hence, generating an adequate antiviral response using RNAi-mediated viral gene silencing should involve the proper interaction between the siRNA duplex (guide strand) and the Ago2 protein that is present within the RISC complex [36] . Therefore, the three-dimensional structure of human argonaute 2 (Ago2) protein ( Figure 3 ) was retrieved from RCSB Protein Data Bank (PDB ID: 4F3T) for computer-assisted molecular docking studies [37] . BIOVIA discovery studio [38] was used to prepare the proteins through correction of bonds, removal of J o u r n a l P r e -p r o o f unrelated chemical complexes, and elimination of water molecules and hetatom groups. The missing side-chain atoms in protein have been reconstructed and then the structure was minimized with the partial implementation of the GROMOS96 force-Field using Swiss-PdbViewer [39] . The molecular docking between siRNA (guide strand) and argonaute protein was performed using the HDOCK server [40] . This server uses a hierarchical Fast Fourier Transform (FFT)related docking program. It is a hybrid docking algorithm that uses free docking and template-based modelling and shows high efficacy and robustness. Further docking analysis of protein-siRNA complex was performed using BIOVIA discovery studio [38] . 3.1. Target specific designing of potential siRNA molecules siDirect 2.0 web server predicted 7 siRNA for the envelope gene, 10 siRNA for the membrane gene, 3 siRNA for the nucleocapsid phosphoprotein gene and 137 siRNA for the surface glycoprotein gene (Supplementary Table 6 , 7, 8 and 9) that followed algorithms of Ui-Tei, Amarzguioui and Reynolds ( Table 1 ). The seed target duplex stability (Tm) values for all predicted siRNAs were less than 21.5 °C, which shows that predicted siRNAs may avoid nontarget binding. The consensus targets obtained by siDirect 2.0 were subjected to NCBI-BLAST for searching similarity against the whole human genome, however, no significant matches were found. This shows that our predicted siRNA may not interact in any place other than the viral target site. The sequence's GC content is a vital parameter that influences siRNA functionality. (Supplementary Tables 10, 11, 12 and 13) . GC content analysis of the predicted siRNAs was ranged from 21% to 38% for the envelope gene, 10% to 38% for the membrane gene, 36% to 40% for the nucleocapsid phosphoprotein gene and 10% to 43% for the surface glycoprotein gene. The calculated free energy of folding ranged from 1.3 to 1.7 kcal/mol for the envelope gene, 1.6 to 2 kcal/mol for the membrane gene, 1.7 to 1.9 kcal/mol for the nucleocapsid phosphoprotein gene and 1.1 to 2.9 kcal/mol for the surface glycoprotein gene (Supplementary Tables 10, 11, 12 and 13 ). The binding free energy between the target and guide strand was calculated. The values spanned from -25.8 to -34.1 kcal/mol for the envelope gene, -20.9 to -33.5 kcal/mol for the membrane gene, -32.5 to -35.8 kcal/mol for the nucleocapsid phosphoprotein gene and -20.5 to -36.1 kcal/mol for the surface glycoprotein gene (Supplementary Tables 10, 11, 12 and 13) . The predicted siRNAs' Tm (Cp) and Tm (Conc) were computed. The greater the values of these two melting temperatures, the more effective the siRNA species. Best siRNAs have been filtered for each gene based on the criteria as follows GC content (33% to 43%) and free energy of binding (>=-31.3 kcal/mol) ( Table 2 ) and subjected to docking studies. The docking score, interaction statistics and interactive residues were provided in Table 3 COVID-19, respiratory illness caused by the SARS-CoV-2 coronavirus, has recently become pandemic. The pandemic is sweeping through India at a pace that has staggered scientists. To combat this pandemic, promising varied techniques such as gene therapy, as well as other therapeutic tactics such as the creation of vaccines, drugs, monoclonal antibodies, peptides, and other therapeutic strategies are being suggested. The present study introduced an in silico approach to design potential siRNA molecules of Indian SARS-CoV-2 strains targeting virus structural genes which included envelope protein (E), membrane protein (M), nucleocapsid phosphoprotein (N) and surface glycoprotein (S). Among the 811 SARS-CoV-2 strains from throughout India, 243 conserved regions (10 envelope protein, 26 membrane protein, 67 nucleocapsid phosphoprotein and 140 surface glycoprotein) (Supplementary Tables 2, 3, 4 and 5) were found. Conserved regions that are less than 21 nucleotides were not taken out of thefor further analysis. The siDirect webserver was used to identify potential targets and generate corresponding siRNAs using conserved sequences. The task is completed by siDirect in three steps: highly functional siRNA selection, seed-dependent off-target effects reduction, and near-perfect matched genes removal. There were 157 target regions, 7 envelope proteinsprotein, 10 membrane proteinsprotein, 3 nucleocapsid phosphoproteinsphosphoprotein, and 137 surface glycoproteinsglycoprotein (Supplementary Tables 6, 7, 8 and 9 ). To improve the results, the U, R, A (Ui-Tei, Amarzguioui, and Reynolds) guidelines (Table 1 ) were used to identify the siRNAs. The formation of non-target siRNA bonds was eliminated by lowering the melting temperature (Tm) below 21.5 degrees Celsius. Despite the fact that the primary goal of siRNA design is to silence specific targets, there is still the possibility of silencing an unknown number of unwanted genes [41] . The off-target activity caused by siRNA can be explained by two mechanisms: mRNA depletion or translational suppression at the protein level [42] . First, siRNA can tolerate several mismatches to the mRNA target while also retaining good silencing capacity with imperfect complementarity [41] . Another possible mechanism for off-target activity is that siRNAs' physical structure, ~21 nucleotides (nt) RNA oligomers, appears to be identical to the related class of microRNAs. MicroRNAs are short, endogenously transcribed RNAs that prevent mRNA from being translated rather than being degraded. MicroRNAs appear to contain RNA oligo and RNA target mismatches built into their structure. Together the mechanisms of siRNA mismatch tolerance and microRNA translation inhibition create a risk of off-target activity when ~21 nt RNAs are introduced into human cells [41, 42] . As a result, off-target effects may occur even when siRNAs are generated using different algorithms. To ensure siRNA target specificity, siRNA design processes are often followed by a BLAST search for cross-reactive 21-bp siRNA sequences. BLAST results show that the SARS-CoV-2 target sequences (E, M, N, and S) had no significant hits against the entire genome, ruling out the possibility of off-target silencing. The GC-content of siRNA influences its functioning, and there is an inverse link between GC-content and siRNA function. siRNA targets with a very high GC content may be more prone to folding, which could weaken target accessibility. For a siRNA to be effective, it needs to have a medium proportion of GC content, which ranges from 31.6 to 57.9% [43] . Therefore, molecules with a GC content of less than ~32% were excluded from the final selection. For each of the 157 predicted species, the GC content ranged from 10% to 43%. The GC content of the filtered siRNAs is greater than or equivalent to 33%. ( Table 2 ). The RNA structure webserver was used to estimate probable folding structures and corresponding minimal free energy using guide strands of predicted siRNAs. The probable folding and minimum free energy of folding of the predicted siRNA guide strands were calculated with the assumption that siRNA with the positive free energy of folding may be more permissible to the target, resulting in better gene silencing [44] . At 37°C, the free energy of folding of the final siRNAs is more than zero ( Figure 4 , Table 2 ), implying that the predicted siRNAs are more accessible for efficient binding. DuplexFold was used to calculate the binding energy of the target mRNA and guide siRNA interactions. Lower binding energy suggests a better interaction, and hence a better likelihood of inhibiting the target. All 157 predicted siRNAs had binding free energy values ranging from -36.1 to -20.5 kcal/mol. (Supplementary Tables 10, 11, 12 and 13) . We chose the best siRNAs for each E, M, S, and N gene based on GC content and free energy of binding criteria for further docking analysis. The selected siRNAs have a free energy of binding equal to or less than -31.3 kcal/mol ( Figure 5 , Table 2 ), implying that predicted siRNAs are more interactive with their targets. The Ago2 and modelled siRNAs were docked with the help of the HDOCK server. To perform traditional global docking, the server employs a hybrid-docking technique that combines template-based modelling and free docking, as well as follows a hierarchical FFTbased docking algorithm. Resultant Ago2-siRNAs docked complexes were retrieved from the server and manually analysed to determine the optimal docked complex based on docking J o u r n a l P r e -p r o o f score, interaction pattern analysis, and placement of siRNA within the Ago2 domains (PAZ, MID, N-terminal and PIWI domain) (Figure 3 ). Based on the interaction pattern analysis of 12 filtered siRNAs guide strand with target mRNA (Table 3) , 4 best final siRNAs (e5, m2, n1 and g80) were sorted and hypothesized for each gene that could mediate the post-transcriptional gene silencing of the targeted gene ( Table 4 ). The docking score of final siRNAs (e5, m2, n1 and g80) showed impressive binding with argonaute protein within the docked complex ( Figure 6 ), ranging from -322.34 to -358.89 kcal/mol. The detailed elucidation about interacting residues is provided in Supplementary Table 14 . Therapeutic applications of siRNAs are quite challenging as they have several issues, such as siRNA instability, limited cellular uptake, and the lack of a competent delivery method [45] . However, an appropriate promoter-controlled vector can aid in the targeting of therapeutic genes to the targeted cell for efficient gene therapy [46] . Vector-based siRNA in plasmid form can also be used to target desired genes within a given cell culture to examine the potentiality of a newly created siRNA [47] . Furthermore, plasmids containing siRNA can be directly delivered into the organs of choice [48] . In this present study, four potential siRNA candidates were hypothesized to be effective at binding and cleaving specific SARS-CoV-2 (Indian strains) mRNA targets of structural proteins ( Table 4 ). The suggested therapeutic molecule may be used to treat the COVID-19 pandemic in India on a broad scale because since the study comprises a large array of 811 SARS-CoV-2 sequences from all over India, however, it needs further validation in in vitro and in vivo studies. siRNA therapy might be a promising tool of the RNAi pathway for controlling viral infections in humans by PTGS of the significant gene in various biological systems. In this study, four siRNA molecules were predicted to be effective against the envelope gene (E), membrane gene (M), nucleocapsid phosphoprotein gene (N) and surface glycoprotein gene (S) of 811 Indian strains of the SARS-CoV-2 virus using various computational tools considering all maximum parameters in prominent conditions molecular modelling and docking analysis. However, future applications will necessitate additional validations in vitro and in vivo animal models. In the battle against the Covid-19 pandemic, these synthetic molecules may be used as novel antiviral therapy and provide a basis for the researchers and pharmaceutical industry to create antiviral therapeutics at the genome level. All necessary data generated or analyzed during this study are included in this article. Any additional data could be available from the corresponding author upon request. Amarzguioui rules Reynolds ruless  A/U at the 5′ end of the antisense strand  G/C at the 5′ end of the  The A/U differential of the duplex end should be >0  The designed siRNA must maintain a GC content between 30% to 52% (1 point) sense strand  At least five A/U residues in the 5′ terminal one-third of the antisense strand  The absence of any GC stretches of more than 9 nt in length  Strong binding of 5 sense strand  Position 1 must have any bases other than U  Position 6 must have A constantly  Weak attachment of 3′ sense/passenger strand. J o u r n a l P r e -p r o o f Virtual Screening Reveals Potential Anti-Parasitic Drugs Inhibiting the Receptor Binding Domain of SARS-CoV-2 Spike protein Clinical features of patients infected with 2019 novel coronavirus in Emerging 2019 novel coronavirus (2019-NCoV) pneumonia COVID-19) -the data -Statistics and Research -Our World in Data Implications of the second wave of COVID-19 in India In silico prediction and structurebased multitargeted molecular docking analysis of selected bioactive compounds against mucormycosis An in silico analysis of effective siRNAs against COVID-19 by targeting the leader sequence of SARS-CoV-2 Discovering small-molecule therapeutics against SARS-CoV-2 Coronavirus envelope (E) protein remains at the site of assembly Genotype and phenotype of COVID-19: Their roles in pathogenesis Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2): An overview of viral structure and host response The SARS-CoV-2 nucleocapsid phosphoprotein forms mutually exclusive condensates with RNA and the membrane-associated M protein COVID-19): A literature review Designing an effective therapeutic siRNA to silence RdRp gene of SARS-CoV-2 Molecular Mechanisms and Biological Functions of siRNA Small RNAs: regulators and guardians of the genome Prophylactic and therapeutic effects of small interfering RNA targeting SARS-coronavirus Optimal choice of functional and off-target effect-reduced siRNAs for RNAi therapeutics Ui-Tei, siDirect 2.0: updated software for designing functional siRNA with reduced seed-dependent off-target effect, BMC Bioinforma The siRNA Non-seed Region and Its Target Sequences Are Auxiliary Determinants of Off-Target Effects CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, positionspecific gap penalties and weight matrix choice Molecular Evolutionary Genetics Analysis across Computing Platforms Guidelines for the selection of highly effective siRNA sequences for mammalian and chick RNA interference An algorithm for selection of functional siRNA sequences Rational siRNA design for RNA interference Basic local alignment search tool Amino acid substitution matrices from protein blocks OligoCalc: an online oligonucleotide properties calculator Improved RNA secondary structure prediction by maximizing expected pair accuracy RNAstructure: web servers for RNA secondary structure prediction and analysis DINAMelt web server for nucleic acid melting prediction 0: An Updated Web Server for RNA 3D Structure Prediction Human Argonaute2 mediates RNA cleavage targeted by miRNAs and siRNAs The Protein Data Bank PdbViewer: an environment for comparative protein modeling The HDOCK server for integrated protein-protein docking Preclinical and clinical development of siRNA-based therapeutics Many commonly used siRNAs risk off-target activity A structural interpretation of the effect of GC-content on efficiency of RNA interference Design of potential siRNA molecules for hepatitis delta virus gene silencing Disulfide crosslinked stearoyl carrier peptides containing arginine and histidine enhance siRNA uptake and gene silencing In Silico Design and Experimental Validation of siRNAs Targeting Conserved Regions of Multiple Hepatitis C Virus Genotypes Small interfering RNA Inhibits Hepatitis B virus replication in mice The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.☐The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:Graphical abstract Highlights  We designed and analyzed siRNA molecules that are predicted to be highly effective considering the wide array of 811 Indian SARS-CoV-2 strains and mediate PTGS of specified mRNA.  4 effective siRNAs were sorted as potential therapeutic agents against E, M, N and S genes of SARS-CoV-2 based on their GC content, free energy of folding, free energy of binding, melting temperature, heat capacity and molecular docking analysis.  The goal of this computational study is to reduce the amount of time and effort required in traditional wet lab methods based on trial and error, and our developed siRNAs have the potential to be an effective RNAi therapeutic.