Skip to main content

ORIGINAL RESEARCH article

Front. Immunol., 13 March 2023
Sec. Vaccines and Molecular Therapeutics
This article is part of the Research Topic Methods in Vaccines and Molecular Therapeutics: 2022 View all 17 articles

Development of multi-epitope vaccines against the monkeypox virus based on envelope proteins using immunoinformatics approaches

Caixia Tan,&#x;Caixia Tan1,2†Fei Zhu,,,,&#x;Fei Zhu2,3,4,5,6†Pinhua Pan,,,,*Pinhua Pan2,3,4,5,6*Anhua Wu,*Anhua Wu1,2*Chunhui Li,*Chunhui Li1,2*
  • 1Department of Infection Control Center of Xiangya Hospital, Central South University, Changsha, Hunan, China
  • 2National Clinical Research Center for Geriatric Disorder, Xiangya Hospital, Changsha, Hunan, China
  • 3Department of Respiratory Medicine, National Key Clinical Specialty, Branch of National Clinical Research Center for Respiratory Disease, Xiangya Hospital, Central South University, Changsha, Hunan, China
  • 4Center of Respiratory Medicine, Xiangya Hospital, Central South University, Changsha, Hunan, China
  • 5Hunan Engineering Research Center for Intelligent Diagnosis and Treatment of Respiratory Disease, Changsha, Hunan, China
  • 6Clinical Research Center for Respiratory Diseases in Hunan Province, Changsha, Hunan, China

Background: Since May 2022, cases of monkeypox, a zoonotic disease caused by the monkeypox virus (MPXV), have been increasingly reported worldwide. There are, however, no proven therapies or vaccines available for monkeypox. In this study, several multi-epitope vaccines were designed against the MPXV using immunoinformatics approaches.

Methods: Three target proteins, A35R and B6R, enveloped virion (EV) form-derived antigens, and H3L, expressed on the mature virion (MV) form, were selected for epitope identification. The shortlisted epitopes were fused with appropriate adjuvants and linkers to vaccine candidates. The biophysical andbiochemical features of vaccine candidates were evaluated. The Molecular docking and molecular dynamics(MD) simulation were run to understand the binding mode and binding stability between the vaccines and Toll-like receptors (TLRs) and major histocompatibility complexes (MHCs). The immunogenicity of the designed vaccines was evaluated via immune simulation.

Results: Five vaccine constructs (MPXV-1-5) were formed. After the evaluation of various immunological and physicochemical parameters, MPXV-2 and MPXV-5 were selected for further analysis. The results of molecular docking showed that the MPXV-2 and MPXV-5 had a stronger affinity to TLRs (TLR2 and TLR4) and MHC (HLA-A*02:01 and HLA-DRB1*02:01) molecules, and the analyses of molecular dynamics (MD) simulation have further confirmed the strong binding stability of MPXV-2 and MPXV-5 with TLRs and MHC molecules. The results of the immune simulation indicated that both MPXV-2 and MPXV-5 could effectively induce robust protective immune responses in the human body.

Conclusion: The MPXV-2 and MPXV-5 have good efficacy against the MPXV in theory, but further studies are required to validate their safety and efficacy.

1 Introduction

Monkeypox is a zoonotic disease caused by an enveloped double-stranded DNA (dsDNA) virus known as monkeypox virus (MPXV), which belongs to the Orthopoxvirus genus of the Poxviridae family (1). Most cases of monkeypox present with mild disease symptoms, such as fever, rash, lymphadenopathy and intense asthenia, but certain population groups (children, pregnant women, and immunocompromised individuals) may suffer from severe disease or complications, including secondary bacterial infections, pneumonitis, respiratory distress, sepsis, encephalitis and even loss of vision following cornea infection (2). The MPXV was first identified in humans in 1970 (3). Since then, sporadic cases have been detected from time to time, but most of them were confined to West and Central African nations, or have a history of travel to African countries or exposure to imported animals (3). However, since May 2022, there has been a sharp increase in the number of MPXV infections worldwide, with the majority of confirmed cases occurring in countries that have not historically reported human monkeypox, and there are no obvious links between the documented cases in different countries (4). As of 13 November 2022, a total of 79,411 laboratory-confirmed monkeypox cases and 50 deaths have been reported in 109 countries, according to WHO (5). However, aside from several antiviral drugs for smallpox (Tecovirimat, cidofovir, and brincidofovir) that can be used to treat monkeypox under certain conditions, there is no proven post-exposure treatment (6). Hence, an effective vaccination against MPXV is warranted.

To date, there is no specialized vaccine against MPXV, but a small study based on previous African data estimated that the smallpox vaccine may provide approximately 85% cross-protection against MPXV (7). There are several available smallpox vaccines, including ACAM20, MVA-BN and LC16. ACAM2000, a live vaccinia virus-based preparation, and MVA-BN, a non-replicating smallpox vaccine, were both licensed by the US Food and Drug Administration (FDA) for active immunization against smallpox (8). However, the former may lead to secondary vaccinia infection in vaccine recipients and those in close contact with them, and some vaccine recipients may develop myocarditis or pericarditis (8). Although data from human immunogenicity trials and pre-clinical studies suggested that the later might offer some protection against MPXV in humans, no large-scale clinical study has been conducted to determine its actual level of protection against monkeypox infection in humans (8). Similar to MVA-BN, LC16 is a live, non-replicating attenuated vaccine, that was approved for the prevention of monkeypox in August 2022 in Japan (8). Also, its efficacy against monkeypox is extrapolated from data from animal studies rather than direct human trials. So, in short, there are research gaps in the monkeypox vaccine research that needs to be filled.

The multi-epitope vaccine, a novel type of vaccine developed in recent decades, has shown good immune effects against microbial infections and has fewer side effects than traditional vaccines in animals and early clinical trials (9). The application of immunoinformatics in vaccine development offers a more effective, cost-effective, and time-saving approach to developing vaccines for infectious diseases (10). The selection of targets is crucial for developing effective multi-epitope vaccines. Since no specific receptor for MPXV on the host cell membrane has been found so far, several envelope proteins, playing a key role in the invasion of host cells by MPXV, may be attractive targets for MPXV vaccine development (11, 12). A35R, homologous to vaccinia virus A33R, is an envelope glycoprotein of EV and is required for the formation of actin-containing microvilli and virion transmission between cells (13). B6R, located both on the membranes of infected cells and on the enveloped virion (EV) envelope, is similar to the complement control protein-like B5R of the vaccinia virus. Mice immunized with a subunit vaccine containing B5 could produce anti-B5 antibodies capable of neutralizing extracellular viruses, and destroying infected cells, according to Paran et al (14). Of note, existing neutralization tests demonstrate only the antibody reaction to B5 is available for neutralizing EV (14). Studies have shown that mice immunized with H3L can produce anti-H3L antibodies that protect them from a lethal dose of vaccinia virus (15). In addition, H3L can be recognized as an antigen by CD8+ T cells (16). Given the above, A35R, B6R, and H3L are attractive targets for MPXV vaccine development.

A recent study found that MPXV has an high mutation rate, which was 6 to 12 times higher than previously believed, and some of these mutations led to stronger MPXV infectivity (17).To provide broader protection against different strains of MPXV, three reference strains from three different epidemic phases of monkeypox were selected as the research objects, and the A35R, B6R and H3L of virus strains were selected as the protein targets to develop a promising multi-epitope vaccine using the immunoinformatics approaches.

2 Methodology

2.1 Data retrieval and sequence alignment

