key: cord-1025337-vybv2srj authors: Akbay, Burkitkan; Abidi, Syed Hani; Ibrahim, Mahmoud A. A.; Mukhatayev, Zhussipbek; Ali, Syed title: Multi-Subunit SARS-CoV-2 Vaccine Design Using Evolutionarily Conserved T- and B- Cell Epitopes date: 2021-06-26 journal: Vaccines (Basel) DOI: 10.3390/vaccines9070702 sha: a8fde4a7eb9fb7487edaa82973d496af2cfb502e doc_id: 1025337 cord_uid: vybv2srj The SARS-CoV-2 pandemic has created a public health crisis worldwide. Although vaccines against the virus are efficiently being rolled out, they are proving to be ineffective against certain emerging SARS-CoV-2 variants. The high degree of sequence similarity between SARS-CoV-2 and other human coronaviruses (HCoV) presents the opportunity for designing vaccines that may offer protection against SARS-CoV-2 and its emerging variants, with cross-protection against other HCoVs. In this study, we performed bioinformatics analyses to identify T and B cell epitopes originating from spike, membrane, nucleocapsid, and envelope protein sequences found to be evolutionarily conserved among seven major HCoVs. Evolutionary conservation of these epitopes indicates that they may have critical roles in viral fitness and are, therefore, unlikely to mutate during viral replication thus making such epitopes attractive candidates for a vaccine. Our designed vaccine construct comprises of twelve T and six B cell epitopes that are conserved among HCoVs. The vaccine is predicted to be soluble in water, stable, have a relatively long half-life, and exhibit low allergenicity and toxicity. Our docking results showed that the vaccine forms stable complex with toll-like receptor 4, while the immune simulations predicted that the vaccine may elicit strong IgG, IgM, and cytotoxic T cell responses. Therefore, from multiple perspectives, our multi-subunit vaccine design shows the potential to elicit a strong immune-protective response against SARS-CoV-2 and its emerging variants while carrying minimal risk for causing adverse effects. The novel severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), has caused a worldwide pandemic that continues to negatively impact the global economy, social dynamics, and health care systems [1] . SARS-CoV-2 has been identified as the seventh coronavirus that can infect humans [2] . Although coronaviruses have infected the human population repeatedly over the past century [3] , only a few have caused severe large-scale outbreaks. Prior to SARS-CoV-2, the SARS-CoV outbreak occurred in southern China in November 2002, followed by the Middle East respiratory syndrome (MERS-CoV) outbreak in Saudi Arabia and South Korea, respectively, in 2012 and 2015. SARS-CoV and MERS-CoV have an estimated mortality rate of~9.6% and 35%, respectively [4, 5] . Moreover, four other human coronaviruses (HCoVs) associated with mild upper respiratory symptoms, have previously been identified, i.e., HCoV-NL63, HCoV-229E, HCoV-OC43, and HCoV-HKU1 [6] [7] [8] . for each HCoV exhibiting >90% sequence similarity to the reference sequences. The retrieved sequences were aligned using the online MAFFT alignment tool (https://mafft.cbrc. jp/alignment/server/ (accessed on 23 May 2020)). In order to analyze the evolutionary relationship between SARS-CoV-2 with other HCoVs (SARS-CoV, MERS-CoV, HCoV-NL63, HCoV-229E, HCoV-OC43, and HKU-1), nucleotide alignment was used to construct a maximum-likelihood (ML) phylogenetic tree using IQ-TREE v1.6.12 [25] . The ML tree was constructed using the general time reversible (GTR) model of nucleotide substitution and gamma correction for among site rate variation. The ultrafast bootstrap method was used with 1000 replicates to infer support branching in the tree topology. The Newick file for the ML tree was visualized using FigTree v1.4.4. In the next step, E, N, M, and S reference sequences for all seven HCoVs (described above) were used in BLAST search (pBLAST), with filters of 90-100% percent identity and 100% query coverage, to retrieve E, N, M, and S protein sequences from HCoVs that are over 90% conserved. The total number of protein sequences from each HCoV used in our study are shown in Table 1 . Subsequently, the E, N, M, and S protein sequences for each HCoV were separately aligned using an online MAFFT alignment tool [26] , and each alignment was used to generate protein specific consensus sequences for each HCoV with the Consensus Maker tool (https://www.hiv.lanl.gov/content/sequence/CONSENSUS/SimpCon.html (accessed on 25 May 2020)). In the next step, protein-specific consensus sequences from the HCoVs were aligned again, and the conserved regions within the E, N, M, and S proteins from the seven HCoVs were identified using ESPript 3.x [27] and Jalview 2.11.1.3 [28] software. Additionally, to take into account the recently emerging mutations among SARS-CoV-2 variants, approximately 5000 whole genome sequences of three SARS-CoV-2 variants:B.1.1.7 (VOC-202012/01) identified in the U.K., 501Y.V2 (B.1.351) identified in the South Africa, and 501Y.V3 (P.1) identified in Brazil were retrieved from GSAID (https://www.gisaid.org/ (accessed on 5 April 2021), and aligned with MAFFT to generate consensus sequences for each of the variants using BioEdit [29] . Nucleotide consensus sequences were then translated into amino acid sequences using the Expasy tool (https://web.expasy.org/translate/ (accessed on 5 April 2021)). The resulting amino acid sequences were trimmed with reference to E, N, M, and S reference protein sequences in Aliview [30] . The protein-specific consensus sequences for each HCoV, including the SARS-CoV-2 variants, were used to predict B and T cell epitopes and to design a multi-subunit vaccine based on the selected epitopes as described below (Figure 1 ). Using E, N, M, and S protein consensus sequences, HCoV-derived cytotoxic T cell epitopes were predicted using the web-based CTLPred tool [31] and a combined support vector machine and artificial neural network approach. Human leukocyte antigen (HLA) restriction was identified using the nHLAPred tool for each predicted T cell epitope [32] . The prediction of linear B cell epitopes was performed using the web-based ABCpred tool [33] . Epitopes located in the conserved regions between the 7 HCoVs or 3 closely related viruses, namely, SARS-CoV-2, SARS-CoV, and MERS-CoV, were selected. Predicted CD8+ T and B cell epitopes were compared to the experimentally validated SARS-CoV epitopes available in IEDB (http://www.iedb.org/ (accessed on 10 April 2021)) database [34] and that had been reported by other studies. The predicted epitopes were shortlisted after screening for immunogenicity and toxicity. VaxiJen v2.0 [35] was used to evaluate immunogenicity (a score >0.4 (default value) was considered as antigenic), while ToxinPred [36] was used to evaluate the toxicity of the epitopes. Epitopes predicted to be immunogenic and non-toxic were selected to design the multi-subunit vaccine. As indicated earlier, the conserved epitopes predicted to be immunogenic and nontoxic were used to design a multi-subunit vaccine (Supplementary Figure S1 ). For the vaccine design, the CD8+ T cell epitopes were linked using an AAY linker, and linear B cell epitopes were linked by an GPGPG linker [37] . At the C-terminal of the vaccine, a six-histidine tail was added while at the N-terminus, β-defensin amino acid sequence (UniProt ID: Q5U7J2) was added to the first epitope via an EAAAK linker. β-defensin acts as an adjuvant and mediates recruitment of naïve T cells and immature dendritic cells at the site of infection/immune activity [38] . Physico-chemical properties, such as half-life, instability index, aliphatic index, and grand average of hydrophobicity of the vaccine were determined using the ProtParam server (http://web.expasy.org/protparam (accessed on April 12 2021)), while the basic property of allergenicity of the multi-subunit vaccine was assessed using AllerTOP v. 2.0 server [39] . The immunogenicity/antigenicity of the multi-subunit vaccine was evaluated using the VaxiJen v2.0 server [35] . The secondary structure of the designed vaccine was predicted using the PSIPRED server [40] , while the 3D structure was predicted using both the I-TASSER [41] and the Phyre2 server (using an intensive modelling approach) [42] . The best model was selected and refined/optimized using the Galaxyrefine server [43] , while PROSA [44] and Ramachandran plot analyses were performed to analyze the quality of the refined structure. Solubility of the protein (vaccine) was determined using the Protein-Sol server [45] . Since Toll-like receptor 4 (TLR4) is capable of recognizing various exogenous and endogenous ligands such as viral protein, lipopolysaccharide, heat shock proteins, βdefensin, etc., which leads to dendritic cell maturation, and up-regulation of costimulatory molecules resulting in type 1 polarized adaptive immune response [46, 47] , we performed molecular docking between TLR4 and our vaccine-construct (containing β-defensin as an adjuvant). For molecular docking calculations, the crystal structure of the human TLR4 deposited in the Protein Data Bank with the PDB accession code of 4G8A [48] was retrieved and utilized as a template. For the preparation of the TLR4 crystal structure, all crystallographic water molecules and heteroatoms were removed. Moreover, mutated residues were reversed and refined. The protonation state of the TLR4 was then examined, and hydrogen atoms were consequently added [49] . Binding modes of the developed vaccine with the human TLR4 were then predicted using the PatchDock server [50] . All predicted docked complexes (exactly 1000 solutions) were then refined and rescored using the FireDock server [51, 52] . The interactions between the vaccine and the TLR4 were resolved using the LigPlot + software [53] . The humoral and cellular immune response against the multi-subunit vaccine construct was simulated using the C-ImmSim tool [54] . For simulations, the following parameters were used: three injections at 4-week intervals; simulation steps = 1000; time to simulation = default; injection volume = vaccine only and no LPS. Additionally, we used NetCTL [55] and NetChop [56, 57] to predict proteasomal sites (using default a value of 0.5) and the possible proteolytic release of the epitopes from the vaccine construct. To evaluate the genetic similarity among the seven human coronaviruses, we constructed ML phylogenetic trees that revealed alpha-HCoVs, NL63 and 229E, clustered distantly from Beta-HCoVs, OC43 and HKU1. Importantly, SARS-CoV-2, SARS-CoV, and MERS-CoV clustered together, indicating a close phylogenetic relatedness between these viruses (Supplementary Figure S2) . In the next step, we identified conserved regions in the consensus sequences of the E, N, M, and S proteins from seven HCoVs (Figure 2A and Supplementary Figure S3 ) and three phylogenetically closely related HCoVs-SARS-CoV-2, SARS-CoV, and MERS-CoV ( Figure 2B and Supplementary Figure S4) using ESPript 3.x and Jalview 2.11.1.3 software [28] . The results revealed several regions within N, M, and S proteins exhibiting more than 40% homology among the seven HCoVs, while no conserved regions were found for the E protein (Supplementary Figure S3) , corroborating the findings of the phylogenetic analysis. Additionally, while the conservation was analyzed for SARS-CoV-2, SARS-CoV, and MERS-CoV, the results showed greater sequence conservation for E, N, M, and S compared to those in the seven hCoVs (Supplementary Figures S3 and S4 ). In the next steps, focusing on these conserved regions, we predicted CD8+ T and B cell epitopes for each SARS-CoV-2 protein. CD8+ T cell epitopes were analyzed for E, N, M, and S proteins, focusing on the conserved regions of the seven HCoVs protein sequences. No CD8+ T cell epitopes were identified in the conserved regions of E protein among all seven HCoVs (Supplementary Figure S3 and Table 2 ). However, when the analysis was carried out using E sequences from only SARS-CoV-2, SARS-CoV and MERS-CoV, two conserved CD8+ T cell epitopes were identified, each exhibiting~44% and~77% conservation (Table 3 and Supplementary Figure S4 ). These epitopes were found to be restricted by HLA-C*0401 HLA-A*0205, HLA-A*0301, HLA-A2, HLA-A3, and HLA-B*5301 (Table 3) . Similarly, analysis of the M protein sequences showed the presence of one CD8+ T cell epitope with~55% conservation among all seven HCoVs (Table 2) , which was found to be restricted by HLA-A24, HLA-A*2402, and HLA-Cw*0401 ( Table 2) . Conservation of this epitope between SARS-CoV-2, SARS-CoV, and MERS-CoV was found to be 66% (Table 2 and Supplementary Figure S4 ). Further analysis of SARS-CoV-2, SARS-CoV, and MERS-CoV epitopes revealed three conserved CD8+ T cell epitopes exhibiting~44%,~66%, and~88% conservation between the three viruses (Table 3 and Supplementary Figure S4 ). These epitopes were found to be restricted by HLA-C*0401 and then by HLA -A*0203, -A*24, -B*2703, and -B*35 (Table 3) . Analysis of the N protein sequences showed one CD8+ T cell epitope with 44% conservation among all HCoVs (Table 2 and Supplementary Figure S3 ). This epitope was found to be restricted by HLA-Cw*0401, HLA-A*0202, HLA-A*0203, HLA-A*0205, HLA-A*0301, HLA-A*0401, HLA-B*35, and HLA-B*3501 (Table 2) . Again, conservation of this epitope was increased to~77% when the analysis was focused only on SARS-CoV-2, SARS-CoV, and MERS-CoV (Table 3 and Supplementary Figure S4 ). Seven additional conserved CD8+ T cell epitopes were identified for the N protein with~77%,~66%, and 55% conservation (Table 3 and Supplementary Figure S4 ). Most of the conserved epitopes were found to be restricted by HLA-C*0401 and then by HLA-A*0201, HLA-A*2402, HLA-A*A2, HLA-A*6801, HLA-B*3501, HLA-B*51, and HLA-B*5301 (Table 3 ). Table 2 . Predicted CD8+ T cell epitopes found to be conserved in the seven HCoVs. Conserved amino acids are shown in bold letters. HLA restriction for each epitope is also given. Table 3 . Predicted CD8+ T cell epitopes found to be conserved SARS-CoV-2, SARS-CoV, and MERS-CoV. Conserved amino acids are shown in bold letters. HLA restriction for each epitope is also given. Table 3 . Cont. Table 3 . Cont. Table 3 . Cont. Unlike E, N, and M proteins, two conserved CD8+ T cell epitopes were found for the S protein among all the HCoVs, with~44% and~66% conservation, respectively (Table 2 and Supplementary Figure S3 ). These epitopes were found to be mainly restricted by HLA-C*0401, HLA-A*0206, and HLA-A*3301 ( Table 2 ). Conservation of these two epitopes increased to~66% and~77% in the SARS-CoV-2, SARS-CoV, and MERS-CoV sequences (Table 3 and Supplementary Figure S4 ). Further analysis of epitopes in the S sequences from these three HCoVs revealed 12 conserved CD8+ T cell epitopes, with~44%,~55%, 66%,~77%, and~88% conservation (Table 3 and Supplementary Figure S4 ). Most of the conserved epitopes were found to be restricted by HLA-C*0401, HLA-A1, HLA-A2, HLA-A*0206, HLA-A24, HLA-A*2402, HLA-A*3301, HLA-A*6801, HLA-B*3501, HLA-B8, HLA-B27, HLA-B*2706, and HLA-B51 (Table 3) . It is important to note that all CD8+ T epitopes predicted in the E, N, M, and S sequences were found to be 100% conserved in the latest deposited 5000 SARS-CoV-2 sequences as well as in the three variants of SARS-CoV-2 (Supplementary Figure S5 ). E, N, M, and S consensus protein sequences of seven HCoVs were used to predict linear B cell epitopes of varying length (between 10 to 20 amino acids). Even though the protein sequences showed several regions exhibiting more than 40% homology among the seven HCoVs, no B cell epitopes were found in the conserved regions. Analysis of SARS-CoV-2, SARS-CoV, and MERS-CoV, however, revealed conserved B cell epitopes among these three viruses. Two B cell epitopes, with~42% and~50% conservation, were predicted for the E protein, while only one epitope, with~86% conservation, was predicted for the M protein (Table 4 , Supplementary Figure S4) . Similarly, for the N proteins, four epitopes were found to be conserved in the N protein with~50%,~56%, and~69%, while five epitopes were found to be conserved in the S protein, showing~50%,~56%, and~63% conservation (Table 4 , Supplementary Figure S4 ). Table 4 . Predicted B cell epitopes found to be conserved in E, N, M, and S protein sequences of SARS-CoV-2, SARS-CoV, and MERS-CoV. Conserved amino acids are in bold letters. Envelope № It is again important to note that all the predicted B cell epitopes for E, N, M, and S were found to be 100% conserved in the latest deposited 5000 SARS-CoV-2 sequences as well as in the three variants (B.1.1.7, B.1.351, and P.1) of SARS-CoV-2 (Table 5 , Supplementary Figure S5 ). Table 5 . List of immunogenic and non-toxic T and B cell epitopes from E, N, M, and S proteins. The table shows sequences of the shortlisted epitopes, their immunogenicity and toxicity score, and conservation (including level of conservation) between seven HCoVs and the three closely related SARS-CoV-2, SARS-CoV, and MERS-CoV as well as in three variants of SARS-CoV-2. The last column shows whether these epitopes have previously been experimentally validated or not. In the next step, the immunogenicity and toxicity profile of all the predicted epitopes found to be conserved in the E, N, M, and S proteins of SARS-CoV-2 and other HCoVs was evaluated. Among the CD8+ T cell epitopes, one out of two epitopes for E, two out of four epitopes for M, three out of eight epitopes for N, and six out of fourteen epitopes for S were predicted to be immunogenic and non-toxic, with the highest immunogenic score (1.5) for the peptides of S, followed by M, N, and E (Table 5 ). From B cell epitopes, two out of two epitopes for E, three out of four epitopes for N, and one of five epitopes for the S protein were predicted to be immunogenic and non-toxic, with the highest immunogenic score (1.8) for the peptides of N, followed by S and E (Table 5 ). B cell epitope for M, however, was predicted to be non-immunogenic. Therefore, the B cell epitope for M was not selected for further analysis. A total of 18 (12 CD8+ T cell and 6 linear B cell) immunogenic, non-toxic, and conserved epitopes were subsequently used to construct a multi-subunit vaccine. Physico-chemical assessment of the vaccine construct revealed the estimated half-life of the vaccine to be 30 h in mammalian reticulocytes (in vitro), >20 h in yeast (in vivo), and >10 h in Escherichia coli (in vivo). The instability index was predicted to be 33.06, classifying the vaccine as stable. The aliphatic index and grand average of hydropathicity (GRAVY) of the vaccine was found to be −0.250, making the construct hydrophilic and capacitating it to interact better with the neighboring water molecules [37] . The molecular weight of the construct was found to be 34,336.60 daltons. The solubility index of the vaccine construct, based on Protein-Sol calculation, was found to be 0.34, indicating an acceptable solubility profile for the vaccine. Allergenicity and immunogenicity prediction revealed the construct to be a probable antigen (VaxiJen score 0.56) and non-allergen. Secondary structure assessment of the vaccine construct showed that the construct had 11 helices and 11 strands (Figure 3 ). The 3D structure of the vaccine was constructed using both the Phyre2 and I-TASSER tools. The structure constructed using I-TASSER had an overall low-quality score after refinement, and on a Ramachandran plot, only 88% residues were in the favorable region (data not shown). To the contrary, the 3D structure predicted by the Phyre2 tool after refinement ( Figure 4A ) exhibited a high confidence score, ranging from 84-96% for most regions, and 66-74% for some regions. The Ramachandran analysis of the structure showed that the predicted structure had 94.3% of the residues in a highly preferred region, 4.5% in a preferred region, and only 1.1% residues in a nonpreferred/disallowed region ( Figure 4B ). The PROSA score of the structure predicted by Phyre2 was −2.94, which lies in the acceptable range ( Figure 4C) . Overall, the assessment of the different protein quality parameters indicates the structure predicted by Phyre2 to be of high quality. To reveal the binding modes and energies of the developed vaccine with human TLR4, the molecular docking technique was utilized. All possible vaccine-TLR4 solutions (exactly 1000 solutions) were first predicted using the PatchDock server [50] . The predicted vaccine-TLR4 complex solutions were then refined, and the global binding energies were rescored using the FireDock server [51, 52] . The 3D representation of the predicted vaccine-TLR4 complex structure is depicted in Figure 5 , and the estimated global binding energy and energy contributions for the top-ranked vaccine-TLR4 complex are listed in Table 6 . Interactions between the vaccine-TLR4 complex were resolved using LigPlot+ software. The global binding energy estimates suggest that our vaccine construct formed a stable complex with TLR4, with global binding energy of −47.76 kcal/mol (Table 6 and Figure 5A ). Assessment of the energy contributions revealed that attractive van der Waals (vdW) was the dominant force in vaccine-TLR4 binding energy, with a value of −36.33 kcal/mol ( Table 6 ). Analysis of the protein-protein interaction suggested that the construct formed multiple interactions with TLR4, especially with several residues of β-defensin (adjuvant used in our construct, and a known agonist of TLR4 [46, 47] were also found to interact with TLR4 ( Figure 5B ). showing T cell epitopes in red, while B cell epitopes in orange, (B) Ramachandran plot for the 3D structure of the vaccine showing 97.8% of the residues in the favorable region, (C) ProSA z-score for the 3D structure of the vaccine, determined using NMR spectroscopy (dark blue) and X-ray crystallography (light blue). We performed in silico immune simulations to examine the potential of our vaccine construct to elicit potent and durable B (including antibody response) and T cell responses. Results from proteasomal processing showed that the proteasomal sites in the vaccine are strategically located to release all 9-mer T cell epitopes and variable-length B cell epitopes (Supplementary Figure S2) . Results from the immune simulation showed that the vaccine construct was predicted to elicit significant titers of IgG and IgM ( Figure 6A ), significant effector and memory B-cell response ( Figure 6B ) as well as cytotoxic T cell and interferon-gamma response ( Figure 6C) . Overall, the vaccine construct was predicted to be stable, non-toxic and highly immunogenic, and capable of mounting a sustained protective immune response. Using immunoinformatic tools, we have designed, constructed, evaluated, and assessed the performance of a non-allergenic and antigenic chimeric multi-epitope vaccine that combines T and B cell epitopes from conserved regions of the E, M, N, and S proteins of seven human coronaviruses. The vaccine shows potential to elicit a strong response against SARS-CoV-2 while also providing immune protection against emerging SARS-CoV-2 variants and other HCoVs. All of the currently available vaccines target the SARS-CoV-2 S protein, which is the most antigenic, highly variable, and prone to mutations [19, 20] . Due to the high substitution/mutation rate and rapid transmissions, new variants of SARS-CoV-2 have been continuously emerging globally, posing threats to circumvent the efficacy of current vaccines, possibly leading to a novel-SARS-CoV-2-variant-assciated surge in transmissions [15] . In fact, there have already been reports of the current vaccines losing effectiveness against the emerging SARS-CoV-2 variants [21, 22] . Variants of high concern, B.1.1.7 (VOC-202012/01) identified in the U.K.; 501Y.V2 (B.1.351) identified in the South Africa; and P.1 (B.1.1.28.1) identified in Brazil, all have more than 20 mutations in the S protein, raising the concern that both natural infection-derived immunity as well as vaccine-derived immunity may fail to confer protection against a SARS-CoV-2 variant (re)infection [15, 21] . This hypothesis has been partly supported by a recent study showing that the 501Y.V2 (B.1.351) variant presents a complete immune escape from natural immunity [67] . Nonetheless, efforts to design new vaccines or to improve the existing vaccines that can confer immunity against SARS-CoV-2 and its variants are underway. In this regard, numerous in silico studies have also been conducted to design multi-epitope vaccine candidates against SARS-CoV-2. However, most of these studies focused on epitopes from SARS-CoV-2 S protein as the antigenic targets [37, [68] [69] [70] [71] , while others included the N and M proteins [72, 73] . The primary structure of the protein can provide information regarding mutations and the conservation of residues, while the secondary structure, depicting the spatial orientation, such as β-folds and α-helices, can give insights into how the protein will fold into its final form. Finally, the tertiary structure depicts the three-dimensional orientation of the protein in space [74] . The nature of the amino acids and protein folding can affect the physio-chemical properties of the protein, such as solubility, stability and interaction with other molecules, etc. [74] . The physico-chemical analyses revealed the molecular weight of our vaccine, 34.34 kDa, indicating that it is an ideal vaccine candidate, as proteins with a molecular weight less than 100 kDa are considered to be efficient vaccine candidates [70] . The instability index of our vaccine construct was shown to be 33.06, indicating that our construct is stable, as proteins with an index greater than 40 are generally considered unstable [38] . The vaccine construct was also found to be thermo-stable, with an aliphatic index of 72.00. A high aliphatic index is an indication of a protein to be thermo-stable over a wide temperature range [75] . The GRAVY index was predicted to be −0.250, reflecting high solubility of the vaccine and its effective interaction with water molecules [76] . The half-life of the vaccine was found to be 30 h in mammalian reticulocytes (in vitro), >20 h in yeast (in vivo), and >10 h in Escherichia coli (in vivo). Vaccine candidates with similar a half-life have been found to induce strong humoral and cellular responses in mice [77] . Finally, a refined and validated vaccine model showed that the designed vaccine is of high quality, with >97% of residues in the favored region [78] . In our study, we analyzed all the structural proteins, E, N, M, and S, not only in SARS-CoV-2 and its known variants, but also from six other human coronaviruses, and selected highly conserved epitopes with the aim of designing an effective multi-subunit vaccine that can be effective against SARS-CoV-2 and its variants and may also confer cross-immunoprotection against other HCoVs. Yazdani et al. [79] conducted a similar study using the E, N, M, and S proteins, however, this work investigated the conservation of peptides for only SARS-CoV-2 and SARS-CoV. In addition, the final epitopes selected in their study presented a lower antigenicity score (~0.7), while the epitopes we used in our design exhibited a high immunogenicity score of 1.5 for CD8+ T epitopes residing in S protein, and of 1.8 for linear B cell epitopes residing in the N protein. As our multi-subunit vaccine contains peptides from multiple proteins which are 100% conserved in the B.1.1.7, 501Y.V2, and P.1 variants SARS-CoV-2, mutations that potentially reduce immunogenicity are not likely to occur. The existence of cross-immune reactivity has been documented in several studies [80] [81] [82] [83] [84] . It has been reported that the antigenic domains of the S and N proteins are highly cross-reactive across coronaviruses [80, 82, 83] , and it has been shown that T cell immunity previously induced by circulating human alpha-and beta-HCoVs in young adults is protective against severe clinical outcomes [81] . An ideal vaccine candidate, on the one hand, should be immunogenic/antigenic, and should be able to elicit targeted humoral and cellular immune responses, while, on the other hand, should also be non-toxic and non-allergenic [85] , hence providing protection against adverse effects. Vaccines using a complete agent (virus, for example) or large proteins carry an increased chance of allergenic response due to the unnecessary antigenic load introduced in the vaccine design [86] . This drawback is avoided in our construct, as we used short immunogenic peptides that were predicted to be non-allergenic and non-toxic. Moreover, the predicted conserved T and B cell epitopes used for our multisubunit vaccine can elicit cellular and humoral immune responses simultaneously, further increasing the chance of triggering an efficient overall immune response. Additionally, our vaccine was constructed with β-defensin as an adjuvant at the N-terminal, which allows for a long-lasting immune response. β-defensin has been shown to be an efficient adjuvant capable of enhancing the immunogenicity of candidate vaccines against viruses, including MERS-CoV [14, 87] . In the immune simulation analyses, these properties were predicted for our vaccine construct as well, showing that the vaccine may be able to elicit a sustained IgG, IgM, Band cytotoxic T-response, and interferon-gamma response for well over 300 days. The activation/release of interferon-gamma is a hallmark of antiviral activity by T cells, as this cytokine is capable of killing the virus-infected cells [88] . Similarly, the analysis of the immunological response in individuals enrolled in phase 1/2 trial of the ChAdOx1 nCoV-19 vaccine showed that a single dose of the vaccine in adults was able to induce a potent immune response, primarily composed of neutralization antibodies of IgG1 and IgG3 subclasses, TH-1 response with interferon-γ, and tumor necrosis factor-α cytokine secretion by CD4+ T cells as well as mono-and poly-functional and cytotoxic CD8+ T cells up to 8 weeks after vaccination [89] . Similar results were observed for the Pfizer/BioNTech and Moderna vaccines, where robust T cell response against spike protein of SARS-CoV-2 was observed in vaccinated individuals, which was more profound than the response observed in convalescent patients [90] . Moreover, the CD8+ T cell epitopes in our vaccine design were linked with the AYY linker that provides binding sites for the transporter associated with antigen processing (TAP) to facilitate epitope presentation [91] , while B cell epitopes were linked with the glycine rich GPGPG linker, which gives the vaccine structural flexibility, prevents junctional epitope formation, and helps immune processing and epitope presentation [91] . Lastly, in molecular docking studies, using TLR4 showed that our vaccine construct, and especially several of β-defensin residues, a known TLR4 agonist [46, 47] used as an adjuvant in our vaccine, can form stable complexes with the TLR4. TLR-4 is capable of recognizing various exogenous and endogenous ligands such as viral protein, lipopolysaccharide, heat shock proteins, β-defensin, etc., which leads to dendritic cell maturation and up-regulation of costimulatory molecules, resulting in a type 1 polarized adaptive immune response [46, 47] . Recently, Zhao et al. showed that SARS-CoV-2 is able to activate TLR4 to activate an anti-bacterial like immune response [92] . Similarly, other in silico studies have shown that the ability of a vaccine construct to form stable complexes with these TLRs may lead to the induction of antiviral immunity [37, 93] . Multi-epitope vaccines designed using immune-informatics have been widely studied against viruses as well as against cancer [94] [95] [96] [97] . Some of these vaccines have shown promising results in vivo, while some entered clinical trials [98] [99] [100] [101] [102] . Analyzing a vaccine construct in silico before proceeding to in vitro evaluations and clinical trials is an efficient way to determine the vaccine's potential with confidence, economizing on time, effort, and money. Here we present in silico analyses of a multi-subunit vaccine design that, from multiple perspectives, shows potential to elicit a strong immune-protective response against SARS-CoV-2 and its emerging variants, while carrying minimal risk for causing adverse effects. Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/vaccines9070702/s1, Figure S1 : Multisubunit vaccine construct, Figure S2 : Phylogenetic analysis between human coronaviruses, Figure S3 : Conservation of envelope, membrane, nucleocapsid and spike protein sequences in seven HCoVs, Figure S4 : Conservation of envelope, membrane, nucleocapsid, and spike protein sequences in SARS-CoV-2, SARS-CoV and MERS-CoV, Figure S5 : Conservation of envelope, membrane, nucleocapsid, and spike protein sequence in 5000 SARS-CoV-2 sequences and in B.1.1.7, B.1.351, and P.1 variants of SARS-CoV-2, Table S1 : Consensus sequences used to predict the epitopes for vaccine design, Table S2 : List of all epitopes predicted with CTLpred. A pneumonia outbreak associated with a new coronavirus of probable bat origin Characteristics of and Important Lessons From the Coronavirus Disease 2019 (COVID-19) Outbreak in China: Summary of a Report of 72314 Cases From the Chinese Center for Disease Control and Prevention Hosts and Sources of Endemic Human Coronaviruses Middle East respiratory syndrome coronavirus (MERS-CoV): Announcement of the Coronavirus Study Group Coronaviruses post-SARS: Update on replication and pathogenesis A previously undescribed coronavirus associated with respiratory disease in humans Human Coronaviruses: A Review of Virus-Host Interactions Characterization and complete genome sequence of a novel coronavirus, coronavirus HKU1, from patients with pneumonia Characteristics of SARS-CoV-2 and COVID-19 SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor Functional assessment of cell entry and receptor usage for lineage B beta-coronaviruses, including 2019-nCoV The Novel Coronavirus: A Bird's Eye View The outbreak of SARS-CoV-2 pneumonia calls for viral vaccines Looking beyond COVID-19 vaccine phase 3 trials Vaccines and SARS-CoV-2 variants: The urgent need for a correlate of protection SARS-CoV-2 Vaccines Sputnik V COVID-19 vaccine candidate appears safe and effective Safety and efficacy of an rAd26 and rAd5 vector-based heterologous prime-boost COVID-19 vaccine: An interim analysis of a randomised controlled phase 3 trial in Russia The SARS-CoV-2 Spike Glycoprotein Biosynthesis, Structure, Function, and Antigenicity: Implications for the Design of Spike-Based Vaccine Immunogens The Impact of Mutations in SARS-CoV-2 Spike on Viral Infectivity and Antigenicity Efficacy of ChAdOx1 nCoV-19 (AZD1222) vaccine against SARS-CoV-2 variant of concern 202012/01 (B.1.1.7): An exploratory analysis of a randomised controlled trial Understanding the T cell immune response in SARS coronavirus infection Potential impact of individual exposure histories to endemic human coronaviruses on age-dependent severity of COVID-19 W-IQ-TREE: A fast online phylogenetic tool for maximum likelihood analysis MAFFT multiple sequence alignment software version 7: Improvements in performance and usability Deciphering key features in protein structures with the new ENDscript server Jalview Version 2-A multiple sequence alignment editor and analysis workbench BioEdit: A User-Friendly Biological Sequence Alignment Editor and Analysis Program for Windows 95/98/NT AliView: A fast and lightweight alignment viewer and editor for large datasets Prediction of CTL epitopes using QM, SVM and ANN techniques A hybrid approach for predicting promiscuous MHC class I restricted T cell epitopes Prediction of continuous B-cell epitopes in an antigen using recurrent neural network The Immune Epitope Database (IEDB): 2018 update A server for prediction of protective antigens, tumour antigens and subunit vaccines Open Source Drug Discovery, C.; Raghava, G.P. In silico approach for predicting toxicity of peptides and proteins Design of multi-epitope vaccine candidate against SARS-CoV-2: A in-silico study Exploring NS3/4A, NS5A and NS5B proteins to design conserved subunit multi-epitope vaccine against HCV utilizing immunoinformatics approaches AllerTOP v.2-A server for in silico prediction of allergens The PSIPRED protein structure prediction server The I-TASSER Suite: Protein structure and function prediction The Phyre2 web portal for protein modeling, prediction and analysis Protein structure refinement driven by side-chain repacking ProSA-web: Interactive web service for the recognition of errors in three-dimensional structures of proteins Protein-Sol: A web tool for predicting protein solubility from sequence Toll-like receptor 4-dependent activation of dendritic cells by beta-defensin 2 Innate immune pattern recognition: A cell biological perspective Structural Analyses of Human Toll-like Receptor 4 Polymorphisms D299G and T399I* H++: A server for estimating pKas and adding missing hydrogens to macromolecules PatchDock and SymmDock: Servers for rigid and symmetric docking Fast interaction refinement in molecular docking A web server for fast interaction refinement in molecular docking Multiple ligand-protein interaction diagrams for drug discovery Immune system simulation online Large-scale validation of methods for cytotoxic T-lymphocyte epitope prediction Prediction of proteasome cleavage motifs by neural networks The role of the proteasome in generating cytotoxic T-cell epitopes: Insights obtained from improved predictions of proteasomal cleavage SARS-CoV-2-specific CD8+ T cell responses in convalescent COVID-19 individuals Magnitude and Dynamics of the T-Cell Response to SARS-CoV-2 Infection at Both Individual and Population Levels T cell and antibody kinetics delineate SARS-CoV-2 peptides mediating long-term immune responses in COVID-19 convalescent individuals SARS-CoV-2 Proteome-Wide Analysis Revealed Significant Epitope Signatures in COVID-19 Patients. Front Comprehensive analysis of T cell immunodominance and immunoprevalence of SARS-CoV-2 epitopes in COVID-19 cases Unbiased Screens Show CD8(+) T Cells of COVID-19 Patients Recognize Shared Epitopes in SARS-CoV-2 that Largely Reside outside the Spike Protein Selective and cross-reactive SARS-CoV-2 T cell epitopes in unexposed humans The immunodominant and neutralization linear epitopes for SARS-CoV-2 SARS-CoV-2 501Y.V2 escapes neutralization by South African COVID-19 donor plasma Design of a multi-epitope vaccine against SARS-CoV-2 using immunoinformatics approach A candidate multi-epitope vaccine against SARS-CoV-2 Identification of putative vaccine candidates against Helicobacter pylori exploiting exoproteome and secretome: A reverse vaccinology based approach A rational design of a multi-epitope vaccine against SARS-CoV-2 which accounts for the glycan shield of the spike glycoprotein Immunoinformatics-guided designing of epitope-based subunit vaccines against the SARS Coronavirus-2 (SARS-CoV-2) Design of a peptide-based subunit vaccine against novel coronavirus SARS-CoV-2. Microb. Pathog. 2020, 145, 104236 Structural and functional analysis of hypothetical and conserved proteins of Clostridium tetani A simple method for displaying the hydropathic character of a protein Vaccination with a novel multi-epitope ROP8 DNA vaccine against acute Toxoplasma gondii infection induces strong B and T cell responses in mice Stereochemistry and Validation of Macromolecular Structures Design an Efficient Multi-Epitope Peptide Vaccine Candidate Against SARS-CoV-2: An in silico High prevalence of pre-existing serological cross-reactivity against severe acute respiratory syndrome coronavirus-2 (SARS-CoV-2) in sub-Saharan Africa Older adults lack SARS CoV-2 cross-reactive T lymphocytes directed to human coronaviruses OC43 and NL63 Epitope-resolved profiling of the SARS-CoV-2 antibody response identifies cross-reactivity with an endemic human CoV Analysis of Serologic Cross-Reactivity Between Common Human Coronaviruses and SARS-CoV-2 Using Coronavirus Antigen Microarray Antigenic crossreactivity between severe acute respiratory syndrome-associated coronavirus and human coronaviruses 229E and OC43 Understanding the immunogenicity and antigenicity of nanomaterials: Past, present and future Peptide Vaccine: Progress and Challenges The potent adjuvant effects of chicken beta-defensin-1 when genetically fused with infectious bursal disease virus VP2 gene Gamma-interferon exerts a critical early restriction on replication and dissemination of yellow fever virus vaccine strain 17D-204 T cell and antibody responses induced by a single dose of ChAdOx1 nCoV-19 (AZD1222) vaccine in a phase 1/2 clinical trial SARS -CoV-2 T-cell immunity to variants of concern following vaccination. bioRxiv 2021 Contriving Multi-Epitope Subunit of Vaccine for COVID-19: Immunoinformatics Approaches. Front SARS-CoV-2 spike protein interacts with and activates TLR41 Designing of a next generation multiepitope based vaccine (MEV) against SARS-COV-2: Immunoinformatics and in silico approaches Immunoinformatic identification of B cell and T cell epitopes in the SARS-CoV-2 proteome Bioinformatic prediction of potential T cell epitopes for SARS-Cov-2 Immunoinformatics and Structural Analysis for Identification of Immunodominant Epitopes in SARS-CoV-2 as Potential Vaccine Targets Immunoinformatic Analysis of T-and B-Cell Epitopes for SARS-CoV-2 Vaccine Design A phase I clinical trial of a multi-epitope polypeptide TAB9 combined with Montanide ISA 720 adjuvant in non-HIV-1 infected human volunteers A randomized phase II trial of multiepitope vaccination with melanoma peptides for cytotoxic T cells and helper T cells for patients with metastatic melanoma (E1602) Immunologic response to the survivin-derived multi-epitope vaccine EMD640744 in patients with advanced solid tumors Immunological features and efficacy of a multi-epitope vaccine CTB-UE against H. pylori in BALB/c mice model Rational design and efficacy of a multiepitope recombinant protein vaccine against foot-and-mouth disease virus serotype A in pigs