key: cord-0905576-f6o1aynz
authors: Samad, Abdus; Ahammad, Foysal; Nain, Zulkar; Alam, Rahat; Imon, Raihan Rahman; Hasan, Mahadi; Rahman, Md. Shahedur
title: Designing a multi-epitope vaccine against SARS-CoV-2: an immunoinformatics approach
date: 2020-07-17
journal: Journal of biomolecular structure & dynamics
DOI: 10.1080/07391102.2020.1792347
sha: 36409011e8f566194cbdc5fd00cfe90bdcd8eab4
doc_id: 905576
cord_uid: f6o1aynz

Ongoing COVID-19 outbreak has raised a drastic challenge to global public health security. Most of the patients with COVID-19 suffer from mild flu-like illnesses such as cold and fever; however, few percentages of the patients progress from severe illness to death, mostly in an immunocompromised individual. The causative agent of COVID-19 is an RNA virus known as severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2). Despite these debilitating conditions, no medication to stop the disease progression or vaccination is available till now. Therefore, we aimed to formulate a multi-epitope vaccine against SARS-CoV-2 by utilizing an immunoinformatics approach. For this purpose, we used the SARS-CoV-2 spike glycoprotein to determine the immunodominant T- and B-cell epitopes. After rigorous assessment, we designed a vaccine construct using four potential epitopes from each of the three epitope classes such as cytotoxic T-lymphocytes, helper T-lymphocyte, and linear B-lymphocyte epitopes. The designed vaccine was antigenic, immunogenic, and non-allergenic with suitable physicochemical properties and has higher solubility. More importantly, the predicted vaccine structure was similar to the native protein. Further investigations indicated a strong and stable binding interaction between the vaccine and the toll-like receptor (TLR4). Strong binding stability and structural compactness were also evident in molecular dynamics simulation. Furthermore, the computer-generated immune simulation showed that the vaccine could trigger real-life-like immune responses upon administration into humans. Finally, codon optimization based on Escherichia coli K12 resulted in optimal GC content and higher CAI value followed by incorporating it into the cloning vector pET28+(a). Overall, these results suggest that the designed peptide vaccine can serve as an excellent prophylactic candidate against SARS-CoV-2. Communicated by Ramaswamy H. Sarma

Coronavirus disease 2019 is an acute highly infectious disease. Patients with COVID-19 mostly feel like flu-like symptoms including cold and fever; a few percentages of the patients suffer from the respiratory tract infection leading to severe atypical pneumonia that eventually ends up in a case of fatality (Boopathi et al., 2020; Joshi et al., 2020; Vankadari & Wilce, 2020) . Moreover, patients admitted to the intensive care unit were likely to report cardiovascular, respiratory disease, cerebrovascular, abdominal pain, endocrine, anorexia and digestive diseases (Chan et al., 2020; . Acute cardiac injury and acute respiratory distress syndrome are commonly observed in severe cases and is strongly associated with the mortality rate (Abdelli et al., 2020; Wu et al., 2020) . The infected patients can transmit the virus through coughs, sneezes, exhales and many other ways, hence, playing an essential role in human to human transmission (Chan et al., 2020) . Infection with the virus is sometimes asymptotic, which also plays a vital role in the transmission process (Shen et al., 2020) .

The COVID-19 outbreak has already been taken place all over the world with a total of 7,732,952 confirmed cases and 428,248 death cases (13 June 2020, 03:33 GMT) over 213 countries and territories around the world (www.worldometers.info). As of 13 June 2020, the outbreak in Bangladesh includes 81,523 confirmed and 1095 death cases, and the number of cases is increasing drastically day by day. The causative agent of the COVID-19 outbreak is severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) which is a positive-sense single-stranded RNA virus that belongs to the family of Coronaviridae (Boopathi et al., 2020; Umesh et al., 2020) , and genus beta-coronavirus Li & De Clercq, 2020) . This virus was first identified in patients with a cluster of pneumonia in the province Wuhan of China on 29 December 2019 (Das et al., 2020; Kumar et al., 2020; Zhu et al., 2020) . The World Health Organization (WHO) country office in China announced the official declaration of this virus on 31 December 2019 (Calisher et al., 2020; Heymann & Shindo, 2020) . The new coronavirus was named 'SARS-CoV-2' by the International Committee on Taxonomy of Viruses (ICTV), and the disease caused by the pathogen was announced as COVID-19 by the WHO on 11 February 2020 (Gupta et al., 2020; Muralidharan et al., 2020; Wu et al., 2020) . Epidemiological investigations of the Wuhan zoonotic virus revealed 89.1% nucleotide similarity between SARS-CoV-2 and previously originated group of SARS-like coronavirus Vankadari & Wilce, 2020; Wu et al., 2020) .