The amino acid sequences of A35R, B6R and H3L of three strains (Zaire-96-I-16, 2001, GenBank: GCA_000857045.1; MPXV-UK_P2, 2018, GenBank: MT903344.1; 2022 reference strain, GenBank: GCA_014621545.1) were obtained in FASTA format from National Center for Biotechnology Information (NCBI) (https://www.ncbi.nlm.nih.gov/) database. To identify the common and specific sites of the single protein from different strains, we performed alignment of the target protein sequences from three objective MPXV strains using the Clustal Omega program (https://www.ebi.ac.uk/Tools/msa/clustalo/) available on the EMBL-EBI server. The visualization of the results of sequence alignment performed by Jalview v 2.11 software (18).

2.2 Prediction and selection of epitope

2.2.1 Prediction and selection of cytotoxic T lymphocyte epitope

The NetCTL 1.2 server (https://services.healthtech.dtu.dk/service.php?NetCTL-1.2) (19) was utilized to predict the Cytotoxic T Lymphocyte (CTL) epitopes restricted to 12 major histocompatibility complex (MHC) class I supertypes (A1, A2, A3, A24, A26, B7, B8, B27, B39, B44, B58, and B62). This method combines the prediction of peptide MHC class I binding, proteasomal C terminal cleavage and Transporter associated with Antigen processing (TAP) transport efficiency (19), and only epitopes with a combined score ≥0.75 were chosen for antigenicity testing using the VaxiJen v2.0 server (http://www.ddg-pharmfac.net/vaxijen/VaxiJen/VaxiJen.html) (20). VaxiJen allows antigen classification entirely based on protein physicochemical properties, overcoming the limitations of alignment-dependent methods. The accuracy of the VaxiJen in predicting the protective antigen of the virus is 70% (20). To minimize the toxicity and allergenicity of the vaccine and improve the safety of the vaccine, ToxinPred2 (https://webs.iiitd.edu.in/raghava/toxinpred2/algo.html), with an accuracy of 95.54% (21), and AllerTOP v. 2.0 (https://www.ddg-pharmfac.net/AllerTOP/index.html) (22), with an accuracy of 85.3%, were used to predict the toxicity and allergenicity of the epitopes, respectively.

The strength of immunogenicity reflects the ability of the antigen to induce an immune response. Class-I Immunogenicity tool hosted in IEDB (http://tools.iedb.org/immunogenicity/) (23) is considered to be able to predict the immunogenicity of peptides to some extent. So, the antigenic epitopes were further submitted to the Class-I Immunogenicity tool to check their immunogenicity. T cells require the detection of epitopes on MHC molecules to be activated and to exert their effector activity. The Tepitool tool (http://tools.iedb.org/tepitool/) (24) provided on Immune Epitope Database (IEDB) helps to predict the affinity between the peptide and class I human leukocyte antigen (HLA) alleles, and the IC50 is an indicator to evaluate the binding affinity of the HLA molecule to a peptide. Hence, The CTL epitopes with antigenicity, immunogenicity, non-allergenicity and non-toxicity were submitted to the Tepitool tool, and the default threshold IC50 (affinity) ≤ 500nM was used as the screening criteria.

2.2.2 Prediction and selection of helper T lymphocyte epitope

The NetMHCIIpan-4.1 server (https://services.healthtech.dtu.dk/service.php?NetMHCpan-4.1) (25), appearing stronger predictive power for the prediction of HTL epitopes and class II MHC molecules compared to other servers, was employed to predict the (Helper T Lymphocyte) HTL epitopes that possess a good affinity to 14 HLA class II alleles (HLA-DRB1*01:01, HLA-DRB1*03:01, HLA-DRB1*04:01, HLA-DRB1*04:04, HLA-DRB1*04:05, HLA-DRB1*07:01, HLA-DRB3*02:02, HLA-DRB1*09:01, HLA-DRB1*11:01, HLA-DRB1*13:02, HLA-DRB1*15:01, HLA-DRB3*01:01, HLA-DRB4*01:01, HLA-DRB5*01:01) covering 95% of the world’s the population. For each HTL epitope, the IC50 and percentile rank were determined, and the epitopes with an IC50 ≤500 and a percentile rank ≤ 5% were selected for antigenicity, toxicity and allergenicity prediction via VaxiJen v2.0, ToxinPred2 and AllerTOP v. 2.0 servers, respectively. Interferon-gamma (IFN-γ) plays an essential role in viral clearance and the activation of the host immune response (26). In addition to regulating B cell growth and immunoglobulin secretion, interleukin 4(IL-4) also influences T cell differentiation (27). Interleukin 2(IL-2) is also a multifunctional cytokine that helps T cells grow, proliferate, and differentiate and potently enhances the function of B-cells (28). Hence, the selected HTL epitopes were further submitted to IFNepitope (http://crdd.osdd.net/raghava/ifnepitope/) (29), IL4Pred (http://crdd.osdd.net/raghava/il4pred/) (30) and IL2Pred (https://webs.iiitd.edu.in/raghava/il2pred/) (31) servers to assess their potential to induce IFN-γ, IL-4, and IL-2 secretion, and only epitopes that induce IFN-γ production and at least IL-2 or IL-4 production would be ultimately considered for vaccine construction.

2.2.3 Prediction and selection of B-cell epitope

The ABCpred (https://webs.iiitd.edu.in/raghava/abcpred/) was applied to predict linear B-cell epitopes based on recurrent neural network (RNN). The 16 amino acids were selected as the window length (32). As with T-cell epitopes, the antigenicity was checked by VaxiJen v2.0, and only epitopes with an antigenicity over 0.8 were further submitted to AllerTOP v. 2.0 and ToxinPred2 servers to check their allergenicity and toxicity, respectively.

2.2.4 Population coverage prediction of T cell epitopes

Population coverage in different regions is one of the important indicators in measuring the effectiveness of the vaccine. And its level size is mainly affected by the distribution and expression frequency of the HLA alleles that bind to the multi-epitope vaccine. To increase the population coverage of limited length multi-epitope vaccine, we calculated the world population coverage for each selected CTL and HTL epitope using the IEDB’s population coverage tool (http://tools.iedb.org/population/) (33), and only epitopes with a world population coverage greater than 5% were subjected to further analysis.

2.2.5 Homology analysis of selected epitopes with human proteomes

The Peptide match server can quickly and accurately match the submitted peptides with the UniProtKB Human complete proteome (34). To avoid inducing host autoimmune diseases and cross-reactions, each selected epitope was further submitted to the Peptide match server (https://research.bioinformatics.udel.edu/peptidematch/index.jsp) (34) to check their homology to the human proteome, and only the epitopes are non-homologous to the human proteome were selected for vaccine construction.

2.3 Molecular docking of T-cell epitopes with major histocompatibility complex molecules

MHC molecules are essential for presenting the pathogen-derived peptides to the respective T cells. Thus, it is necessary to perform molecular docking between T-cell epitopes and MHC molecules to understand the possible binding mode between them. First, the PEP-FOLD 3.5 server (https://bioserv.rpbs.univ-paris-diderot.fr/services/PEP-FOLD3/) (35)was used to model the three-dimensional (3D) structure of T-cell epitopes as the ligand, and RCSB Protein Data Bank (RCSB PDB) website (https://www2.rcsb.org/search/advanced/sequence) (36) was employed to acquire the 3D structure of MHC molecules as for the receptor. For some class I and class II MHC molecules for which the 3D structures are unavailable on the RCSB PDB website, we use the 3D structures of HLA-A*02:01 and HLA-DRB1*01:01 as the substitutes, respectively. Then, the energy minimization of receptors and ligands was carried out by Swiss-PdbViewer PDB viewer software (37), and following that, the molecular docking was run by Autodock Vina software hosted in PyRx v0.8 (38). AutoDock Vina, a classic molecular docking tool, has better performance for the prediction of docking between small molecules and proteins (39). The docking results were further analyzed using LigPlot+ v.2.2 software (40). Furthermore, to comparatively evaluate of binding affinity between the selected epitopes and MHC molecules, the molecular dockings between seven experimental verified antigen epitopes of poxvirus and HLA-A*02:01 molecules were selected as the positive controls.

2.4 Design of multi-epitope vaccine constructs

All T-cell and B-cell epitopes that satisfied the filtration criteria were incorporated into the construction of the multi-epitope vaccine. The linker is an essential part of recombinant fusion proteins, which is capable of forming effective epitope conjugation between epitopes and allowing the independent immunological activities of the epitopes (41). Low immunogenicity is a major deficiency of multi-epitope vaccines. According to the previous study, a suitable adjuvant molecule can significantly boost the peptide vaccine’s immunogenicity and longevity (42). 50s ribosomal protein L7/L12 from Mycobacterium tuberculosis is capable of promoting DC maturation and diverse pro-inflammatory cytokine production via interacting with TLR4, which helps enhance cellular immunity (43). In addition to direct antiviral activity, β-defensins can regulate adaptive immunity to viral infection by recruiting naive T cells and immature dendritic cells (DCs) to the infected site (44). LT-IIC, a type II heat-intolerant enterotoxin (HLT) produced by Escherichia coli (E.coli), is an attractive adjuvant with the unique property of enhancing the body’s humoral and cellular immune responses to co-administered antigens (45). The Cholera Toxin B subunit (CTB) is also considered an attractive adjuvant due to its ability to efficiently bind to antigen-presenting cells and enhance the stability and half-life of the entire vaccine molecule (46). RS09, an LPS peptide mimic, is an agonist of TLR4. When used as an adjuvant in certain vaccines, it is capable of resulting in stronger immune activation and higher antibody production (47). All five adjuvants have been widely used in the design of multi-epitope vaccines in silico (4851) and have been shown to stimulate the host’s immune response in animal experiments (43, 45, 5254). So, in this study, each of the five adjuvants was added to the N-terminus of the vaccine using EAAAK liners, resulting in the formation of five vaccine constructs. Afterward, the PADRE sequence was linked to the adjuvant with the help of the EAAAK linker. When added to the vaccine construct, the PADRE, which can bind to most common HLA-DR molecules with high affinity, has been shown to optimize antibody responses and deliver significant helper T-cell activity in vivo (55). The CTL epitopes were connected by AAY linkers, the HTL epitopes were linked using GPGPG linkers, and the B-cell epitopes were spaced with the help of KK linkers. The last CTL epitope was linked to the first HTL epitope and the first HTL epitope was connected to the first B-cell-epitopes using the HEYGAEALERAG linker. Finally, the “6xHis tag” was appended to the C-terminal of the vaccine construct via the RVRR linker to facilitate protein purification but not affect the functionality of the fusion protein.

2.5 Prediction of biophysical and biochemical features of the multi-epitope vaccine constructs

2.5.1 Prediction of antigenicity, allergenicity, toxicity and physicochemical properties of the vaccine

To determine the ability of the vaccine constructs to bind to immune cell receptors (antigenicity), both the VaxiJen and the ANTIGENpro (http://scratch.proteomics.ics.uci.edu/) servers were employed. The allergenicity of the vaccine was assessed using AllerTOP v. 2.0 server. The vaccine’s toxicity was predicted by the ToxinPred2 server. The physicochemical properties of vaccine constructs, including molecular weight, the number of amino acids, theoretical isoelectric point (pI), estimated half-life, instability index, aliphatic index, and grand average of hydropathicity (GRAVY) index, were predicted using ProtParam tool (https://web.expasy.org/protparam/) (56).

2.5.2 The prediction of transmembrane domains, probability of solubility, and signal peptide sequence of the vaccine

Since most amino acids constituting transmembrane region proteins are hydrophobic amino acids and membrane proteins cannot be expressed in prokaryotic expression systems, it is necessary to predict the transmembrane sequence of vaccines. The DeepTMHMM (https://dtu.biolib.com/DeepTMHMM) server, which is the most complete and best-performing method for the prediction of the topology of both alpha-helical and beta-barrel transmembrane proteins (57), was employed to predict the transmembrane domains of the vaccine. Many experimental studies are hampered by protein insolubility. The SOLpro server (http://scratch.proteomics.ics.uci.edu/) (58) can accurately predict the probability of solubility in E. coli of the protein upon overexpression. So, all designed vaccines were submitted to the SOLpro server to predict their solubility. To avoid the amino acid sequence of the vaccine from being cut off by signal peptidase as a signal peptide, the signalP-6.0 online server(https://services.healthtech.dtu.dk/service.php?SignalP-6.0) (59) was utilized to predict whether a vaccine has signal peptides, as well as the location of their cleavage sites.

2.6 Secondary structure prediction, tertiary structure modeling, refinement, and validation

Out of all vaccine constructs, only antigenic, non-allergenic, non-toxic, stable, hydrophilic, thermostable, without transmembrane domains and signal peptide sequence, and highly soluble constructs were further submitted to the PSIPRED 4.0 server (http://bioinf.cs.ucl.ac.uk/psipred/) to obtain its secondary structure. The PSIPRED 4.0 server performed protein secondary structure prediction based on position-specific scoring matrices with a prediction accuracy of over 84% (60). The tertiary structures of the vaccine constructs were modeled using the Robetta server (http://robetta.bakerlab.org). Based on the fact that a confident template for the putative domains of a submitted sequence is available or unavailable, this server performs comparative modeling or de novo modeling, respectively (61). And this server will generate five models of each submitted sequence. Then, the optimal crude 3D models of each vaccine sequence were submitted to the GalaxyRefine module placed in the GalaxyWEB server(https://galaxy.seoklab.org/), which can rebuild unreliable loops or termini of the initial model structures using an optimization-based refinement method and further generate five refined models of each crude model (62). Referring to the criterion that the better model with a lower MolProbity score, we picked the best 3D refined model for each vaccine construct and submitted it to the ERRAT tool, PROCHECK tool of structure validation service SAVES v6.0 (https://saves.mbi.ucla.edu/) and ProSA-Web server (https://prosa.services.came.sbg.ac.at/prosa.php) for quality verification. The overall quality factor generated by the ERRAT program is the number of non-bonded interactions (side chains) formed between pairs of different atomic types within 0.35nm, and models with an overall quality factor value greater than 85 are generally considered to have good quality (63). The Ramachandran plot generated by the PROCHECK tool checks the stereochemical quality of the protein structure, and a better model has more residues in the Ramachandran-favored region and fewer residues in the Ramachandran-disallowed region. The Z-score predicted by the ProSA-Web server reflects the overall quality of the model, and a positive Z-score indicates the input structure has problematic or erroneous parts (64).

2.7 Molecular docking

Stable interaction between the vaccine candidate and immune receptor is an important precondition to generating successful immune reactions. The TLRs not only hold a vital position in the first line of defense against pathogens but also play a vital role in linking innate immunity with adaptive immunity (65). Previous studies revealed that TLR2 and TLR4 can induce antiviral immune responses by detecting the virus coat proteins (66, 67). Class I and class II MHC molecules help the CD4+ T and CD8+ T cells recognize foreign antigens(vaccine). Therefore, the molecular docking of the vaccine candidates with TLRs (TLR2, TLR4) and MHCs (HLA-A*02:01, HLA-DRB1*01:01) was performed using the ClusPro 2.0 server (https://cluspro.bu.edu/login.php) (68). Since 2004, ClusPro has consistently been one of the best docking servers as demonstrated in Critical Assessment of Predicted Interactions (CAPRI), providing high predictive performance for the docking of protein-protein complexes (68). In each docking case, a total of thirty models were generated, and ten models with highly populated clusters were presented. Then, the best model in each docking case was further submitted to the HADDOCK 2.4 server (https://wenmr.science.uu.nl/haddock2.4/) (69), which is capable of refining the complex structures by rearranging the side-chain in the restricted interface and optimizing soft rigid-body, for refinement. After refinement, the docked models were submitted to the PRODIGY server (https://bianca.science.uu.nl/prodigy/) (70) for predicting binding energy. Furthermore, the PDBsum server (http://www.ebi.ac.uk/thornton-srv/databases/cgi-bin/pdbsum/GetPage.pl?pdbcode=index.html) (71) was utilized to investigate interacting residues between docked chains (i.e., the vaccines, the TLRs and MHCs). The PyMOL v2.5 was used to visualize the vaccine-receptor complex (72), and the Blender v3.3 software was employed to render the final image (73).

2.8 Molecular dynamics simulation and MM-PBSA analysis

2.8.1 Molecular dynamics simulation

The Molecular dynamics (MD) simulations of several complexes with better docking results were conducted using the Gromacs v2022.1 software (74), which is helpful in better understanding the dynamics and structural stability of the docking complexes. In the preprocessing stage, the Amber14SB _parmbsc1 force field parameters were added to the complex system to obtain the coordinate file and topology file of the vaccine-receptor complexes (75). Subsequently, each system was solvated using the transferable intermolecular potential 3P (TIP3P) water model, following being neutralized with Cl- ions, leaving the topological and structural coordinates in a steady state. Following this preprocessing, energy minimization was performed for the entire protein complex system; Then, under the NVT ensemble, each system was gradually heated from 0K to 310 K, and the modified Berendsen thermostat was used to equilibrium temperature (76). Then, in the NPT ensemble, the pressure of the system was adjusted to 1atm with regulating of Parinello-Rahman barostat (77). Finally, a 50ns MD simulation was run for all docking complex systems, adopting constant temperature and pressure. The LINCS algorithm (78) and the PME (Particle mesh Ewald) method (79) were used to constrain the hydrogen bond and calculate the long-range electrostatic energy, respectively. After 50ns MD simulation, the post-MD simulation analyses of the MD trajectories were performed, including the root mean square deviation (RMSD) of each system, each system chain’s root mean square fluctuation (RMSF), the radius of gyration (Rg) and hydrogen bonds of each system.

Then, the trajectories of the last 20ns were extracted for further analysis. In addition to classic MD analysis, the conformational free energy of each complex was studied at the transient stages via studying free energy landscapes (FELs) (80). The PCA was run using gmx covar and gmx anaeig programs, and the PC1 and PC2 were further used to calculate and analyze Gibb’s free energy landscape (Gibb’s FEL) with the gmx sham program (80). Gibb’s FEL plots were generated using the DuIvyTools (DOI:10.5281/zenodo.6340263). In Gibb’s FEL plots, the stable energy states were shown blue. In addition, the correlations between the position fluctuations of each residue in the vaccine chain and the position fluctuations of each residue in the receptor chain were analyzed with the dynamical cross-correlation (DCC) analysis by the Bio3D package in the R v4.22 software.

2.8.2 MM-PBSA calculations

The Gmx_MMPBSA v1.52 tool was applied to calculate the binding free energy between the vaccine and the receptors using Poisson–Boltzmann and surface area continuum solvation (MM-PBSA) method (81), and the calculation formula is as follows:

ΔGbind,slov=ΔGbind,vaccum+ΔGslov,complex(ΔGsolv, ligand+ΔGsolv,receptor)

2.9 Population coverage

More than a thousand distinct HLA alleles have been identified in humans, and the frequency of distribution of each HLA allele varies across regions and ethnicities. Given that T lymphocytes can only recognize molecular complexes composed of epitopes bound to MHC molecules, a particular epitope triggers immune responses only in individuals expressing HLA alleles with an affinity for this epitope. The IEDB population coverage tool (http://tools.iedb.org/population/) was employed to calculate the population coverage of the designed vaccine in this work based on the data of TCR specificity, HLA restriction, and HLA allele frequencies.

2.10 Immune simulation

For assessing the potential immune efficacy of the designed vaccine, we submitted it to the C-IMMSIM v10.1 tool (https://kraken.iac.rm.cnr.it/C-IMMSIM/) (82). The C-IMMSIM server recognizes epitopes by Position Specific Scoring Matrix (PSSM) methods, and can simultaneously stimulate the immune response of three compartments that represent three separate anatomical regions found in mammals (the bone marrow, the thymus and lymph nodes) to antigen epitopes. During the immune simulation, all parameters except for time steps were left at their default settings. In this study, the simulation time steps were set to 1050 (1-time step equals 8 h), and each injection point was set at time steps 1, 84, and 168, respectively. The default dose of each antigen injection was 1,000, and the time interval between each dose was four weeks, which is the recommended interval between injections for most commercial vaccines. To comparatively evaluate the candidate vaccine’s immunogenicity, a novel multi-epitope vaccine against Echinococcus granulosus, whose antigenicity and immunogenicity had been verified in vitro and in vivo experiments (83), was chosen as a positive control.

2.11 Codon optimization and in silico cloning

Due to the degeneracy of codons, each amino acid may correspond to multiple codons. However, there are great differences in codon preference among different species and organisms. To acquire the codon that is capable of efficiently encoding the targeted amino acid in the selected expression host, the codon optimization of the vaccine was conducted using the JCat server (http://www.jcat.de/CAICalculation.jsp) (84). The JCat server is an easy-to-use program, which can perform real-time calculations, eliminate the need for the definition of highly expressed genes, and avoid cleavage sites for specific restriction enzymes (84). The Escherichia coli (E. coli) K12 strain was selected as the expression host. Two parameters, codon adaptation index (CAI), and the percentage of GC content, were used to evaluate the expression efficiency of the vaccine in E. coli. The CAI, with a value range of 0 to 1, refers to the coincidence degree between synonymous codons and optimal codons, and a higher CAI indicates a higher expression level of the foreign gene in the host. The percentage of GC content represents the guanine (G)+cytosine (C) ratio in the DNA sequence, and its optimal range is 30%–70%. Otherwise, undesirable gene expression levels and transcriptional efficiency will result. Afterward, the BamHI and XhoI restriction endonuclease sites were integrated into the N- and C-terminal sites of optimized codon sequences, respectively. The entire sequence was cloned into the pET28a (+) vector using SnapGene software. Finally, the MPXV-5 vaccine protein was synthesized by Sangon Biotech co ltd (Shanghai, China). The flowchart of vaccine design is illustrated in Figure 1.

FIGURE 1
www.frontiersin.org

Figure 1 Flowchart of the multi-epitope vaccine design.

3 Result

3.1 Protein sequences retrieval and alignment

The NCBI accession numbers and amino acid sequences of A35R, B6R and H3L derived from three objective MPXV strains were presented in Supplementary Table 1. Protein sequence alignment results (Figure 2) demonstrated that there was no difference in three target protein sequences (A35R, B6R and H3L) between MPXV-UK_P2 and 2022 reference strain. Compared with the Zaire-96-I-16 MPXV strain, there were 3 mutation sites in A35R and 2 mutation sites in H3L, but no mutation sites in B6R. As shown in Figure 2, there were two mutation sites falling into the epitope regions.

FIGURE 2
www.frontiersin.org

Figure 2 The sequence alignment of three target proteins in three selected reference strains. (A) The sequences of A35R; (B) The sequences of B6R; (C) The sequences of H3L; The Zaire-96-I-16 strain is selected as the reference sequence in the sequence alignment, i.e., the first-row sequences were the sequences of Zaire-96-I-16 strain, the second-row sequences were the sequences of MPXV-UK_P2 strain, and the third-row sequences were the sequences of 2022 reference strain. The red arrows refer to mutation sites; The CTL, HTL, and B cell epitopes were highlighted with red, blue and purple box lines in the sequence, respectively.

3.2 Epitopes selection

The immune system is comprised of many specialized cell types, all of which work together to keep people healthy. CTLs, the main force of cellular immunity, are immune cells that can directly mediate the death of target cells and eliminate circulating antigens in the host body (85). Stimulated B cells can differentiate into antibody-producing plasma cells that play a vital role in the generation of humoral immune responses (86). HTLs, another group of immune cells, can assist both humoral and cellular immunity by secreting different cytokines (87). Hence, to construct an effective multi-epitope vaccine, both T-cell epitopes (CTL and HTL epitopes) and B-cell epitopes should be incorporated. After the layer-by-layer screening (Figure 3), a total of 7 CTL epitopes (Table 1), 6 HTL epitopes (Table 2), and 7 liner B-cell epitopes (Table 3) were finally included in the vaccine constructs. The binding affinity (IC50) of the experimentally determined poxvirus-specific T cell epitopes to HLA molecules was shown in Supplementary Table 2.

FIGURE 3
www.frontiersin.org

Figure 3 Flowchart of the epitope screening.

TABLE 1
www.frontiersin.org

Table 1 List of the CTL epitopes selected for vaccine construction.

TABLE 2
www.frontiersin.org

Table 2 List of the HTL epitopes selected for vaccine construction.

TABLE 3
www.frontiersin.org

Table 3 List of the B-cell epitopes selected for vaccine construction.

3.3 Molecular docking of the T-cell epitopes with MHC molecules

Except for only one epitope (STETSFNDK), all other epitopes had multiple binding alleles. Some had even as high as 8 (TMSAFLIVR) or 13 (FTYTGGYDV) alleles (Table 1). In this study, we performed molecular docking between each T-cell epitope and one of their corresponding alleles for a representative docking analysis (Supplementary Table 3). As shown in Supplementary Table 3, Autodock Vina docking data indicated that the binding energy between CTL epitopes and MHC class I alleles ranged from -10.2 to -5.9 Kcal/mol, and the binding energy range between HTL epitopes and MHC class II alleles was -7 to -5.8 Kcal/mol, the binding energy range between seven experimental verified antigen epitopes of poxvirus and HLA-A*02:01 molecules was -9 to -5.8 Kcal/mol, indicating the T-cell epitopes had a good affinity for MHC molecules. The predicted results of Ligplus v2.25 demonstrated that many hydrogen bonds were formed between T-cell epitopes and MHC molecules, and salt bridges were also detected in certain docking cases (Supplementary Table 3), which once again proved that T-cell epitopes had a strong affinity to MHC molecules.

3.4 Multi-epitope vaccine features

3.4.1 Prediction of antigenicity, allergenicity, and toxicity of the vaccine construct

The final multi-epitope vaccine construct was built with the following six parts: CTL epitopes, HTL epitopes, B-cell epitopes, PADRE, adjuvant and His-tag sequence. Five different adjuvants (50S ribosomal protein L7/L12, β-defensin, LT-IIC, CTB, RS09) resulted in five different vaccine constructs (Supplementary Table 4) named MPXV-1, MPXV-2, MPXV-3, MPXV-4 and MPXV-5, respectively. Except for the adjuvant, all other components of these constructs are identical. The antigenicity ranges of five vaccine constructs predicted by VaxiJen and ANTIGENpro servers were 0.62~0.69 and 0.56~0.84, respectively, which demonstrated that the five vaccine constructs can bind or interact with the immune effector cells or soluble antibodies (Table 4). Moreover, none of the five vaccine constructs were allergenic. However, the ToxinPred2 predicted results showed that MPXV-3 and MPXV-4 were toxic. So, the vaccine bodies of MPXV-3-4 and their corresponding adjuvants were separately submitted to the ToxinPred2 server, and its results revealed that the vaccine bodies were not toxic, while both adjuvants were toxic.

TABLE 4
www.frontiersin.org

Table 4 Properties of the vaccine constructs.

3.4.2 Physicochemical properties analysis of the vaccine construct

The physicochemical properties of the vaccine constructs are presented in Table 4. The molecular weight of five vaccine constructs (MPXV-1–5) was below 110 kDa, which was considered favorable for vaccine development. The theoretical pI for all vaccine constructs ranged from 8.61 to 9.56. The range of instability index (II) values for all constructs was between 34.19 and 36.62, suggesting that the constructs were stable in a test tube. In general, proteins with an instability index of less than 40 were considered to be stable. The aliphatic index is regarded as a positive factor for the increase in thermostability of globular proteins. In all developed constructs, the aliphatic index ranged from 82.25 to 89.06, indicating that these constructs were thermostable. All designed vaccines have an estimated half-life of 30 hours in mammalian reticulocytes, >20 hours in yeast, and >10 hours in E. coli. The range of the GRAVY index of MPXV-1-5 was -0.19 to -0.051. The GRAVY index is correlated with protein solubility, with negative GRAVY indicating a hydrophilic molecule.

Predictions from the DeepTMHMM server showed that none of the constructed vaccines, except MPXV-1, has a transmembrane helix (Supplementary Figure 1). This may be due to the adjuvant of MPXV-1 is located inside the membrane while the vaccine body is located outside the membrane, resulting in the formation of transmembrane sequences at the bridging site between the two parts. The probabilities of being soluble upon overexpression in E. coli of five constructs were 0.975118 (MPXV-1), 0.835770 (MPXV-2), 0.610239 (MPXV-3), and 0.867569 (MPXV-4), 0.869869 (MPXV-5), respectively, which demonstrated that all vaccine constructs had good solubility when heterologous expression in E. coli. Moreover, the predicted results of the signalP-6.0 server revealed that none of the vaccines had signal peptide sequences (Supplementary Figure 2).

Based on all the above predictions about the properties of all vaccine constructs, MPXV-2 and MPXV-5 were regarded as the most ideal vaccine constructs for immunological applications. Therefore, MPXV-2 and MPXV-5 were selected for further analysis.

3.5 Secondary structure prediction of the vaccine construct

The secondary structures of the final selected vaccine constructs were predicted by the PSIPRED server and predicted results (Supplementary Figure 3) showed that the MPXV-2 contained 58.22% α-helix, 31.22% coils and 10.56% β-strands. Compared to MPXV-2, the MPXV-5 secondary structure comprised a greater proportion of α-helix (61.10%), a smaller proportion of β-strands (6.79%), with the similar proportion of coils (32.11%) (Supplementary Figure 3).

3.6 Three-dimensional structure prediction, refinement, and validation of the vaccine construct

The 3D structure of the selected vaccine constructs was predicted by the Robetta server and then refined by the GalaxyWEB server. Based on the model quality of the refined models (Supplementary Table 5), the refined model 5 (Figure 4A) for MPXV-2 and the refined model 4 (Figure 4A) for MPXV-5 were chosen for further quality validation. The corresponding Overall Quality Factors generated by the ERRAT program were 91.8465 (MPXV-2-5) and 89.2105 (MPXV-5-4) (Supplementary Figure 4), and the corresponding Z-scores predicted by the ProSA-Web server were -6.94 (MPXV-2-5) and -6.36 (MPXV-5-4) (Figure 4B). The Ramachandran plot analysis of the refined model MPXV-2-5 revealed that 92.3% of residues located in the most favorable region, 6.4% in the allowed region, and only 1.3% in the disallowed region (Figure 4C). For the refined model MPXV-5-4, 92.4% of residues suited in the most favored region, 6.4% in the allowed region, and only 1.2% of residues in the disallowed regions (Figure 4C). The above analysis confirmed the excellent quality of the refined 3D models of the MPXV-2 and MPXV-5 vaccine constructs.

FIGURE 4
www.frontiersin.org

Figure 4 Refinement and quality validation of the three dimensional (3D) structure. (A) The refined 3D structure. (C) Ramachandran plots of the refined 3D structure generated by PROCHECK server; Red regions represent the most favored regions, dark yellow regions represent the additional allowed region, light yellow regions represent the generally allowed regions, and white regions denote the disallowed regions. (B) Z-Score plot of the refined 3D structure predicted by ProSA-Web server.

3.7 Molecular docking

The ClusPro generated thirty clusters in each docking case, and the best one with the largest cluster size in each docking case was selected to present (Table 5). Table 5 displayed the detailed docking results for all refined docking complexes. Compared with the MPXV-2-immune receptor complexes, all corresponding MPXV-5-immune receptor complexes had lower docking scores and binding energy, indicating that MPXV-5 had stronger affinity for all immune receptors (TLR2, TLR4, HLA-A*02:01 and HLA-DRB1*01:01). The predicted results of the PDBsum server (Figure 5; Supplementary Figure 5 and Table 5) showed that 6 hydrogen bonds, 11 hydrogen bonds and 2 salt bridges, 12 hydrogen bonds and 1 salt bridges, 16 hydrogen bonds and 6 salt bridges were formed in MPXV-2-TLR2, MPXV-2-TLR4, MPXV-2-HLA-A*02:01 and MPXV-2-HLA-DRB1*02:01 complexes, respectively. There are 12 hydrogen bonds with 3 salt bridges, 18 hydrogen bonds with 3 salt bridges, 12 hydrogen bonds with 3 salt bridges, and 12 hydrogen bonds with 5 salt bridges in MPXV-5-TLR2, MPXV-5-TLR4 and MPXV-5-HLA-A*02:01 and MPXV-HLA-DRB1*02:01 complexes, respectively (Figure 5; Supplementary Figure 5 and Table 5). Docking studies showed that there were significant interactions between the designed vaccines (MPXV-2 and MPXV-5) and the immune receptors.

TABLE 5
www.frontiersin.org

Table 5 Docking analysis of vaccine-receptors complexes.

FIGURE 5
www.frontiersin.org

Figure 5 Vaccine candidates (MPXV-2 and MPXV-5) docked with TLR2 and TLR4. (A) The binding mode between the MPXV-2 and TLR2 and the interacting residues of the MPXV-2 with TLR2; (B)The binding mode between the MPXV-2 and TLR4 and the interacting residues of the MPXV-2 with TLR4; (C) The binding mode between the MPXV-5 and TLR2 and the interacting residues of the MPXV-5 with TLR2; (D) The binding mode between the MPXV-5 and TLR4 and the interacting residues of the MPXV-5 with TLR4.

3.8 Molecular dynamics simulation

3.8.1 Root mean square deviation, root mean square fluctuation radius of gyration and hydrogen bonds analyses

RMSD can be used to measure changes in protein structure during MD simulations compared to the starting point. A leveling off of the RMSD curve can also indicate that the protein structure has equilibrated. The RMSD value of all MPXV-2-receptors complexes (MPXV-2-TLR2, MPXV-2-TLR4, MPXV-2-HLA-A*01:01) increased gradually in the range of 1-20ns and then maintained the level of about 1.3nm, 0.9nm, 1.0nm, respectively (Figure 6A). The average RMSD value of the MPXV-5-TLR2 was 0.82nm (Figure 6B), and several small fluctuations in RMSD value were observed. The average RMSD value of MPXV-5-TLR4 complex was the highest among all complexes, which was 1.51nm (Figure 6B). The RMSD value of the MPXV-5-HLA-A*01:01 complex stabilized at about 20 ns of the MD simulation with an average value of 1.05 nm (Figure 6B).

FIGURE 6
www.frontiersin.org

Figure 6 The plot of Root means square deviation (RMSD), Root means square fluctuation (RMSF), Radius of Gyration (Rg) and hydrogen bonds of vaccine-receptor complexes during the Molecular dynamics (MD) simulation. (A) RMSD of MPXV-2-TLR2, MPXV-2-TLR4 and MPXV-2-HLA-A*02:01 complexes; (B) RMSD of MPXV-5-TLR2, MPXV-5-TLR4 and MPXV-5-HLA-A*02:01 complexes; (C) Rg of MPXV-2-TLR2, MPXV-2-TLR4 and MPXV-2-HLA-A*02:01 complexes; (D) Rg of MPXV-5-TLR2, MPXV-5-TLR4 and MPXV-5-HLA-A*02:01 complexes; (E) RMSF of MPXV-2 chain in MPXV-5-TLR2, MPXV-5-TLR4 and MPXV-5-HLA-A*02:01 complexes; (F) RMSF of MPXV-5 chain in MPXV-5-TLR2, MPXV-5-TLR4 and MPXV-5-HLA-A*02:01; (G) Hydrogen bonds in MPXV-2-TLR2, MPXV-2-TLR4 and MPXV-2-HLA-A*02:01 complexes; (H) Hydrogen bonds in MPXV-5-TLR2, MPXV-5-TLR4 and MPXV-5-HLA-A*02:01 complexes.

The radius of gyration (Rg) is defined as the distribution of atoms of a protein around its axis and can be used to characterize the compactness of molecules. To some extent, it can characterize the overall changes in the ligand structure upon binding to the receptor protein. In this study, the RG values of MPXV-2-TLR4 and MPXV-2-HLA-A*02:01 complex systems were all kept at a relatively constant level (Figure 6C), indicating that the protein complex system was relatively stable with no obvious conformational change. The RG value of the MPXV-2-TLR2 complex initially maintained a stable level, began to decline at about 9 ns and reached a stable level again at about 22 ns (Figure 6C). The RG values of all MPXV-5-receptors complexes (MPXV-5-TLR2, MPXV-5-TLR4, MPXV-5-HLA-A*01:01) fluctuated slightly in the initial phase, and keep stable at about 20ns of the MD simulation (Figure 6D).

RMSF is a measurement of individual residue flexibility and the higher the RMSF value indicates better flexibility. The RMSF value of MPXV-2 chain in complex with TLR2 and TLR4 had less magnitude of fluctuation with an average value of 0.40 ± 0.11nm and 0.35 ± 0.12 nm, respectively (Figure 6E). In the case of MPXV-2 chain complexed with HLA-A*02:01, the RMSF value of residues 228 to 317 showed significant fluctuations, reaching a peak value of 1.51nm (Figure 6E). The RMSF value of MPXV-5 chain in complex with TLR2 had little magnitude of fluctuation with an average value of 0.33 ± 0.09 nm (Figure 6F). In the case of MPXV-5 chain in the MPXV-5-TLR4 and MPXV-5-HLA-A-02:01 complexes, residues 386-389 and residues 235-241 showed the largest magnitude of fluctuations, respectively, and the average RMSF value of them was 0.33 ± 0.13nm and 0.29 ± 0.09 nm, respectively (Figure 6F). For MPXV-2-receptors complexes, the RMSF values of almost all residues of the TLR2 chain were less than 0.5 nm, except for a few residues, such as residues 34-36 and 377-381(Supplementary Figure 6); for the TLR4 chain, the maximum and minimum RMSF values were 0.5838 nm and 0.1759 nm, respectively, and the RMSF values of bulk amino acids were less than 0.5nm (Supplementary Figure 6); for both A chain and B chain of HLA-A-02:01, the RMSF values showed major fluctuations in almost all residues. For MPXV-5-receptor complexes, the RMSF values of the TLR2 chain and HLA-A*02:01 chain had less magnitude of fluctuation, with an average value of 0.24nm ± 0.06 and 0.24 ± 0.05/0.30 ± 0.07(A chain/B chain), respectively (Supplementary Figure 6). In the MPXV-5-TLR4 complex, the maximum and minimum RMSF values in the TLR4 chain were 0.66 and 0,16, with a higher average value of 0.32 ± 0.09 nm (Supplementary Figure 6).

Hydrogen bonds is a nonbonded interaction, which is essential in maintaining the overall stability of the vaccine-receptor complex. In general, the number of hydrogen bonds in all complexes except the MPXV-5-TLR4 complex increased with time. The number of hydrogen bonds in MPXV-2-TLR2, MPXV-2-HLA-A02:01, MPXV-2-TLR4, MPXV-5-TLR2, MPXV-5-HLA-A02:01, MPXV-5-TLR4 complexes finally stabilized at approximately 7, 12, 8, 11,15 and 8, respectively (Figures 6G, H). The calculations of RMSD, RG, RMSF and hydrogen bonds demonstrated the stability and stiffness of the vaccine-TLRs and vaccine-HLA-A-02:01 complexes.

3.8.2 Gibb’s free energy analysis

The low Gibb’s free energy value corresponds to the stable conformation and system. In this study, the calculation of Gibb’s free energy value was based on the two principal components: the two-dimensional projections of PC1 vs PC2 were analyzed. As shown in Figure 7 and Supplementary Figure 7, among the six complexes, MPXV-5-TLR4, MPXV-2-HLA-A*02:01 and MPXV-2-TLR4 complexes have more lowest energy basins, followed by MPXV-5-TLR3 and MPXV-5-HLA-A*02:01, MPXV-5-TLR2 complexes, while MPXV-5-TLR2 complex has the fewer lowest energy basins. The lowest energy conformations of the MPXV-2-TLR2 complex were limited to the narrow range of energy basins -8.7 to 10.8 kJ/mol on PC1 and -5.6 to 7.5 kJ/mol on PC2(Figure 7A). For the MPXV-5-TLR2 complex, the energy range for PC1 was -11.9 to 10.7 kJ/mol and -8.6 to 10.8 kJ/mol for PC2 (Figure 7B). For the MPXV-2-TLR4 complex, the energy range for PC1 was -16.1 to 14.1 kJ/mol and -11.6 to 12.7 kJ/mol for PC2 (Figure 7C). For the MPXV-5-TLR4 complex, the energy range was -8.8 to 6.6 kJ/mol for PC1 and -5.8 to 7.4 kJ/mol for PC2 (Figure 7D). For the MPXV-2- HLA-A*02:01 complex, the energy range was wider, with PC1 having an energy range of -16.1 to 15.2 kJ/mol and PC2 having an energy range of -11.4 to 12.7 kJ/mol (Supplementary Figure 7A). For the MPXV-5-HLA-A*02:01 complex, the energy range for PC1 was -8.9 to 16.0 kJ/mol and -6.2 to 8.2 kJ/mol for PC2 (Supplementary Figure 7B).

FIGURE 7
www.frontiersin.org

Figure 7 Gibb’s free energy landscape. (A) MPXV-2-TLR2 complex. (B) MPXV-5-TLR2 complex. (C) MPXV-2-TLR4 complex. (D) MPXV-5-TLR4 complex.

3.8.3 MM-PBSA calculation

The trajectories extracted at 100 ps during the 40-50-ns simulation periods were submitted to MM-PBSA calculation. As shown in Table 6, the total binding free energies (ΔTOTAL) of MPXV-2-TLR2, MPXV-2-TLR4, MPXV-2-HLA-A*02:01, MPXV-5-TLR2, MPXV-5-TLR4, MPXV-5-HLA- A*02:01 complexes were -143.12 kcal/mol, -129.11 kcal/mol, -126.44 kcal/mol, -103.64 kcal/mol, -184.82kcal/mol, and -100.55 kcal/mol, respectively, indicating that both MPXV-2 and MPXV-5 have good affinity to TLR2, TLR4 and HLA-A*02:01. The analysis of detailed interaction energies was shown in Table 6.

TABLE 6
www.frontiersin.org

Table 6 MM-PBSA calculations of vaccines-receptors complexes.

3.8.4 Dynamic cross-correlation analysis

The correlations between the movements of different residues of the two chains in the vaccine-receptor complex systems were shown with the DCCM generated by DDC analysis. In the figure of DCCM, the direction and strength of the correlation between different residues were shown in different colors, and the red, white and blue corresponded to the correlation coefficients of 1 0, and -1, respectively. As shown in Supplementary Figure 8, both positively and negatively correlated residue contacts were detected between all vaccine chains and all receptor chains. In the MPXV-2-TLR2 complex, the residues around 51-151 (residues 600-700 in the figure) in the vaccine chain showed a slightly positive correlation with the residues 240-500 (residues 240-500 in the figure) from the TLR2 chain (Supplementary Figure 8A). In the MPXV-2-TLR4 complex, the relatively strong positive correlations between the residues 49-149 and 299-339 (residues 650-750, 900-940 in the figure) from the vaccine chain and the TLR4 chain were observed (Supplementary Figure 8B). In the MPXV-2-HLA-A*02:01 complex, the residues 28-128 (residues 400-500 in the figure) from the vaccine chain showed a positive correlation with the residues 1-200 (residues 1-200 in the figure) from the A chain of HLA-A*02:01 (Supplementary Figure 8C). In the MPXV-5-TLR2 complex, the residues around 26-46, 61-111 and 251-331(residues 575-595, 610-650 and 800-880 in the figure) from the vaccine chain showed positive correlations with the residues 1-300 (residues 1-300 in the figure) from the TLR2 chain (Supplementary Figure 8D). In the MPXV-5-TLR4 complex, the residues around 9-109 and 249-299 (residues 610-710 and 850-900 in the figure) from the vaccine chain showed positive correlations with the residues 1-250 (residues 1-250 in the figure) from the TLR4 chain (Supplementary Figure 8E). In the MPXV-5-HLA-A*02:01 complex, the residues around 28-118 (residues 400-490 in the figure) from the vaccine chain showed positive correlations with the residues 1-190 (residues 1-190 in the figure) from the HLA-A*02:01 chain (Supplementary Figure 8F).

3.9 Population coverage

Based on the HLA alleles binding to T-cell epitopes in the multi-epitope vaccine, the population coverage predictions were performed to determine the proportion of the population in different parts of the world that would respond to the designed vaccines. As shown in Figure 8, the population coverage of the multi-epitope vaccine was high in countries suffering relatively severe MPXV attacks in this outbreak, such as Europe (Spain, France, Germany, the United Kingdom, the Netherlands, Portugal), North America (Mexico, Canada, the United States), and South America (Brazil, Peru). The analyses of population coverage were suggestive of the fact that the designed vaccine candidates (MPXV-2 and MPXV-5) have the potential to combat the MPXV infection globally.

FIGURE 8
www.frontiersin.org

Figure 8 Population coverage of the designed vaccine in different regions suffering relatively severe MPXV attacks is calculated by the IEDB population coverage tool.

3.10 Immune simulation

The immune simulation results of MPXV-2, MPXV-5 and control group (C1) were shown in Figure 9; Supplementary Figure 9. Similar immune response patterns were observed between the simulations of candidate vaccines (MPXV-2 and MPXV-5) and C1. The primary response in the three simulations were both characterized by elevated immunoglobulin IgM (Figure 9A; Supplementary Figure 9A), B cells (isotype: lgM) (Figure 9B; Supplementary Figure 9B) and plasma cells (isotype: lgM) (Supplementary Figure 9C). In secondary and tertiary reactions, these elevations were elevated to a greater degree, and the elevation of different immunoglobulins and the immunocomplexes (IgM+IgG, IgG1+IgG2, IgG1, memory B cells and plasma cells (isotype: lgG1, lgG2) were also observed (Figures 9B, C; Supplementary Figure 9B, C). As shown in Figures 9A–C; Supplementary Figure 9A–C the titers of antibody (IgM + IgG, IgG1 + IgG2) and the level of plasma B cell (isotype: IgG1) in MPXV-2 group were higher than those in C1, while the titer of antibody (IgM + IgG, IgG1 + IgG2) and the level of plasma B cell (isotype: IgG1) in MPXV-5 group were lower than those in C1. Additionally, the helper T-cells along with corresponding memory cells were observed for two vaccine candidates and C1 (Figure 9D; Supplementary Figure 9D). The levels of active cytotoxic T cells and resting cytotoxic T cells in MPXV-2 and MPXV-5 groups changed significantly with time, while the resting cytotoxic T cell level in C1 changed little, and no inactive cytotoxic t cells were observed in the C1 group (Figures 9E; Supplementary Figure 9E). High levels of IFN-γ and IL-2 were observed in vaccine candidate groups and C1 group in the cytokines’ simulation plot (Figure 9F; Supplementary Figure 9F), while the concentration level of IL-2 was slightly higher for MPXV-2. During the immune simulation, the increased level of NK cells, macrophage (MA) and DC cells were also detected in three groups (Figure 9G–I; Supplementary Figure 9G–I). Together, the simulation results suggest that the MPXV-2 and MPXV-5 both can induce stronger immunity in the host body, thereby effectively combating the MPXV. Moreover, we can know that the humoral immune response provoked by MPXV-2 was much stronger than that of MPXV-5.

FIGURE 9
www.frontiersin.org

Figure 9 The plots relative to the immune stimulation of MPXV-5. (A) Titers of immunoglobulins and the immunocomplexes after vaccination. (B) Levels of B-cell population after vaccination. (C) Levels of plasma B-cell after vaccination. (D) Levels of helper T-cell cell population after vaccination. (E) Levels of cytotoxic T cells population after vaccination. (F) Concentration of cytokines and interleukins after vaccination. (G) Levels of NK cell population after vaccination. (H) Levels of MA population after vaccination. (I) Levels of DC population after vaccination.

3.11 Codon optimization and in silico cloning

For the optimized codon sequences of MPXV-2 and MPXV-5 (Supplementary Data), the CAI were 0.97 and 0.94 respectively, and the GC content was 48.79% and 48.18%, respectively, showing that both the MPXV-2 and MPXV-5 might well express in E. coli K12 strain and have the potential for mass production. The optimized codon sequences of MPXV-2 and MPXV-5 were inserted into the pET28a (+) vector (Supplementary Figure 10). Electrophoretogram of enzyme digestion of pET28a (+)-MPXV-5 was shown in (Supplementary Figure 11).

4 Discussion

Since May 2022, monkeypox has broken the traditional prevalent situation, beginning to spread across several countries at the same time (4). However, no available therapeutics and preventions are available for MPXV infections highlighting the need for the development of a novel vaccine against MPXV. The multi-epitope vaccine can elicit specific and strong immune responses based on multiple conserved epitopes in several antigenic proteins, avoiding adverse reactions against undesirable epitopes (88). The application of immunoinformatics methods to design multi-epitope vaccines is becoming increasingly popular, which can significantly save time and expense in vaccine development. The selection of the target protein will deeply affect the vaccine’s efficacy. A35R, B6R and H3L, which were needed for viral attachment and entry and targeted by neutralizing antibodies (nAbs), were potential targets for diagnostic tests and vaccine/drug development (89). It is well known that mutations in the strain will weaken the effectiveness of the vaccine. Hence, to design a vaccine as effective as possible against the current epidemic MPXV strains and most MPXV strains, we contrived a multi-epitope vaccine from three reference strains of monkeypox at different epidemic stages as objects and A35R, B6R and H3L as targets.

Polypeptide-combined T-cell and B-cell epitopes often equip much stronger immunogenicity than T-cell or B-cell epitopes alone. So, to maximize the immunogenicity of the designed vaccine, every vaccine construct involved both T-cell epitopes (CTL and HTL epitopes) and B-cell epitopes of each target protein. Moreover, all selected epitopes were antigenic, non-allergenic, non-toxic and non-homologous with the human proteome, indicating that the multi-epitope vaccine can induce safe antigenic response but no cross-reaction with proteins in humans. Cytokines are a category of signaling molecules that play an important role in mediating and regulating immunity and inflammation. IFN-γ, IL-4 and IL-2 are several important cytokines produced by the immune system, which contribute to the differentiation of T cells and mediate cytotoxicity and antibody production (2628). Hence, all antigenic HTL epitopes were assessed for the potential to induce the IFN-γ, IL-2 and IL-4 secretion, and only HTL epitopes that were capable of mediating the release of IFN-γ and IL-4 or IL-2 were considered for vaccine construction. Notably, the molecular docking between the T-cell epitopes and MHC molecules was performed, and the docking results revealed that the T-cell epitopes could successfully bind to MHC molecules, thereby being recognized by the antigen-presenting cell.

To acquire the optional vaccine candidates, five adjuvants were applied in this study. Although five selected adjuvants have been widely used in the design of multi-epitope vaccines in silico (4851) and have been shown to stimulate the host’s immune response in animal experiments (43, 45, 5254), any of these adjuvants have not been used for human vaccination. In the future, we can perform experiments to compare the immune enhancement effects of the adjuvants used in human vaccination, such as aluminum hydroxide and MF59 (90), and these adjuvants on multi-epitope vaccines. After the analysis of predicted biophysical and biochemical features of each vaccine construct (MPXV-1–5), two antigenic, safe, stable, thermostable and hydrophilic and soluble upon expression vaccine constructs (MPXV-2 and MPXV-5) were regarded as excellent vaccine candidates and subjected to further analysis.

Molecular docking analysis revealed strong binding affinities between MPXV-2 and MPXV-5 with TLRs (TLR2 and TLR4) and MHC molecules, which confirmed the capacity of the multi-epitope vaccine to be recognized by the human immune system, inducing stable and strong immune responses. MD simulation analysis indicated the MPXV-2-receptors and MPXV-5-receptors complexes’ stability through varying environmental conditions, including changing pressure, temperature, and motion. The initial trajectory evaluations, including the calculation of RMSD, RG, RMSF and hydrogen bonds, all correspond to the high stability of the vaccine-TLR and vaccine-MHC complexes in a biological environment. Gibb’s free energy analysis indicated that the MPXV-5-TLR4 had more stable conformations than other complexes. As shown in MM-PBSA calculations, the values of binding free energy of all complexes were negative, showing that all vaccine-receptor complexes could remain stable under natural conditions, especially the MPXV-5-TLR4 complex. Also, compared with other complexes, the binding mode of MPXV-5-TLR4 might be the most stable under natural conditions.

Notably, comparative immune stimulation was conducted between vaccine candidates and positive control. Similar immune response patterns were observed between the simulations of vaccine candidates (MPXV-2 and MPXV-5) and control group (C1), and the results of immune stimulation revealed that the MPXV-2 and MPXV-5 both can induce protective immune responses in the host. However, immune simulation performed by the C-IMMSIM server has simplified the process of antigen presentation. In general, the CTL and HTL epitopes are processed and presented by two distinct, proteasomal degradation (endogenous) and endo-lysosomal degradation (exogenous) pathways, respectively (9193). The multi-epitope vaccine is considered an extracellular antigen. After being injected into the host body, the vaccine mainly activates CD4+T cells through the exogenous pathway to induce adaptive immunity and then activates B cells to generate humoral immunity. In the field of synthetic peptide vaccines, because the clinical benefits of early direct injection of synthetic peptides or the synthetic peptides loaded on DC cells are relatively inefficient, the use of complete proteins or long peptides to activate cellular immunity by cross-presentation had attracted more attention (94). The designed multi-epitope vaccines (MPXV-2 and MPXV-5) belong to long peptides. Therefore, the CTL epitopes from MPXV-2 and MPXV-5 can be cross-presented by DCs, thus activating CD8+T cells and generating cellular immunity. However, the cellular immunity induced by multi-epitope vaccines cross-presented by DC cells may not be strong enough. The high salt for mulation of Al(OH)3 can be selected as an adjuvant, which can significantly enhance the cross-presentation of extracellular antigens by DC cells (95). However, in the future, in vitro and in vivo experiments still need to verify that the multi-epitope vaccine with a suitable adjuvant can trigger an effective cellular immune response

5 Conclusion

Since May 2022, the United Kingdom, the United States, France, Spain, Brazil and other countries have successively reported their “first” cases of monkeypox, followed by a continuous increase in the number of cases. However, no specific treatment or vaccine for MPXV is currently available on the market. In this study, five multi-epitope vaccine constructs were developed against MPXV using immunoinformatics approaches, and the MPXV-2 and MPXV-5 were identified as the most promising vaccine constructs after the immunological and physiochemical analyses. The molecular docking and MD simulation confirmed that the MPXV-2 and MPXV-5 can form stable bindings to immune receptors and thereby trigger effective immune responses against MPXV. Additionally, the results obtained from immune simulation suggested that the MPXV-2 and MPXV-5 equip the capacity to elicit immune responses in the human body, but the antibody response induced by the MPXV-2 was stronger. Having said all this, the MPXV-2 and MPXV-5 could, in theory, effectively combat the MPXV.

6 Limitation

First, the immunogenicity of multi-epitope vaccines is weaker than that of inactivated or attenuated vaccines. The application of other available adjuvants, such as aluminum hydroxide, may address this issue (95, 96). Second, the immunogenicity of multi-epitope vaccines can be affected by the expression system, necessitating careful selection of the expression system during vaccine preparation. Although animal experiments using immunoinformatics to construct multi-epitope vaccines have shown promising results, and several vaccines designed using this approach are currently in clinical studies, as an emerging technology, more in-depth studies are still needed to verify its effectiveness. Finally, its protective effect needed to be confirmed by in vitro and in vivo experiments.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

Author contributions

CT, FZ, PP, AW and CL designed the research and analyzed the data; CT and FZ charged the drawing; CT, FZ, PP, AW and CL wrote the paper. PP, AW and CL reviewed and edited the paper. CL, AW and PP acquired the funding. All authors contributed to the article and approved the submitted version.

Funding

This research was supported by the National Key Research and Development Program of China (No. 2022YFC2009801), the Natural Science Foundation of Hunan Province (No. 2021JJ31071), the National Clinical Research Center for Geriatric Disorders (Xiangya Hospital, No. 2021KFJJ05) and the Research Project on Education and Teaching Reform of Central South University (No. 2021 jy139-2). Health Development Research Center of the National Health Commission, "Evidence-based Evaluation and Demonstration Base Construction Project of Infection Control Measures in Healthcare Institutions"(CNHDRC-KJ-L-2020-53-04375).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2023.1112816/full#supplementary-material

References

1. Weaver JR, Isaacs SN. Monkeypox virus and insights into its immunomodulatory proteins. Immunol Rev (2008) 225:96–113. doi: 10.1111/j.1600-065X.2008.00691.x

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Pourriyahi H, Aryanian Z, Afshar ZM, Goodarzi A. A systematic review and clinical atlas on mucocutaneous presentations of monkeypox: With a comprehensive approach to all aspects of the new and previous monkeypox outbreaks. J Med Virol (2022) 95(2):e28230. doi: 10.1002/jmv.28230

CrossRef Full Text | Google Scholar

3. Petersen E, Kantele A, Koopmans M, Asogun D, Yinka-Ogunleye A, Ihekweazu C, et al. Human monkeypox: Epidemiologic and clinical characteristics, diagnosis, and prevention. Infect Dis Clin North Am (2019) 33(4):1027–43. doi: 10.1016/j.idc.2019.03.001

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Thornhill JP, Barkati S, Walmsley S, Rockstroh J, Antinori A, Harrison LB, et al. Monkeypox virus infection in humans across 16 countries - April-June 2022. N Engl J Med (2022) 387(8):679–91. doi: 10.1056/NEJMoa2207323

PubMed Abstract | CrossRef Full Text | Google Scholar

5. World Health Organization. Multi-country outbreak of monkeypox, external situation report 10 (Accessed 27 November 2022).

Google Scholar

6. McCarthy MW. Therapeutic strategies to address monkeypox. Expert Rev Anti Infect Ther (2022) 20(10):1249–52. doi: 10.1080/14787210.2022.2113058

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Fine PE, Jezek Z, Grab B, Dixon H. The transmission potential of monkeypox virus in human populations. Int J Epidemiol (1988) 17(3):643–50. doi: 10.1093/ije/17.3.643

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Gruber MF. Current status of monkeypox vaccines. NPJ Vaccines (2022) 7(1):94. doi: 10.1038/s41541-022-00527-4

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Lennerz V, Gross S, Gallerani E, Sessa C, Mach N, Boehm S, et al. Immunologic response to the survivin-derived multi-epitope vaccine Emd640744 in patients with advanced solid tumors. Cancer Immunol Immunother (2014) 63(4):381–94. doi: 10.1007/s00262-013-1516-5

CrossRef Full Text | Google Scholar

10. Oli AN, Obialor WO, Ifeanyichukwu MO, Odimegwu DC, Okoyeh JN, Emechebe GO, et al. Immunoinformatics and vaccine development: An overview. Immunotargets Ther (2020) 9:13–30. doi: 10.2147/itt.S241064

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Hatmal MM, Al-Hatamleh MAI, Olaimat AN, Ahmad S, Hasan H, Ahmad Suhaimi NA, et al. Comprehensive literature review of monkeypox. Emerg Microbes Infect (2022) 11(1):2600–31. doi: 10.1080/22221751.2022.2132882

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Hooper JW, Custer Dm, Thompson E. Four-Gene-Combination DNA vaccine protects mice against a lethal vaccinia virus challenge and elicits appropriate antibody responses in nonhuman primates. Virology (2003) 306(1):181–95. doi: 10.1016/s0042-6822(02)00038-7

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Roper RL, Wolffe Ej, Weisberg A, Moss B. The envelope protein encoded by the A33r gene is required for formation of actin-containing microvilli and efficient cell-to-Cell spread of vaccinia virus. J Virol (1998) 72(5):4192–204. doi: 10.1128/jvi.72.5.4192-4204.1998

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Paran N, Lustig S. Complement-bound human antibodies to vaccinia virus B5 antigen protect mice from virus challenge. Expert Rev Vaccines (2010) 9(3):255–9. doi: 10.1586/erv.10.5

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Davies DH, McCausland Mm, Valdez C, Huynh D, Hernandez JE, Mu Y, et al. Vaccinia virus H3l envelope protein is a major target of neutralizing antibodies in humans and elicits protection against lethal challenge in mice. J Virol (2005) 79(18):11724–33. doi: 10.1128/jvi.79.18.11724-11733.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Drexler I, Staib C, Kastenmuller W, Stevanović S, Schmidt B, Lemonnier FA, et al. Identification of vaccinia virus epitope-specific hla-a*0201-Restricted T cells and comparative analysis of smallpox vaccines. Proc Natl Acad Sci USA (2003) 100(1):217–22. doi: 10.1073/pnas.262668999

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Isidro J, Borges V, Pinto M, Sobral D, Santos JD, Nunes A, et al. Phylogenomic characterization and signs of microevolution in the 2022 multi-country outbreak of monkeypox virus. Nat Med (2022) 28(8):1569–72. doi: 10.1038/s41591-022-01907-y

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Waterhouse AM, Procter JB, Martin DM, Clamp M, Barton GJ. Jalview version 2–a multiple sequence alignment Editor and analysis workbench. Bioinformatics (2009) 25(9):1189–91. doi: 10.1093/bioinformatics/btp033

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Larsen MV, Lundegaard C, Lamberth K, Buus S, Lund O, Nielsen M. Large-Scale validation of methods for cytotoxic T-lymphocyte epitope prediction. BMC Bioinf (2007) 8:424. doi: 10.1186/1471-2105-8-424

CrossRef Full Text | Google Scholar

20. Doytchinova IA, Flower DR. Vaxijen: A server for prediction of protective antigens, tumour antigens and subunit vaccines. BMC Bioinf (2007) 8:4. doi: 10.1186/1471-2105-8-4

CrossRef Full Text | Google Scholar

21. Sharma NA-O, Naorem LA-O, Jain SA-O, Raghava GA-O. Toxinpred2: An improved method for predicting toxicity of proteins. Brief Bioinform (2022) 23(5):1–12. doi: 10.1093/bib/bbac174

CrossRef Full Text | Google Scholar

22. Magnan CN, Zeller M, Kayala MA, Vigil A, Randall A, Felgner PL, et al. High-throughput prediction of protein antigenicity using protein microarray data. Bioinformatics (2010) 26(23):2936–43. doi: 10.1093/bioinformatics/btq551

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Calis JJ, Maybeno M, Greenbaum JA, Weiskopf D, De Silva AD, Sette A, et al. Properties of mhc class I presented peptides that enhance immunogenicity. PloS Comput Biol (2013) 9(10):e1003266. doi: 10.1371/journal.pcbi.1003266

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Paul S, Sidney J, Sette A, Peters B. Tepitool: A pipeline for computational prediction of T cell epitope candidates. Curr Protoc Immunol (2016) 114:18.9.1-.9.24. doi: 10.1002/cpim.12

CrossRef Full Text | Google Scholar

25. Reynisson B, Barra C, Peters B, Nielsen M. Improved prediction of mhc ii antigen presentation through integration and motif deconvolution of mass spectrometry mhc eluted ligand data. J Proteome Res (2020) 19(6):2304–15. doi: 10.1021/acs.jproteome.9b00874

CrossRef Full Text | Google Scholar

26. Kak G, Raza M, Tiwari BK. Interferon-gamma (Ifn-Γ): Exploring its implications in infectious diseases. Biomol Concepts (2018) 9(1):64–79. doi: 10.1515/bmc-2018-0007

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Smiley ST, Grusby MJ. Interleukin 4. In: Delves PJ, editor. Encyclopedia of immunology, 2nd ed., vol. . p . Oxford: Elsevier (1998). p. 1451–3.

Google Scholar

28. Olejniczak K, Kasprzak A. Biological properties of interleukin 2 and its role in pathogenesis of selected diseases–a review. Med Sci Monit (2008) 14(10):Ra179–89.

PubMed Abstract | Google Scholar

29. Dhanda SK, Vir P, Raghava GPS. Designing of interferon-gamma inducing mhc class-ii binders. Biol Direct (2013) 8(1):30. doi: 10.1186/1745-6150-8-30

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Dhanda SK, Gupta S, Vir P, Raghava GPS. Prediction of Il4 inducing peptides. Clin Dev Immunol (2013) 2013:263952. doi: 10.1155/2013/263952

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Lathwal A, Kumar R, Kaur D, Raghava GPS. In silico model for predicting il-2 inducing peptides in human. bioRxiv (2021). doi: 10.1101/2021.06.20.449146

CrossRef Full Text | Google Scholar

32. Saha S, Raghava GP. Prediction of continuous b-cell epitopes in an antigen using recurrent neural network. Proteins (2006) 65(1):40–8. doi: 10.1002/prot.21078

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Bui H-H, Sidney J, Dinh K, Southwood S, Newman MJ, Sette A. Predicting population coverage of T-cell epitope-based diagnostics and vaccines. BMC Bioinf (2006) 7(1):153. doi: 10.1186/1471-2105-7-153

CrossRef Full Text | Google Scholar

34. Chen C, Li Z, Huang H, Suzek BE, Wu CH. A fast peptide match service for uniprot knowledgebase. Bioinformatics (2013) 29(21):2808–9. doi: 10.1093/bioinformatics/btt484

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Lamiable A, Thévenet P, Rey J, Vavrusa M, Derreumaux P, Tufféry P. Pep-Fold3: Faster de Novo structure prediction for linear peptides in solution and in complex. Nucleic Acids Res (2016) 44(W1):W449–54. doi: 10.1093/nar/gkw329

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Sussman JL, Lin D, Jiang J, Jiang J, Manning NO, Prilusky J, et al. Protein data bank (Pdb): Database of three-dimensional structural information of biological macromolecules. Acta Crystallogr D Biol Crystallogr (1998) 54(Pt 6 Pt 1):1078–84. doi: 10.1107/s0907444998009378

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Kaplan W, Littlejohn TG. Swiss-Pdb viewer (Deep view). Brief Bioinform (2001) 2(2):195–7. doi: 10.1093/bib/2.2.195

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Eberhardt J, Santos-Martins DA-O, Tillack AA-O, Forli SA-O. Autodock vina 1.2.0: New docking methods, expanded force field, and Python bindings. J Chem Inf Model (2021) 61(8):3891–8. doi: 10.1021/acs.jcim.1c00203

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Trott O, Olson AJ. Autodock vina: Improving the speed and accuracy of docking with a new scoring function, efficient optimization, and multithreading. J Comput Chem (2010) 31(2):455–61. doi: 10.1002/jcc.21334

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Laskowski RA, Swindells MB. Ligplot+: Multiple ligand-protein interaction diagrams for drug discovery. J Chem Inf Model (2011) 51(10):2778–86. doi: 10.1021/ci200227u

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Chen X, Zaro Jl Fau - Shen W-C, Shen WC. Fusion protein linkers: Property, design and functionality. Adv Drug Deliv Rev (2013) 65(10):1357–69. doi: 10.1016/j.addr.2012.09.039

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Marciani DJ. Vaccine adjuvants: Role and mechanisms of action in vaccine immunogenicity. Drug Discov Today (2003) 8(20):934–43. doi: 10.1016/S1359-6446(03)02864-2

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Lee SJ, Shin SJ, Lee MH, Lee MG, Kang TH, Park WS, et al. A potential protein adjuvant derived from mycobacterium tuberculosis Rv0652 enhances dendritic cells-based tumor immunotherapy. PloS One (2014) 9(8):e104351. doi: 10.1371/journal.pone.0104351

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Yang D, Chertov O, Bykovskaia Sn, Chen Q, Buffo MJ, Shogan J, et al. Beta-defensins: Linking innate and adaptive immunity through dendritic and T cell Ccr6. Science (1999) 286(5439):525–8. doi: 10.1126/science.286.5439.525

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Hu JC, Mathias-Santos C, Greene CJ, King-Lyons ND, Rodrigues JF, Hajishengallis G, et al. Intradermal administration of the type ii heat-labile enterotoxins Lt-iib and Lt-iic of enterotoxigenic escherichia coli enhances humoral and Cd8+ T cell immunity to a Co-administered antigen. PloS One (2014) 9(12):e113978. doi: 10.1371/journal.pone.0113978

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Antonio-Herrera L, Badillo-Godinez O, Medina-Contreras O, Tepale-Segura A, García-Lozano A, Gutierrez-Xicotencatl L, et al. The nontoxic cholera b subunit is a potent adjuvant for intradermal dc-targeted vaccination. Front Immunol (2018) 9:2212. doi: 10.3389/fimmu.2018.02212

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Shanmugam A, Rajoria S, George AL, Mittelman A, Suriano R, Tiwari RK. Synthetic toll like receptor-4 (Tlr-4) agonist peptides as a novel class of adjuvants. PLoS One (2012) 7(2):e30839. doi: 10.1371/journal.pone.0030839

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Aziz S, Waqas M, Halim SA, Ali A, Iqbal A, Iqbal M, et al. Exploring whole proteome to contrive multi-Epitope-Based vaccine for neocov: An immunoinformtics and in-silico approach. Front Immunol (2022) 13:956776. doi: 10.3389/fimmu.2022.956776

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Tan C, Zhu F, Xiao Y, Wu Y, Meng X, Liu S, et al. Immunoinformatics approach toward the introduction of a novel multi-epitope vaccine against clostridium difficile. Front Immunol (2022) 13:887061. doi: 10.3389/fimmu.2022.887061

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Guo L, Yin R, Liu K, Lv X, Li Y, Duan X, et al. Immunological features and efficacy of a multi-epitope vaccine ctb-ue against h. pylori in Balb/C mice model. Appl Microbiol Biotechnol (2014) 98(8):3495–507. doi: 10.1007/s00253-013-5408-6

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Ojha R, Pandey RK, Prajapati VK. Vaccinomics strategy to concoct a promising subunit vaccine for visceral leishmaniasis targeting sandfly and leishmania antigens. Int J Biol Macromol (2020) 156:548–57. doi: 10.1016/j.ijbiomac.2020.04.097

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Li M, Jiang Y, Gong T, Zhang Z, Sun X. Intranasal vaccination against hiv-1 with adenoviral vector-based nanocomplex using synthetic tlr-4 agonist peptide as adjuvant. Mol Pharm (2016) 13(3):885–94. doi: 10.1021/acs.molpharmaceut.5b00802

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Hou J, Liu Y, Hsi J, Wang H, Tao R, Shao Y. Cholera toxin b subunit acts as a potent systemic adjuvant for hiv-1 DNA vaccination intramuscularly in mice. Hum Vaccin Immunother (2014) 10(5):1274–83. doi: 10.4161/hv.28371

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Kalenik BM, Góra-Sochacka A, Sirko A. β-defensins - underestimated peptides in influenza combat. Virus Res (2018) 247:10–4. doi: 10.1016/j.virusres.2018.01.008

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Alexander J, del Guercio M-F, Maewal A, Qiao L, Fikes J, Chesnut RW, et al. Linear padre T helper epitope and carbohydrate b cell epitope conjugates induce specific high titer igg antibody responses. J Immunol (2000) 164(3):1625. doi: 10.4049/jimmunol.164.3.1625

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Gasteiger E, Gattiker A, Hoogland C, Ivanyi I, Appel RD, Bairoch A. Expasy: The proteomics server for in-depth protein knowledge and analysis. Nucleic Acids Res (2003) 31(13):3784–8. doi: 10.1093/nar/gkg563

PubMed Abstract | CrossRef Full Text | Google Scholar

57. Hallgren J, Tsirigos KD, Pedersen MD, Almagro Armenteros JJ, Marcatili P, Nielsen H, et al. Deeptmhmm predicts alpha and beta transmembrane proteins using deep neural networks. bioRxiv (2022). doi: 10.1101/2022.04.08.487609

CrossRef Full Text | Google Scholar

58. Magnan CN, Randall A, Baldi P. Solpro: Accurate sequence-based prediction of protein solubility. Bioinformatics (2009) 25(17):2200–7. doi: 10.1093/bioinformatics/btp386

PubMed Abstract | CrossRef Full Text | Google Scholar

59. Nielsen H, Tsirigos KD, Brunak S, von Heijne G. A brief history of protein sorting prediction. Protein J (2019) 38(3):200–16. doi: 10.1007/s10930-019-09838-3

PubMed Abstract | CrossRef Full Text | Google Scholar

60. Buchan DWA, Jones DT. The psipred protein analysis workbench: 20 years on. Nucleic Acids Res (2019) 47(W1):W402–w7. doi: 10.1093/nar/gkz297

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Kim DE, Chivian D Fau - Baker D, Baker D. Protein structure prediction and analysis using the robetta server. Nucleic Acids Res (2004) 32(Web Server issue):W526–31. doi: 10.1093/nar/gkh468

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Ko J, Park H Fau - Heo L, Heo L Fau - Seok C, Seok C. Galaxyweb server for protein structure prediction and refinement. Nucleic Acids Res (2012) 40(Web Server issue):W294–7. doi: 10.1093/nar/gks493

PubMed Abstract | CrossRef Full Text | Google Scholar

63. Colovos C, Yeates TO. Verification of protein structures: Patterns of nonbonded atomic interactions. Protein Sci (1993) 2(9):1511–9. doi: 10.1002/pro.5560020916

PubMed Abstract | CrossRef Full Text | Google Scholar

64. Wiederstein M, Sippl MJ. Prosa-web: Interactive web service for the recognition of errors in three-dimensional structures of proteins. Nucleic Acids Res (2007) 35(Web Server issue):W407–10. doi: 10.1093/nar/gkm290.

PubMed Abstract | CrossRef Full Text | Google Scholar

65. Vijay K. Toll-like receptors in immunity and inflammatory diseases: Past, present, and future. Int Immunopharmacol (2018) 59:391–412. doi: 10.1016/j.intimp.2018.03.002

PubMed Abstract | CrossRef Full Text | Google Scholar

66. Zhou R, Liu L, Wang Y. Viral proteins recognized by different tlrs. J Med Virol (2021) 93(11):6116–23. doi: 10.1002/jmv.27265

PubMed Abstract | CrossRef Full Text | Google Scholar

67. Sartorius R, Trovato M, Manco R, D'Apice L, De Berardinis P. Exploiting viral sensing mediated by toll-like receptors to design innovative vaccines. NPJ Vaccines (2021) 6(1):127. doi: 10.1038/s41541-021-00391-8

PubMed Abstract | CrossRef Full Text | Google Scholar

68. Kozakov D, Hall DR, Xia B, Porter KA, Padhorny D, Yueh C, et al. The cluspro web server for protein-protein docking. Nat Protoc (2017) 12(2):255–78. doi: 10.1038/nprot.2016.169

PubMed Abstract | CrossRef Full Text | Google Scholar

69. van Zundert GCP, Rodrigues JPGLM, Trellet M, Schmitz C, Kastritis PL, Karaca E, et al. The Haddock2.2 web server: User-friendly integrative modeling of biomolecular complexes. J Mol Biol (2016) 428(4):720–5. doi: 10.1016/j.jmb.2015.09.014

PubMed Abstract | CrossRef Full Text | Google Scholar

70. Vangone A, Bonvin AMJJ. Contacts-based prediction of binding affinity in protein–protein complexes. eLife (2015) 4:e07454. doi: 10.7554/eLife.07454

PubMed Abstract | CrossRef Full Text | Google Scholar

71. Laskowski RA-O, Jabłońska J, Pravda L, Vařeková RS, Thornton JM. Pdbsum: Structural summaries of pdb entries. Protein Sci (2018) 27(1):129–34. doi: 10.1002/pro.3289

PubMed Abstract | CrossRef Full Text | Google Scholar

72. Delano WL. The pymol molecular graphics system. Proteins (2002) 30:442–54.

Google Scholar

73. Wickes RD ed. Blender overview. In: Foundation blender compositing. Berkeley, CA: Apress. p. 1–23.

Google Scholar

74. Kutzner C, Páll S, Fechner M, Esztermann A, de Groot BL, Grubmüller H. Best bang for your buck: Gpu nodes for gromacs biomolecular simulations. J Comput Chem (2015) 36(26):1990–2008. doi: 10.1002/jcc.24030

PubMed Abstract | CrossRef Full Text | Google Scholar

75. Maier JA, Martinez C, Kasavajhala K, Wickstrom L, Hauser KE, Simmerling C. Ff14sb: Improving the accuracy of protein side chain and backbone parameters from Ff99sb. J Chem Theory Comput (2015) 11(8):3696–713. doi: 10.1021/acs.jctc.5b00255

PubMed Abstract | CrossRef Full Text | Google Scholar

76. Berendsen HJ, Postma Jv, Van Gunsteren WF, DiNola A, Haak JR. Molecular dynamics with coupling to an external bath. J Chem Phys (1984) 81(8):3684–90. doi: 10.1063/1.448118

CrossRef Full Text | Google Scholar

77. Parrinello M, Rahman A. Strain fluctuations and elastic constants. J Chem Phys (1982) 76(5):2662–6. doi: 10.1063/1.443248

CrossRef Full Text | Google Scholar

78. Hess B, Bekker H, Berendsen HJC, Fraaije JGEM. Lincs: A linear constraint solver for molecular simulations. J Comput Chem (1997) 18(12):1463–72. doi: 10.1002/(SICI)1096-987X(199709)18:12<1463::AID-JCC4>3.0.CO;2-H

CrossRef Full Text | Google Scholar

79. Darden T, York D, Pedersen L. Particle mesh ewald: An N·Log(N) method for ewald sums in Large systems. J Chem Phys (1993) 98(12):10089–92. doi: 10.1063/1.464397

CrossRef Full Text | Google Scholar

80. Maisuradze GG, Liwo A, Scheraga HA. Relation between free energy landscapes of proteins and dynamics. J Chem Theory Comput (2010) 6(2):583–95. doi: 10.1021/ct9005745

PubMed Abstract | CrossRef Full Text | Google Scholar

81. Valdés-Tresanco MS, Valdés-Tresanco ME, Valiente PA, Moreno E. Gmx_Mmpbsa: A new tool to perform end-state free energy calculations with gromacs. J Chem Theory Comput (2021) 17(10):6281–91. doi: 10.1021/acs.jctc.1c00645

PubMed Abstract | CrossRef Full Text | Google Scholar

82. Rapin N, Lund O, Bernaschi M, Castiglione F. Computational immunology meets bioinformatics: The use of prediction tools for molecular binding in the simulation of the immune system. PloS One (2010) 5(4):e9862. doi: 10.1371/journal.pone.0009862

PubMed Abstract | CrossRef Full Text | Google Scholar

83. Yu M, Zhu Y, Li Y, Chen Z, Sha T, Li Z, et al. Design of a novel multi-epitope vaccine against echinococcus granulosus in immunoinformatics. Front Immunol (2021) 12:668492. doi: 10.3389/fimmu.2021.668492

PubMed Abstract | CrossRef Full Text | Google Scholar

84. Grote A, Hiller K, Scheer M, Münch R, Nörtemann B, Hempel DC, et al. Jcat: A novel tool to adapt codon usage of a target gene to its potential expression host. Nucleic Acids Res (2005) 33(Web Server issue):W526–31. doi: 10.1093/nar/gki376

PubMed Abstract | CrossRef Full Text | Google Scholar

85. Chaudhri G, Quah BJ, Wang Y, Tan AH, Zhou J, Karupiah G, et al. T Cell receptor sharing by cytotoxic T lymphocytes facilitates efficient virus control. Proc Natl Acad Sci USA (2009) 106(35):14984–9. doi: 10.1073/pnas.0906554106

PubMed Abstract | CrossRef Full Text | Google Scholar

86. Cancro MP, Tomayko MM. Memory b cells and plasma cells: The differentiative continuum of humoral immunity. Immunol Rev (2021) 303(1):72–82. doi: 10.1111/imr.13016

PubMed Abstract | CrossRef Full Text | Google Scholar

87. Bacchetta R, Gregori S, Roncarolo MG. Cd4+ regulatory T cells: Mechanisms of induction and effector function. Autoimmun Rev (2005) 4(8):491–6. doi: 10.1016/j.autrev.2005.04.005

PubMed Abstract | CrossRef Full Text | Google Scholar

88. Vartak A, Sucheck SJ. Recent advances in subunit vaccine carriers. Vaccines (Basel) (2016) 4(2):12–19. doi: 10.3390/vaccines4020012

PubMed Abstract | CrossRef Full Text | Google Scholar

89. Pacchioni SM, Bissa M, Zanotto C, Morghen CDG, Illiano E, Radaelli A. L1r, A27l, A33r and B5r vaccinia virus genes expressed by fowlpox recombinants as putative novel orthopoxvirus vaccines. J Transl Med (2013) 11:95. doi: 10.1186/1479-5876-11-95

PubMed Abstract | CrossRef Full Text | Google Scholar

90. Pulendran B, Arunachalam PS, O'Hagan DT. Emerging concepts in the science of vaccine adjuvants. Nat Rev Drug Discov (2021) 20(6):454–75. doi: 10.1038/s41573-021-00163-y

PubMed Abstract | CrossRef Full Text | Google Scholar

91. Tellam J, Connolly G, Green KJ, Miles JJ, Moss DJ, Burrows SR, et al. Endogenous presentation of Cd8+ T cell epitopes from Epstein-Barr virus-encoded nuclear antigen 1. J Exp Med (2004) 199(10):1421–31. doi: 10.1084/jem.20040191

PubMed Abstract | CrossRef Full Text | Google Scholar

92. Vyas JM, van der Veen AG, Ploegh HL. The known unknowns of antigen processing and presentation. Nat Rev Immunol (2008) 8(8):607–18. doi: 10.1038/nri2368

PubMed Abstract | CrossRef Full Text | Google Scholar

93. Pishesha N, Harmand TJ, Ploegh HL. A guide to antigen processing and presentation. Nat Rev Immunol (2022) 22(12):751–64. doi: 10.1038/s41577-022-00707-2

PubMed Abstract | CrossRef Full Text | Google Scholar

94. Joffre OP, Segura E, Savina A, Amigorena S. Cross-presentation by dendritic cells. Nat Rev Immunol (2012) 12(8):557–69. doi: 10.1038/nri3254

PubMed Abstract | CrossRef Full Text | Google Scholar

95. Luo M, Shao B, Yu JY, Liu T, Liang X, Lu L, et al. Simultaneous enhancement of cellular and humoral immunity by the high salt formulation of Al(Oh)(3) adjuvant. Cell Res (2017) 27(4):586–9. doi: 10.1038/cr.2017.14

PubMed Abstract | CrossRef Full Text | Google Scholar

96. Güven E, Duus K, Laursen I, Højrup P, Houen G. Aluminum hydroxide adjuvant differentially activates the three complement pathways with major involvement of the alternative pathway. PloS One (2013) 8(9):e74445. doi: 10.1371/journal.pone.0074445

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: monkeypox virus (MPXV), multi-epitope vaccine, molecular docking, molecular dynamics (MD) simulation, immunoinformatics

Citation: Tan C, Zhu F, Pan P, Wu A and Li C (2023) Development of multi-epitope vaccines against the monkeypox virus based on envelope proteins using immunoinformatics approaches. Front. Immunol. 14:1112816. doi: 10.3389/fimmu.2023.1112816

Received: 30 November 2022; Accepted: 21 February 2023;
Published: 13 March 2023.

Edited by:

P Bernard Fourie, University of Pretoria, South Africa

Reviewed by:

Xiancai Rao, Army Medical University, China
Hardik Patel, Seattle Children’s Research Institute, United States

Copyright © 2023 Tan, Zhu, Pan, Wu and Li. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Chunhui Li, lichunhui@csu.edu.cn; Anhua Wu, xywuanhua@csu.edu.cn; Pinhua Pan, pinhuapan668@csu.edu.cn

These authors have contributed equally to this work

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.