During the SARS-CoV-2 infection, the counts of CD4þ and CD8þ T-cells are increased in the peripheral blood, and cytotoxic granules are upturned with high concentration . The over-activation of T-cells causes injury to the immune system of the infected patients resulting in the characteristic feature of 'Lymphopenia' increasing the disease severity. In contrast, less effective T-cell responses may allow the progression of viral pathology and thus increased mortality in SARS-CoV-2 infected patients (Ahammad et al., 2019; . The CD4þ and CD8þ responses provide a long-lasting protection against COVID-19 (Enayatkhani et al., 2020) . Moreover, antibody-mediated immune response along with cellular immunity plays a critical role to induce protectivity against these infections (Enayatkhani et al., 2020) . In recent studies, it has been illustrated that the nucleotide structure of this virus particle has a similarity with SARS-like coronavirus. Their genome was encoded with 16 different non-structural proteins and four main structural proteins, including spike (S), envelope (E), nucleocapsid (N), and membrane (M) proteins (Hasan et al., 2020; Sarma et al., 2020; Wahedi et al., 2020; Wu et al., 2020) . The S-protein formed the viral outer layer; the N-protein helps in the viral replication, genome construction and host cellular response; the Mprotein determined the envelope shape and the E-protein functions in production and maturation of the SARS-CoV-2 (Astuti & Ysrafil, 2020) . The S-protein of the virus contains two subunit, S1 and S2; the S1 subunit of the S-protein recognizes the host T-cells while S2 subunit mediates fusion between the viral and host T-cells (Astuti & Ysrafil, 2020) and characterizes as a highly antigenic and surface exposure Wrapp et al., 2020) . The CD8þ and CD4þ T-cells recognize viral epitopes presented by the major histocompatibility complex class I (MHC I) and class II (MHC II), respectively (Abdellrazeq et al., 2020; Borthwick et al., 2020) . The heterogeneity in T-cells responses to SARS-CoV-2 may, in part, be related to the capacity to recognize the viral antigens in the context of MHC I and MHC II proteins (Astuti & Ysrafil, 2020) . It has been found that T-cell epitopes of SARS-CoV-2 spike protein elicit a T-cells immune response in patients who recover from the disease, and most of these immunogenic epitopes were localized to the S protein of the virus (Astuti & Ysrafil, 2020) . The S-protein has strong interactions and binding affinity to the human angiotensin-converting enzyme 2 (ACE2) receptor and facilitates viral entry into the target cell (Sinha et al., 2020; . The S-protein of the SARS-CoV-2 is the major host interacting protein, which causes cell adhesion and virulence to the human host (Vankadari & Wilce, 2020; Wu et al., 2020) . The virus S-protein entry is mediated by ACE2, and results in an inflammatory cascade initiation by the innate immune system of the host (Astuti & Ysrafil, 2020) . So, targeting S-protein can provide an immunogenic response in the human host, and has been chosen for designing a multi-epitopes vaccine candidate against the SARS-CoV-2. The structural pattern of SARS-CoV-2 protein can be recognized by the transmembrane toll-like receptor 4 (TLR4) which induces inflammatory cytokines or chemokines reaction. The TLR4 protein plays a vital role in the host pathogenesis. Moreover, the involvement of the receptors has also been reported in various immune protective responses by the host. Immune responses are a crucial step to the pathophysiology of the SARS-CoV-2 virus-related disease, and initiation of immune response targeting TLR4 can trigger the anti-viral host defense mechanisms necessary for the elimination of the COVID-19 related infection (Astuti & Ysrafil, 2020) .

Vaccine is an immune-modulatory preparation that triggers a specific immune response against a foreign particle within the host body. A vaccine is now the primary demand to save millions of people from the COVID-19 pandemic. The current world situation is releasing the necessity of an implausible and effectiveness of different anti-viral drugs or vaccine candidates against the SARS-CoV-2. However, no effective drug or vaccine candidates have been developed that can fight against the SARS-CoV-2 (Elfiky, 2020) . Therefore, a multi-epitope vaccine consisting of potential Tand B-cell epitopes can be an ideal approach for the prevention of COVID-19 (Astuti & Ysrafil, 2020) . The vaccine can produce both cellular and humoral immune responses against specific pathogens without producing any immune complications. Besides, it is very easy to control, cause the effectiveness of the vaccine to be regulated by choosing the specific and desired allelic interactions, which provide robust and diverse immune response over a large group of people (Elfiky, 2020) . In multi-epitope vaccine, the biohazard risk is lower as compared to other types of immunizations. In this research, a multi-epitope vaccine has been constructed using the immunoinformatics approach. Epitopes used for the vaccine construction were non-toxic, non-allergenic, highly immunogenic and antigenic. A sufficient number of linkers were used to combine those selected epitopes resulting in busting the immunogenic activity of SARS-CoV-2 vaccine (Gaafar et al., 2019; Li et al., 2014) . A flow chart representing the overall procedure from the antigen selection to vaccine construction and evaluation is illustrated in Figure 1 .

For antigen selection, we collected available SARS-CoV-2 proteomes from the ViPR (https://www.viprbrc.org/) database (Pickett et al., 2012) . The outer membrane of the SARS-CoV-2 is formed by the spike glycoproteins. With the help of these glycoproteins, they adhere to the human host and enter into the host immune system . Due to the direct involvement of glycoproteins in pathogenesis, we considered the spike glycoprotein of the SARS-CoV-2 for multiepitope vaccine design. Initially, we isolated all the spike glycoprotein, and the selected protein sequences of the virus were downloaded in FASTA format. The protective antigens of the surface glycoprotein were checked by VaxiJen v2.0 (http://www.ddg-pharmfac.net/vaxijen/) server (Doytchinova & Flower, 2007) and ANTIGENpro (http://scratch.proteomics. ics.uci.edu/) server with a threshold value 0.4 was set for both of them (Magnan et al., 2010) . Finally, we selected the spike glycoprotein with the highest antigenic score for further investigations.

Cytotoxic T-lymphocytes (CTLs) represent one of several types of cells of the immune system that have the capacity to kill other infectious cells directly . They go right away inside the virus-cell and play an important role in the host defense mechanism. For the prediction of CTLs epitope, the sequence of the selected protein was submitted into the NetCTL v1.2 server available at http://www.cbs.dtu. dk/services/NetCTL/ (Larsen et al., 2007) . The predicted epitopes were further assessed through the VaxiJen v2.0 (Doytchinova & Flower, 2007) , MHC class I immunogenicity (http://tools.iedb.org/immunogenicity/) (Calis et al., 2013) , ToxinPred (http://crdd.osdd.net/raghava/toxinpred/) , and AllerTop v2.0 (https://ddg-pharmfac.net/ AllerTOP/) (Dimitrov et al., 2013) servers. The default parameters of those servers were used for all the predictions.

Helper T-cells (HTLs) are an integral part of adaptive immunity that recognizes foreign antigens and activates B and cytotoxic T-cells resulting in destruction of the infectious pathogen . To determine the HTL epitopes, we used the IEDB's MHC class II binding allele prediction tool, available at http://tools.iedb.org/mhcii/. The HTL epitopes were selected based on a percentile rank of 5% using the CONSENSUS method (Wang et al., 2010) . The predicted epitopes were further evaluated based on their antigenicity and cytokine, i.e. interferon-c (IFNc), interleukin-4 (IL4) and interleukin-10 (IL10) inducing abilities. The antigenicity was anticipated with the VaxiJen v2.0 server while IFNc, IL4, and IL10 features were predicted using IFNepitope (http://crdd. osdd.net/raghava/ifnepitope/) , IL4pred (http://crdd.osdd.net/raghava/il4pred/) and IL10pred (http://crdd.osdd.net/raghava/IL-10pred/) (Nagpal et al., 2017) servers, respectively, with default parameters.

B-cell epitopes are essential to induce humoral or antibodymediated immunity. B-cells consist of groups of amino acids that interact with the secreted antibodies and activate the immune system to destroy the pathogens (Nain et al., 2019) . Therefore, we predicted the linear B-lymphocyte (LBL) epitopes using the iBCE-EL server, available at http://www.thegleelab.org/iBCE-EL/ with default parameters (Manavalan et al., 2018) . The predicted LBL epitopes were also evaluated using the VaxiJen v2.0, ToxinPred, and AllerTop v2.0 servers.

In computational vaccine design, the population coverage directly indicates the worldwide effectiveness of the vaccine by evaluating the prevalence of HLA (Human Leukocyte Antigen) alleles related to the epitope of interest. Therefore, the population coverage was calculated using the T-cell epitopes with their respective HLA binding alleles. To accomplish this, selected epitopes along their allelic information was submitted to the IEDB population coverage tool (Bui et al., 2006) .

For modeling, the selected epitopes of CTL and HTL were submitted into PEP-FOLD v3.0 (https://bioserv.rpbs.univ-parisdiderot.fr/services/PEP-FOLD3/) server. The sOPEP sorting scheme with 200 simulations was selected for the operation (Latysheva & Babu, 2016) . By analyzing the epitope-wise HLA binding alleles, allele HLA-B Ã 15:01 and HLA-C Ã 06:02 were considered for selected CTL epitopes, while DRB1 Ã 01:01 and DRB1 Ã 15:01 were selected for HTL epitopes. The crystal structures of the HLA alleles were retrieved from the Protein Data Bank (PDB) (https://www.rcsb.org/) (Berman et al., 2000) followed by processing with BIOVIA Discovery Studio 2017. For molecular docking, a grid-box around the active site of each HLA allele was defined by the AutoDock tool. Finally, molecular docking was performed between the epitopes and respective HLA alleles using the AutoDock Vina script (Trott & Olson, 2010) . The respective co-crystal ligands were used as the positive control to compare the epitope binding efficiency. The docked complex was visualized in BIOVIA Discovery Studio 2017.

The vaccine construct was designed by using the selected CTL, HTL, and LBL epitopes as well as a suitable adjuvant that was linked by the appropriate linkers (Dorosti et al., 2019; Nain et al., 2019) . Here, we used TLR4 agonist as the adjuvant since TLR4 was recognized by viral glycoproteins, and the adjuvant is required for optimal translation and maximal rate of synthesis of the target vaccine candidate (Olejnik et al., 2018; Pandey et al., 2018) . Therefore, 50S ribosomal protein L7/L12 (NCBI ID: P9WHE3) was considered as the adjuvant to improve the immunogenicity of the vaccine candidate. The adjuvant was linked to the vaccine front with a bi-functional linker EAAAK that has the ability of several lengths of helixforming peptides to separate two weakly interacting b domains. In contrast, the selected CTL was linked with the help of Ala-Ala-Tyr (AAY) linkers, the HTL was linked with Gly-Pro-Gly-Pro-Gly (GPGPG) linkers and the LBL was linked with Lys-Lys (KK) linker (Dorosti et al., 2019; Nain et al., 2019) . The AAY linker is a type of cleavage site of proteasomes that was used to influence protein stability, reduce less immunogenicity and enhance epitope presentation (Abdellrazeq et al., 2020; Borthwick et al., 2020) . The GPGPG, known as the glycine-proline linker, prevents the formation of 'junctional epitopes' and facilitates the immune processing, where the bi-lysine KK linker helps to preserve their independent immunogenic activities of the vaccine construct.

The physiochemistry indicates the basic properties of a protein. The physicochemical features of the vaccine were anticipated using the ProtParam server available at https://web. expasy.org/protparam/ to understand the fundamental nature of the vaccine (Gasteiger et al., 2005) . We also evaluated the immunological properties through VaxiJen v2.0 (Doytchinova & Flower, 2007) , MHC-I immunogenicity (Calis et al., 2013) , AllerTop (Dimitrov et al., 2013) , and SOLpro (Magnan et al., 2009 ) servers.

The two-dimensional (2D) structural features such as alphahelix, beta-turn, and random coils of the construct were identified by SOPMA (Self-Optimized Prediction Method with Alignment) server at https://npsa-prabi.ibcp.fr/NPSA/npsa_ seccons.html (Geourjon & Del eage, 1995) and PSIPRED v4.0 (PSI-blast based secondary structure prediction) server at http://bioinf.cs.ucl.ac.uk/psipred/ (Buchan et al., 2013) with default parameters. SOPMA has more than 80% prediction accuracy (Geourjon & Del eage, 1995) . The 2D structural features were retrieved and evaluated to understand the composition quality of the vaccine.

2.10. Homology modeling, 3D structure refinement and validation

The constructed vaccine was submitted into I-TASSER (Iterative Threading Assembly Refinement) online web portal (https://zhanglab.ccmb.med.umich.edu/I-TASSER/) for threedimensional (3D) structure prediction (Roy et al., 2010) . The I-TASSER web produces the structure of the protein and its functions most accurately using a state-of-the-art algorithm in the form of a 3D structure (Roy et al., 2010) . This web server can predict and determine the C-score, TM-score value, RMSD and top five models of the given protein sequence. The produced 3D structure was downloaded into the PDB format, which was chosen based on the C-score value. The server contains a C-score ranging from -5 to 2, where a higher value indicates a protein model with high confidence. The identified 3D structure was submitted into the GalaxyRefine (http://galaxy.seoklab.org/refine) online web-based server for the refinement of the vaccine structure. This webserver was run by the CASP10 refine technique (Nugent et al., 2014) . The GalaxyRefine website provides the RMSD, energy score and overall quality score. The refined structure was downloaded, and the selected structure was identified depending on the energy score of the lowest and highest RMSD value. The refined and identified structure was visualized using the PyMOL v2.3.4 software (DeLano, 2002) . The resulted 3D structure was evaluated depending on the Ramachandran plot score (vaccine structure validity) and Z-score value that determine the standard deviations from the mean value (Z-score within the known native protein range indicating the good quality of the prepared model). The Ramachandran plot was analyzed by the Rampage server (http://mordred.bioc.cam.ac. uk/$rapper/rampage.php), which runs considering allowed and disallowed regions of amino acid (Lovell et al., 2003; Ramachandran et al., 1963) ; and Z-score plot was analyzed by the ProSA-web (https://prosa.services.came.sbg.ac.at/prosa. php) tool (Wiederstein & Sippl, 2007) .

Molecular docking studies can reveal the binding interactions between modeled protein and receptor molecules. For this purpose, we submitted the refined vaccine model as ligand and TLR4 protein as immunological receptor into the ClusPro v2.0 server, available at https://cluspro.bu.edu/, for molecular docking (Kozakov et al., 2017) . The TLR4 receptor (PDB ID: 4G8A) was selected and downloaded from the PDB server (Sussman et al., 1998) . Initially, the receptor was prepared by separating the attached ligand from the protein, followed by the removal of waters and other chemicals. All these processes were performed in PyMOL v2.3.4 software (DeLano, 2002) . Binding interactions and residues involved in the interacting plane were analyzed with Discovery Studio 2017.

For molecular dynamics (MD) simulation, we used both software and server-based tools to evaluate the dynamics and stability of the vaccine-receptor complex critically. The stability of the docked complex was evaluated by a highly intuitive and accurate molecular dynamic simulation tool YASARA, where parameters for macromolecules to facilitate simulations were generated using the AMBER14 force field . We evaluated the stability, fluctuation and compactness of the vaccine-receptor complex in terms of root mean square deviation (RMSD), root mean square fluctuation (RMSF), and radius of gyration (R g ) values, respectively. Also, the complex was submitted to the iMODS server, available at http://imods.chaconlab.org/ (L opez- Blanco et al., 2014) . Based on the normal mode analysis (NMA), this server provides eigenvalues, deformability, B-factors, and elastic network model to clarify the aggregate protein movement in the inside directions.

To evaluate the possible immune response of the vaccine, the whole construct was submitted into C-IMMSIM v10.1 server available at http://www.cbs.dtu.dk/services/C-ImmSim-10. 1/, and the generated responses were retrieved for detailed observation (Rapin et al., 2010) . In this case, we considered a minimum interval period of 30 days between two doses, as described earlier (Castiglione et al., 2012) . In silico administration of three injections were given with time steps of 1, 84, and 168, respectively, where one-time step is equal to 8 h in real life. The maximum value for simulation steps was set to 300, while the rest of the stimulation parameters were kept default.

For the expression of a foreign gene in a host organism, codon optimization is necessary according to the specific host organism (Grote et al., 2005) . Therefore, the construct was submitted into the JCat server (http:/jcat.de/) for the codon adaptation. Herein, we considered widely used E. coli K12 as the host, and the whole operation is carried out by avoiding the following three criteria: (1) restriction enzymes cleavage sites, (2) binding sites of the prokaryotic ribosome and (3) rho-independent termination of transcription. The adapted sequence was evaluated based on the codon adaptation index (CAI) value and guanine-cytosine (GC) content (Grote et al., 2005) . Finally, the adapted nucleotide sequence was used for in silico cloning into the pET28a (þ) expression vector. The whole in silico cloning operation was executed in SnapGene v4.2 software (Goldberg et al., 2018) .

We found 250 S-proteins from all the retrieved SARS-CoV-2 proteomes. Based on the antigenicity, we selected a spike protein with an antigenic score of 0.4646 (VaxiJen) and 0.717 (ANTIGENpro), which was the highest among all tested proteins. The length of the selected S-protein was 1273 amino acids long while the GenBank accession was QIC53213. The primary sequence of the selected protein was used for further analysis.

A total of 270 CTL epitopes, each with a length of 9 amino acids, were predicted from the selected spike protein. Assessment revealed that twenty-nine CTL epitopes were antigenic, immunogenic, non-toxic and non-allergenic (Supplementary  Table S1 ). Due to the high number of potential epitopes, we selected the top four CTL epitopes for the final vaccine construction based on the antigenicity score (Table 1) .

A total of 478 HTL epitopes, each with a length of 15 amino acids, were identified initially using the IEDB server. Among them, only 16 HTL epitopes were able to induce the evaluated three types of cytokines, such as IFNc, IL4, and IL10 (Supplementary Table S2 ). Likewise, we considered the top four HTL epitopes for incorporating into the final vaccine construct based on the antigenic score (Table 2) .

Preliminary analysis revealed a total of 61 LBL epitopes, each with a length of 12 amino acids. Later with further evaluation, 14 epitopes were found as antigenic, non-toxic and non-allergenic (Supplementary Table S3 ). Out of 14 LBL epitopes, we selected the top four LBL epitopes for vaccine design purposes based on the antigenicity score (Table 3) .

The selected CTL and HTL epitopes, along with their binding alleles, were used to evaluate the population coverage, as shown in Figure 2 . Both CTL and HTL epitopes provided a high percentage (93.30%) of population coverage across the world. The selected epitopes showed interactions with a high number of HLA alleles from different countries such as the United States (99.38%), North America (99.35%), South Korea (99.14%), South Asia (99.10%) and India (99.05%). This result suggests that the vaccine designed with these epitopes could be effective on most of the population in the world (Figure 2 and Supplementary Table S4 ).

We used the docking method to validate the efficacy of selected epitopes in binding their respective HLA alleles. The epitopes, along with their respective docking allele, binding affinities, interactions and residues involved in the hydrogen bonds, are described in Table 4 . The binding affinities of CTL epitopes were between -7.1 and -9.0 kcal/mol, while for HTL epitopes, it was between -5.8 and -6.9 kcal/mol. The binding affinities were either very close or even higher than that of the positive control (Table 4 ). In addition to the tabulated details, we presented the best interacting CTL (VVFLHVTYV) and HTL (QYIKWPWYIWLGFIA) epitopes in Figure 3 . Herein, the best CTL epitope produced a total of 12 hydrogen bonds, in which 8 were classical interactions involved with the active site residue Lys80, Tyr84, Lys146, Val2, Thr7, Tyr8, Val9, Lys66, Asn77, and Thr143. On the other hand, the best HTL epitope showed nine hydrogen bonds, including six classical interactions while it interacted with Ser53, Glu55, Asn62, His328, Trp7, Ala15, Phe13, Tyr8, Ile14, and Ile3 residues.

The vaccine construct was formulated using the previously selected 12 epitopes belonging to three different classes (4 CTL, 4HTL, and 4 LBL). The epitopes were added together with AAY, GPGPG and KK linkers, respectively, as shown in Figure 4 . An adjuvant was added ahead of the construct to improve the immunogenicity. The TLR4 agonist 50S ribosomal protein L7/L12 was linked to the first CTL epitope as an adjuvant by using EAAAK linker. The final vaccine construct was 316 amino acids long (Figure 4 ).

The physicochemical properties of the vaccine construct were assessed as shown in Table 5 . The molecular weight of the construct was found to be 33,614.95 Da. At the same time, IFNc stands for interferon-gamma while IL-4 and IL-10 indicate the interleukin-4 and interleukin-10, respectively. other properties such as theoretical isoelectric point (pI) was 8.28, chemical formula was C 1539 H 2432 N 390 O 439 S 6 , instability index was 24.33, the aliphatic index was 93.64, and grand average of hydropathicity was 0.035. In addition, physicochemical features and the immunological potency of the construct were evaluated. For instance, the antigenicity of the construct was 0.6166, while immunogenicity was 1.58298. Furthermore, the vaccine was non-allergenic and soluble, with a score of 0.871753 out of 1 (Table 5 ). The secondary structural features include a-helix, b-strand and random coils that were evaluated using two different servers. The SOPMA server predicted 39.56% a-helix, 23.42% b-strand and 37.03% random coils in the construct (Table 6 , Supplementary Figure S1 ). On the other hand, the PSIPRED server anticipated the features as 42.41% a-helix, 10.44% b-strand, and 47.15% random coils (Table 6 , Supplementary Figure S1 ).

In homology modeling, the I-TASSER server used 1DD4 (PDB ID) as the best template to generate the top five models. Among the five models, we considered the model with the lowest C-score (-4.82) as recommended by the server (Supplementary Figure S4) . A structural representation of the designed vaccine is provided in Figure 5 . After refinement, the vaccine (model 3) showed 84.7% residues in the favorable region in the Ramachandran plot, with GDT-HA score 0.9146, RMSD value 0.526, MolProbity 2.597, Clash score 27.9 and Poor rotamers score 0.8 (Supplementary Table S5 ).

The refined 3D vaccine model was further validated with the RAMPAGE and ProSA-web servers. Before refinement, the Ramachandran plot of the vaccine showed 61.1% residues in the favorable region and 27.4% in allowed regions, while . Graphical map of the formulated multi-epitope vaccine construct. The vaccine constructs included (left to right) an adjuvant, CTL, HTL and LBL epitopes are shown in the dark blue, red, olive green and green rectangular boxes. Herein, the adjuvant and the first CTL epitope were linked by EAAAK linker (blue), CTL epitopes were added together by AYY linkers (off-white), HTL epitopes by GPGPG linkers (orange), and LBL epitopes by KK linkers (black).

11.5% residues in disallowed regions (Supplementary Figure  S3) . The Ramachandran plot of the refined vaccine model showed 86.3% residues in the favorable region and 9.6% in allowed regions, while 4.1% residues in disallowed regions (Figure 6(A) ). Likewise, the crude model showed a Z-score value of -7.17, while the refined model provided a value of -7.4 (Figure 6(B), Supplementary Figure S4 ).

The docking between the vaccine (ligand) and TLR4 (receptor) was performed to anticipate their binding affinity and interactions. In doing so, the ClusPro v2.0 server provided 30 docked complexes with different poses. Among them, we selected the complex with the least energy score and binding pose with functional interactions (Supplementary Table   Table S6). Thus, model 1 fulfilled the inclined criteria. Hence, it was picked as the best vaccine-TLR4 complex, which had an energy score of -964.6 ( Figure 7) . The selected complex was analyzed for binding interactions and involved in active site residues. The number of hydrogen and hydrophobic bonds present in the interaction plane was 35 and 9, respectively. Among the hydrogen bonds, 28 were classical hydrogen bonds (CHB). The interacting residues in the CHB from the vaccine were Ala36, Arg283, Asn281, Gln267, Gln285, Glu30, Glu33, Lys3, Met1, Pro184, Pro265, Thr268, Thr279, Thr31, Thr6, Leu186, Lys275 and Ser280. Moreover, associated TLR4 active site residues were Gln115, Asn137, Gln163, Lys186, Lys20, Gln21, Ser45, Asn47, Cys51, Asn114, Lys130, Gln91, His159, Asn26, Glu111, Leu154, Thr112 and Trp23 (Figure 7) . Other hydrogen bond interactions were as follows: four were electrostatic salt bridges, two were carbon-hydrogen bonds and a single Pi-donor hydrogen bond.

We executed software-based MD simulation where trajectories from 10 ns long simulation showed structural stabilization around 6.2 ns and light fluctuation afterward ( Figure  8(A) ). The calculated average RMSD value was 3.25 Å, while the average RMSF score was 2.65 Å. The fluctuation was higher in the vaccine part from AA 1500 to AA 1800 ( Figure  8(B) ). Furthermore, the average simulation energy was -7,229,534.41 kJ/mol. The average R g score was 37.62 that fluctuates between 30.90 and 44.07 (Figure 8(C) ). The relatively higher pick at 1500 to 1800 amino acid residues was due to the flexible regions of the docked complex, which was the vaccine part. The MD simulation was also carried out in the iMODS server, where NMA assessment was applied to the internal coordinates of the complex. The deformability builds up the independent distortion of each residue portrayed by the method of chain hinges (Figure 8(D) ). The eigenvalue determined for the complex was found to be 1.871 e-05 (Figure 8(E) ). The variance of each typical mode was gradually decreased (Figure 7 (F)). All these results suggest stable binding interactions with compact conformation and minor fluctuations in the vaccine-TLR4 complex.

The simulated immune response showed similar to actual immunological phenomena provoked by specific pathogens as shown in Figure 9 . For instance, secondary and tertiary immune responses were higher than the primary immune response (Figure 9(A) ). Secondary and tertiary responses showed higher levels of antibodies (i.e. IgG1 þ IgG2, IgM, and IgG þ IgM), which coincided with an antigen extenuation indicating the development of memory cells, thus, intensified antigen clearance upon successive exposures (Figure 9(A) ). Additionally, a prolonged period of viability in B-cells, cytotoxic T-cells and helper T-cells were noticed, indicating the class switching between immune cells and IgM memory formation (Figure 9(B-D) ). The elevated levels of IFNc, IL-4 and IL-10 were additionally apparent (Figure 9 (E)). The percentage (%) and amount (cells/mm 3 ) of Th0 type immune reaction were lower than the Th1 type reaction (Figure 9 (F)). During the presentation, expanded macrophage movement was illustrated, while dendritic cell movement was predictable ( Figure 9(G,H) ). 

We optimized the codons present in the vaccine construct according to the E. coli K12 in the JCat server to increase their translation efficiency. The peptide vaccine construct (316 AA residues) produced 948 lengths of nucleotide sequences. Moreover, the adapted nucleotide sequence has GC content and CAI value of 59.04% and 1.0, respectively. To insert the adapted sequence into the pET28a (þ) vector, we selected XhoI and BamHI restriction sites as the start and end cut points, respectively. Thus, the optimized vaccine construct was cloned into the pET28a (þ) cloning vector with the SnapGene software ( Figure 10 ). The final size of the cloning vector was 6281 nucleotide base pairs (bp).

The present demonic appearance of COVID-19 generates a life-threatening situation to the global public health (Vankadari & Wilce, 2020) , which influences us to design this multi-epitope vaccine applying immunoinformatics approach.

The glycoprotein-based vaccine demonstrated an extraordinary significance ordained by immunoinformatics and revealed our attempt trustworthy. A vaccine is a safe and effective way to protect against infectious diseases (Li et al., 2014) . It should have the ability to provide acquired immunity against contagious diseases (Bol et al., 2016) . In this study, we designed an epitope-based vaccine that could provide a strong immune response against SARS-CoV-2, thereby, preventing the COVID-19 pandemic. A vaccine can prevent future outbreaks (Melief et al., 2015) . However, in the absence of an effective vaccine, control and prevention of COVID-19 infection and transmission are very difficult. Besides, effective vaccination is yet to be developed in controlling the current situation. Thus, a new strategy of vaccine development is a prime need that will contribute to finding a solution to solve this present life-threatening public health issue. Since S-protein of SARS-CoV-2 plays a major role in immune invasion as well as the human to human transmission, our purpose was to design an epitope vaccine by targeting the S-protein. The location of the antigenic region of Figure 7 . Molecular docking between the vaccine and the TLR4 receptor. The interacting residues from the vaccine were Ala36, Arg283, Asn281, Gln267, Gln285, Glu30, Glu33, Lys3, Met1, Pro184, Pro265, Thr268, Thr279, Thr31, Thr6, Leu186, Lys275, and Ser280 while associated TLR4 active site residues were Gln115, Asn137, Gln163, Lys186, Lys20, Gln21, Ser45, Asn47, Cys51, Asn114, Lys130, Gln91, His159, Asn26, Glu111, Leu154, Thr112, and Trp23. the surface of S-protein was evaluated so that this protein can be recognized by cellular and humoral immune systems. First, all potential CTL, HTL and LBL epitopes were identified and evaluated. Then, the vaccine was designed with the top four antigenic CTL, HTL and LBL epitopes with their desired linkers. They were incorporated in vaccine construction as a part of the essential element that enhances the stability, folding and expression patterns of our vaccine candidate (Shamriz et al., 2016) . The adjuvant was attached to the CTL epitope by EAAAK linker, which helps to induce high levels of both cellular and immunogenic humoral responses for particular antigens, and amplify the vaccine's stability and longevity (Bonam et al., 2017; Lee & Nguyen, 2015) . Finally, the vaccine construction was found to accumulate 316 amino acid residues long. Solubility, a type of physicochemical property of a vaccine candidate, is counted as a vital characteristic of any recombinant vaccine (Khatoon et al., 2017) . Hence, the solubility of the vaccine construct was predicted by using a solubility assessing tool to determine the quality of being solvability of the construct inside the host E. coli, and the vaccine construct was found to be solvable inside the host E. coli. The nature of the vaccine determined by theoretical PI value was found to be acidic. Instability index suggested by server tools indicate that the protein would remain stable after synthesis. In contrast, the GRAVY value and aliphatic index portrayed the vaccine to be the hydrophilic and thermostable, respectively. A favorable physicochemical property predicted for the vaccine and all the scores on different parameters relies on a high possibility to confer this vaccine as a valid candidate against SARS-CoV-2.

In our approach, we found the best population coverage all-over the world (93.30%) for the combined results that make this vaccine construct a good candidate and significant weapon. The most cases of infection and mortality were found in the United States and the United Kingdom, and we find the parentage of population coverage for those countries in a significant level (United States (99.38%), Europe (98.92%)). After the 3D structure prediction (based on cscore), the identified models were refined and selected the best model (based on the lowest energy score). In the validation test of 3D structure, we found a good number of Zscore (-7.4) and the superior features of most favored, accepted and disallowed regions for the Ramachandran plot. Molecular docking between the peptide vaccine and virus glycoprotein binding favorable receptor of TLR4 with lowest energy score of À964.6 confirmed the possibility of infection inhibitory activity of the vaccine and suggested a possible tight interaction between the modeled vaccine as a ligand and the TLR4 receptor surface. Molecular dynamics simulation is a potential technique for accessing the physical ground of the protein structure and function of biological macromolecules. Protein dynamic simulations can provide certain information regarding individual atomic movement as a function of time. For the dynamics evaluation of the vaccine candidate, 10 ns dynamic simulations have been performed, and results have been analyzed based on the RMSD and RMSF score. RMSD value is used to compare different atomic conformations of a given molecular system. In this study, the RMSD value was used to determine the significant flexibility and departure of vaccine candidates from the receptor structure, where the RMSF of the complex structure was determined to measure the displacement of our particular vaccine candidate's atom relative to the receptor structure. The calculated average RMSD and RMSF value was 3.25 Å and 2.65 Å, respectively. The fluctuation was found higher in the vaccine part, but it smoothly became stable after 600 ps suggesting possible stability of the modeled vaccine and the receptor. Finally, we performed an immune simulation to observe the optimal behavior and cell density parameters for successful target clearance and find the best immunological response against the pathogen. The vaccine doses were upgraded immunological reaction causing memory B-cell (having a half-life of several months) and T-cell. Sustained generation of IFN-c and IL-2 were seen after immunization because of expanded aide T-cell initiation. In this manner, the vaccine effectively simulated a humoral immunological response to increasing immunoglobulin creation. The MD simulation was performed to evaluate the stability of the vaccine candidate with the receptor, where codon optimization was performed to stabilize the construct vaccine within the host for optimum multi-epitope vaccine production. Finally, the codon was optimized, and desired vaccine candidate in silico cloning was performed successfully into the pET28a (þ) cloning vector of E. coli K12 expression host.

In this study, a series of computational approaches led to the discovery of potential T-and B-cell epitopes in S-protein of SARS-CoV-2 that eventually embroidered into a multi-epitope vaccine. The newly designed vaccine has desired immunodominant properties with high population coverage. Importantly, it was able to bind with the immune receptor TLR4 strongly as well as to elicit robust immune response upon SARS-CoV-2 infection. Based on our findings, we believe that the vaccine candidate can be an important starting point for developing a potent vaccine against the etiological agent of COVID-19 outbreak. Moreover, the potential epitopes identified in this study can be used in future studies as well. However, further experimental assessments are required to confirm our formulated vaccine as an effective prophylactic against SARS-CoV-2.

In silico study the inhibition of Angiotensin converting enzyme 2 receptor of COVID-19 by Ammoides verticillata components harvested from

The in silico cloning of the designed vaccine into the pET-28a (þ) vector. Herein, purple color represents the vector DNA, while the red color indicates the adapted DNA sequence of the designed vaccine. western Algeria

Simultaneous cognate epitope recognition by bovine CD4 and CD8 T cells is essential for primary expansion of antigen-specific cytotoxic T-cells following ex vivo stimulation with a candidate Mycobacterium avium subsp. paratuberculosis peptide vaccine

Contemporary strategies and current trends in designing antiviral drugs against dengue fever via targeting hostbased approaches

Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2): An overview of viral structure and host response

The Protein Data Bank

Prophylactic vaccines are potent activators of monocytederived dendritic cells and drive effective anti-tumor responses in melanoma patients at the cost of toxicity

An overview of novel adjuvants designed for improving vaccine efficacy

Novel 2019 coronavirus structure, mechanism of action, antiviral drug promises and rule out against its treatment

Novel nested peptide epitopes recognized by CD4þ T cells induced by HIV-1 conserved-region vaccines. Vaccines

Scalable web services for the PSIPRED Protein Analysis Workbench

Predicting population coverage of T-cell epitope-based diagnostics and vaccines

Properties of MHC class I presented peptides that enhance immunogenicity

Statement in support of the scientists, public health professionals, and medical professionals of China combatting COVID-19

How the interval between prime and boost injection affects the immune response in a computational model of the immune system

A familial cluster of pneumonia associated with the 2019 novel coronavirus indicating person-to-person transmission: A study of a family cluster

An investigation into the identification of potential inhibitors of SARS-CoV-2 main protease using molecular docking study

Pymol: An open-source molecular graphics tool

Prediction of IL4 inducing peptides

Designing of interferongamma inducing MHC class-II binders

AllerTOP -A server for in silico prediction of allergens

Vaccinomics approach for developing multi-epitope peptide pneumococcal vaccine

VaxiJen: A server for prediction of protective antigens, tumour antigens and subunit vaccines

SARS-CoV-2 RNA dependent RNA polymerase (RdRp) targeting: An in silico perspective

Novel guanosine derivatives against MERS CoV polymerase: An in silico perspective

Reverse vaccinology approach to design a novel multi-epitope vaccine candidate against COVID-19: An in silico study

Immunoinformatics approach for multiepitope vaccine prediction from H, M, F, and N proteins of Peste des Petits ruminants virus

Protein identification and analysis tools on the ExPASy server

SOPMA: Significant improvements in protein secondary structure prediction by consensus prediction from multiple alignments

Salmonella persist in activated macrophages in T cell-sparse granulomas but are contained by surrounding CXCR3 ligand-positioned Th1 cells

JCat: A novel tool to adapt codon usage of a target gene to its potential expression host

In-silico approaches to detect inhibitors of the human severe acute respiratory syndrome coronavirus envelope protein ion channel

In silico approach for predicting toxicity of peptides and proteins

A review on the cleavage priming of the spike protein on coronavirus by angiotensin-converting enzyme-2 and furin

COVID-19: What is next for public health? The Lancet

Clinical features of patients infected with 2019 novel coronavirus in Wuhan. The Lancet

Discovery of potential multi-target-directed ligands by targeting host-specific SARS-CoV-2 structurally conserved main protease

Targeting SARS-CoV-2: A systematic drug repurposing approach to identify promising inhibitors against 3C-like proteinase and 2 0 -O-ribose methyltransferase

Exploring Leishmania secretory proteins to design B and T cell multi-epitope subunit vaccine using immunoinformatics approach

The ClusPro web server for protein-protein docking

Understanding the binding affinity of noscapines with protease of SARS-CoV-2 for COVID-19 using MD simulations at different temperatures

Large-scale validation of methods for cytotoxic T-lymphocyte epitope prediction

Discovering and understanding oncogenic gene fusions through data intensive computational approaches

Recent advances of vaccine adjuvants for infectious diseases

Therapeutic options for the 2019 novel coronavirus (2019-nCoV)

Peptide vaccine: Progress and challenges. Vaccines

iMODS: Internal coordinates normal mode analysis server

Structure validation by Calpha geometry: Phi, psi and Cbeta deviation

SOLpro: Accurate sequencebased prediction of protein solubility

High-throughput prediction of protein antigenicity using protein microarray data

iBCE-EL: A new ensemble learning framework for improved linear Bcell epitope prediction

Therapeutic cancer vaccines

Computational studies of drug repurposing and synergism of lopinavir, oseltamivir and ritonavir binding with SARS-CoV-2 protease against COVID-19

Computer-aided designing of immunosuppressive peptides based on IL-10 inducing potential

Proteome-wide screening for designing a multi-epitope vaccine against emerging pathogen Elizabethkingia anophelis using immunoinformatic approaches

Structural basis and designing of peptide vaccine using PE-PGRS family protein of Mycobacterium ulcerans -An integrated vaccinomics approach

Evaluation of predictions in the CASP10 model refinement category

Toll-like receptor 4 in acute viral infection: Too much of a good thing

Novel immunoinformatics approaches to design multi-epitope subunit vaccine for malaria by investigating anopheles salivary protein

Transmission routes of 2019-nCoV and controls in dental practice

ViPR: An open bioinformatics database and analysis resource for virology research

Stereochemistry of polypeptide chain configurations

Computational immunology meets bioinformatics: The use of prediction tools for molecular binding in the simulation of the immune system

I-TASSER: A unified platform for automated protein structure and function prediction

In-silico homology assisted identification of inhibitor of RNA binding against2019-nCoV N-protein (N terminal domain)

Effect of linker length and residues on the structure and stability of a fusion protein with malaria vaccine application

The outbreak of SARS-CoV-2 pneumonia calls for viral vaccines

Diagnosis, treatment, and prevention of 2019 novel coronavirus infection in children: Experts' consensus statement

An in-silico evaluation of different Saikosaponins for their potency against SARS-CoV-2 using NSP15 and fusion spike glycoprotein as targets

Protein Data Bank (PDB): Database of three-dimensional structural information of biological macromolecules

AutoDock Vina: Improving the speed and accuracy of docking with a new scoring function, efficient optimization, and multithreading

Identification of new anti-nCoV drug chemical compounds from Indian spices exploiting SARS-CoV-2 main protease as target

Emerging WuHan (COVID-19) coronavirus: glycan shield and structure prediction of spike glycoprotein and its interaction with human CD26. Emerging Microbes & Infections

Stilbene-based natural compounds as promising drug candidates against COVID-19

A novel coronavirus outbreak of global health concern

Clinical characteristics of 138 hospitalized patients with 2019 novel coronavirus-infected pneumonia in Wuhan

Peptide binding predictions for HLA DR, DP and DQ molecules

ProSA-web: Interactive web service for the recognition of errors in three-dimensional structures of proteins

Cryo-EM structure of the 2019-nCoV spike in the prefusion conformation

A new coronavirus associated with human respiratory disease in China

Evolution of the novel coronavirus from the ongoing Wuhan outbreak and modeling of its spike protein for risk of human transmission

Pathological findings of COVID-19 associated with acute respiratory distress syndrome. The Lancet

A novel coronavirus from patients with pneumonia in China

We extend our gratitude towards the Deanship of Scientific Research (DSR) at King Abdulaziz University and Biological Solution Centre (BioSol Centre) for providing technical support. Special thanks go to Monokesh Kumer Sen, School of Medicine, Western Sydney University, Australia for his contribution in revising the whole manuscript.

The authors declare no conflict of interest.