Development of Multi-epitope Based Subunit Vaccine Against Crimean-Congo Hemorrhagic Fever Virus Using Reverse Vaccinology Approach.

Crimean-Congo hemorrhagic fever (CCHF) is a viral disease caused by the Crimean-Congo hemorrhagic fever virus (CCHFV) of the Nairovirus genus. CCHF has occurred endemically in several regions of Africa, Southern Europe, and Central and Southeast Asia, with a case fatality rate of 5 to 80%. The World health organization enlisted CCHF as one of the top prioritized diseases for research and development in emergency contexts that making it a public health concern as no effective vaccine is available till date. Therefore, the present study aims to develop an effective multi-epitope subunit vaccine using immunoinformatics and reverse vaccinology approach against this virus. The B-cell and T-cell epitopes were predicted from structural and non-structural proteins, and filtered by immunogenicity, allergenicity, toxicity, conservancy, and cross-reactivity. The computational analysis revealed that the epitopes could induce an adequate immune response and had strong associations with their respective human leukocyte antigen (HLA) alleles with 98.94% of total world population coverage. Finally, the vaccine with 427 amino acids was constructed by connecting 8 cytotoxic T-lymphocytes, 4 helper T-lymphocytes, and 10 B-cell epitopes with appropriate linkers and β-defensin as an adjuvant. The antigenicity, allergenicity, solubility, and physiochemical properties of the vaccine were evaluated, followed by structural modelling, refinement, and validation. In addition, molecular docking and molecular dynamic simulations revealed a robust binding affinity and stability of the vaccine-immune receptor complex. Moreover, the codons were optimized for its higher expression in Escherichia coli (E. coli) K12 strain followed by in silico cloning. The proposed subunit vaccine developed in this study could be a potential candidate against CCHFV. However, further experimental validation is required to ensure the immunogenicity and safety profile of the proposed vaccine for combating and eradicating CCHFV. Supplementary Information: The online version contains supplementary material available at 10.1007/s10989-022-10430-0.

Introduction

The first recognized outbreak of CCHF occurred during the World War II in Crimea. The arboviral disease is caused by the CCHFV classified within the Nairovirus genus of the Bunyaviridae family, including over 350 species. Transmission of CCHFV occurs mostly through ticks of the Hyalomma genus in several vertebrates (Hoogstraal 1979; Goldfarb et al. 1980; Zivcec et al. 2016). Infection can occur either from a tick bite or direct contact with the body fluid of the infected individual or livestock. It is asymptomatic in animals but lethal for humans with a case fatality rate of 5 to 80% (Dowall et al. 2017). There are four phases of CCHFV infection including incubation, prehemorrhagic, hemorrhagic, and convalescent (Bente et al. 2013). The symptoms are carried by fever, fatigue, frequent vomiting, and diarrhoea. But in severe cases, CCHFV causes thrombocytopenia, elevated circulating liver enzymes, and hemorrhagic syndrome (Goldfarb et al. 1980; Whitehouse 2004; Vorou et al. 2007). Some studies showed CCHF symptoms resemble Ebola hemorrhagic fever (Bente et al. 2013). Abnormal coagulation, characterized by disseminated intravascular coagulation, was found at the early stage of the illness in the patient with CCHFV (Burt et al. 1997).

CCHFV is a single-stranded negative-sense RNA virus. The virions are spherical in shape and about 90–100 nm in diameter, with approximately 10 nm long spike-like projections from the surface (Zivcec et al. 2016). It contains a tri-segmented (S, M, and L) negative-sense RNA genome encoding the nucleoprotein (NP), glycoprotein (Gc and Gn), and RNA-dependent RNA polymerase (RDRP) (Burt et al. 1997; Whitehouse 2004; Flick and Whitehouse 2005). The RDRP is encoded by the L segment of the viral genome and is responsible for the replication and transcription (Kinsella et al. 2004). The envelopment polyprotein (EPP) is a key structural component of CCHFV that is cleaved into five chains which are mucin-like variable region, glycoprotein38 (Gp38), glycoprotein N (Gn), non-structural protein M (NSm) and glycoprotein C (Gc) (Bergeron et al. 2007). The Gc and Gn are present at the virion surface, interacting with each other. They can attach the virion with the host cell receptor. Previously the disease has been reported in Sub-Saharan Africa, Southern Europe, and Central and Southeast Asia. The recent emergence of this virus has been reported in Spain, Turkey, and India. The existence of the tick vector in different regions may indicate the expansion of its geological range and risk of an outbreak (Hoogstraal 1979; Hawman and Feldmann 2018). The wide distribution and high fatality rate of the virus with no approved vaccine lead to a serious public health risk in many regions (Mertens et al. 2013). A recent study has shown the efficacy of a antiviral drug ribavirin against CCHFV, but its therapeutic benefit is still not clear (Duygu et al. 2012). As the emerging risk of an outbreak with high mortality rate, the development of an effective vaccine is obligatory.

A mouse-brain-derived CCHFV vaccine was developed in the Soviet Union after the outbreak in 1970. The virus was inactivated by chloroform, followed by heating and absorbtion on aluminium hydroxide (Papa et al. 2011). Although after receiving several doses, the vaccine could stimulate a low-level neutralizing antibody response and raised concerns due to its autoimmune responses. Recent advances created a DNA vaccine expressing the M segment of the glycoprotein precursor gene of CCHFV that shows strong humoral immune responses and protective efficacy against CCHF in a mouse model (Garrison et al. 2017). In another study, a cell culture-based vaccine was developed showing to be protective against CCHFV in mice (Canakoglu et al. 2015). A subunit vaccine was introduced with a high level of antibody induction against CCHFV, but no protective immune response was observed in mice (Kortekaas et al. 2015). Similarly, potent antibody and T cell responses were found in a non-human primate disease model (Cynomolgus macaques) after introducing a glycoprotein precursor and NP encoding DNA vaccine (Hawman et al. 2021). However, human immune responses of the constructed vaccines were not investigated yet. Recently, a formalin-inactivated and alum-formulated CCHFV vaccine was developed for human use that showed strong immune induction in the vaccinated mice in a dose dependent manner and could reduce virus infection induced mortality and morbidity (Berber et al. 2021). As a result, there is no available effective and commercially approved vaccine against this virus. Generally, RNA viruses have a high mutation rate in their genome because of the error-prone nature of their polymerase (Holland and Domingo 1998). With this wide variety of genomic diversity of the CCHFV genome, effective live-attenuated or inactivated vaccine development is challenging. Moreover, a potential drawback of the live attenuated vaccine is its possible reversion to its virulent form, which may cause severe problems after immunization.

With the advancement in genome sequencing and recombinant DNA technologies, reverse vaccinology approach can outstrip the conventional vaccine development. In this process, multi-epitope subunit vaccine is constructed by identifying potential epitopes after analyzing genomic information of a target pathogen with the aid of numerous bioinformatics tools. This approach has been embraced widely as it can effectively induce both humoral and cell mediated immune responses with high specificity and cost-efficient production (Sette et al. 2001). In this study, the protein sequences of EPP and RDRP of CCHFV (strain Nigeria/IbAr10200/1970) were used to predict T and B cell epitope-based peptides, analyze protein-ligand interactions, and design protein structure. The proteins of CCHFV were evaluated to predict the T-cell and B-cell epitopes based on their antigenicity, allergenicity, and toxicity, along with their respected major histocompatibility (MHC) alleles. The selected epitopes are then docked with HLA molecules for further screening to confirm their interaction. The screened epitopes were linked with appropriate linkers and an adjuvant was attached at the N terminal to increase the immunogenicity of the proposed vaccine. Antigenicity, allergenicity, solubility, and physicochemical properties for the designed vaccine construct were analyzed. Prediction of the tertiary structure and disulfide engineering was performed for the refinement and validation of the vaccine. The binding affinity and interaction of the refined vaccine with the human immune receptors, e.g., the toll-like receptors (TLR2, TLR3, TLR4), were studied using molecular docking. Furthermore, the interactions and stability of the vaccine-receptor complex were evaluated by molecular dynamics simulation. Finally, codon optimization followed by in silico cloning strategy was developed to ensure higher expression of the designed protein within an ideal host.

Methodology

Retrieval of Protein Sequences

The protein sequences were obtained from the Universal Protein Resource (UniProt) database (https://www.uniprot.org/) which provide a high quality manually annotated and non-redundant protein sequences, together with experimental results, computed features, and scientific conclusions. CCHFV Nigeria/IbAr10200/1970 was the only completely reviewed strain, and thus used as a reference strain for the study. The viral EPP (UniProt Entry: Q8JSZ3), RDRP (UniProt Entry: Q6TQR6), and NP (UniProt Entry: P89522) were retrieved in FASTA format (The UniProt Consortium 2017). The EPP could be cleaved into the five chains: Mucin-like variable region, G38, Gn/G2, NSm, and Gc/G1.

Protein Antigenicity Prediction and Physicochemical Characterization

The antigenic propensity of the retrieved protein sequences was evaluated through VaxiJen v2.0 (http://www.ddg-pharmfac.net/vaxijen/VaxiJen/VaxiJen.html), keeping the threshold at 0.4 (Doytchinova and Flower 2007). Different physical and chemical parameters of the target proteins, including molecular weight, theoretical isoelectric point (pI), extinction coefficient, instability index, aliphatic index, estimated half-life, and grand average of hydropathicity (GRAVY) were analyzed by using a widely used online tool of ExPASy, ProtParam (https://web.expasy.org/protparam/) (Islam et al. 2019). Subcellular localization of the targeted proteins was predicted using the CELLO2GO online prediction tool (http://cello.life.nctu.edu.tw/cello2go/) (Yu et al. 2014).

Prediction of T-Cell Epitopes

The selected protein sequences were run on NetCTL 1.2 server (http://www.cbs.dtu.dk/services/NetCTL/) with default settings (the threshold for epitope identification, weight on C terminal cleavage, and weight on transporter associated with antigen processing (TAP) transport efficiency were 0.75, 0.15, and 0.05 respectively) to predict 9-mer cytotoxic T-lymphocytes (CTL) epitopes for 12 different class-I MHC supertypes (A1, A2, A3, A24, A26, B7, B8, B27, B39, B44, B58, and B62) integrating the prediction of MHC class-I binding, proteasomal C terminal cleavage, and TAP protein transport efficiency (Larsen et al. 2007). The 15 amino acids (aa) long helper T-lymphocytes (HTL) epitopes of the selected proteins were predicted using the NN-align method of Immune Epitope Database (IEDB, http://tools.iedb.org/mhcii/) (Jensen et al. 2018). The complete set of 27 human alleles was used, covering 99% of the world population of different ethnicities and regions (Vita et al. 2015). The epitopes with half-maximal inhibitory concentration (IC50) less than 50 nM are considered high affinity, less than 500 nM as intermediate affinity, and less than 5000 nM as low affinity. In this study, IC50 value lower than 50 nM were selected for further analysis.

Prediction of B-Cell Epitopes

The ABCpred server (https://webs.iiitd.edu.in/raghava/abcpred/ABC_submission.html) was used to identify the 16-mer linear B-cell epitopes of the targeted proteins at a default threshold of 0.51 (Saha and Raghava 2006). Identification of the linear B-cell epitope was based on several parameters, including hydrophilicity, antigenic propensity, flexibility, surface accessibility, and beta-turn that were evaluated using Parker hydrophilicity prediction algorithms (Parker et al. 1986), Kolaskar and Tongaonkar Antigenicity scale (Kolaskar and Tongaonkar 1990), Karplus and Schulz flexibility prediction tool (Karplus and Schulz 1985), Emini surface accessibility prediction method (Emini et al. 1985), and Chou and Fasman beta-turn prediction algorithm (Chou and Fasman 1978), respectively. Structure-based epitopes prediction of the B-cell were performed using DiscoTope 2.0 of the IEDB server (http://tools.iedb.org/discotope/) (Kringelum et al. 2012).

Evaluation of the Predicted Epitopes

The primarily selected T-cell and B-cell epitopes were sorted based on their antigenicity, allergenicity, and toxicity profile. VaxiJen v2.0 server was used to evaluate the antigenicity of the predicted epitopes. To determine the allergic response and toxic activity of the short-listed epitopes, AllerTop v2.0 (https://www.ddg-pharmfac.net/AllerTOP/) and ToxinPred server (https://webs.iiitd.edu.in/raghava/toxinpred/design.php) were used, respectively (Gupta et al. 2013; Dimitrov et al. 2013). Induction of cytokines such as interferon-γ (IFN-γ), interleukin-4 (IL-4), and interleukin-10 (IL-10) by the predicted epitopes was analysed using IFNepitope (https://webs.iiitd.edu.in/raghava/ifnepitope/developer.php) (Dhanda et al. 2013b), IL4pred (https://webs.iiitd.edu.in/raghava/il4pred/design.php), and IL10pred (https://webs.iiitd.edu.in/raghava/il10pred/predict3.php) web servers (Dhanda et al. 2013a).

For conservancy analysis of the selected CTL and HTL epitopes among different viral serotypes, the BLASTp tool (https://blast.ncbi.nlm.nih.gov/Blast.cgi) was used to identify the homologous sequences of the selected antigenic proteins from the NCBI database. Further, the IEDB Epitope Conservancy Analysis server (http://tools.iedb.org/conservancy/) was used to calculate the epitopes conservancy level. To check the cross-reactivity with human proteome, the PIR peptide matching program (https://research.bioinformatics.udel.edu/peptidematch/batchpeptidematch.jsp) was used (Chen et al. 2013).

Prediction of MHC-I and MHC-II Binding Alleles and Population Coverage Analysis

To predict the MHC class-I and MHC class-II binding alleles for each CTL and HTL epitopes, respectively, IEDB recommended method and consensus prediction method (Percentile rank ≤ 5) and was employed (Wang et al. 2008). The population coverage tool from IEDB (http://tools.iedb.org/population/) was used to calculate the percentage of population coverage of the selected CTL and HTL epitopes and their MHC binding alleles, both individually and combined.

Prediction of the 3D Structure and Molecular Docking of the Selected Epitopes

The three-dimensional structures of top-ranked T-cell epitopes were generated using an online peptide structure prediction tool, namely PEP-FOLD (https://bioserv.rpbs.univ-paris-diderot.fr/services/PEP-FOLD/) (Shen et al. 2014). Based on the results of epitope-allele interaction, we considered HLA-A*01:01 (PDB ID: 4nqv) and HLA-DRB1*15:01 (PDB ID: ID: 1BX2) for docking with MHC class-I and class-II epitopes, respectively. First, the protein was prepared using BIOVIA Discovery Studio Visualizer (Biovia et al. 2000) and then, the molecular docking was conducted in the PyRx virtual screening tool through Autodock Vina (Trott and Olson 2010; Dallakyan and Olson 2015). The size of the grid box was set at 41.4218, 62.1016, 47.2505 Å for X, Y, and Z axes, respectively, of HTL epitopes and HLA-DRB1*15:01 interaction and 53.2882, 68.0418, and 61.4209 Å along X, Y, and Z axes, respectively for HLA-A*01:01 and CTL epitopes. The binding pattern and interactions were analyzed using Discovery Studio Visualizer. Moreover, for MHC class-I, a 9 aa long influenza peptide NP44 (“CTELKLSDY”) and MHC class-II, a 15 aa long human myelin basic protein (“ENPVVHFFKNIVTPR”) were used as control peptide for docking study. These control peptides were already bound to the respective receptor PDB structures.

Construction and Evaluation of the Multi-epitope Subunit Vaccine

The highest-scoring epitopes for all cell types were used to construct the final vaccine constructs against CCHFV. The epitopes were joined together by using a different kinds of linkers. Intra CTL epitopes were connected by GGGS linker, whereas GPGPG and KK linkers were used to link HTL and B-cell epitopes, respectively. Different adjuvants including β-adjuvant (accession number: AF295370), HABA protein (accession number: AGV15514.1), 50 S ribosomal protein L7/L12 (accession number: P9WHE3), and TLR4 agonist (APPHALS) was attached at the N terminal end of each vaccine, followed by a pan HLA-DR epitope (PADRE) sequence joined by an EAAAK linker (Wu et al. 2010). The antigenicity, allergenicity, solubility and toxicity of the constructed vaccines were evaluated through VaxiJen 2.0, AllerTop, Protein-Sol (https://protein-sol.manchester.ac.uk/) (Hebditch et al. 2017), and ToxinPred server (Gupta et al. 2013), respectively.

Structure Prediction, Refinement, and Validation of the Vaccine Construct

The secondary structure of the vaccine was constructed using PRISPRED online tool (http://bioinf.cs.ucl.ac.uk/psipred/) (Jones 1999). It predicts the α-helix, β-sheet, and coil structure of the vaccine. Protein structure prediction server RaptorX (http://raptorx.uchicago.edu/ContactMap/) was used to predict the tertiary structure of the vaccine construct (Xu 2019). Further, GalaxyRefine (http://galaxy.seoklab.org/) online tool was used to refine the predicted 3D structure of the vaccine construct (Heo et al. 2013). The accuracy of the refined vaccine model was validated by an online server, PROCHECK v6.0 (https://servicesn.mbi.ucla.edu/PROCHECK/) (Laskowski et al. 1993). The ERRAT server was used to predict non bonded atomic interaction patterns (Colovos and Yeates 1993). Moreover, the protein structure analysis program of the ProSA-web (https://prosa.services.came.sbg.ac.at/prosa.php) was used to calculate a quality score for the input protein by comparing it with known structures obtained from X-ray analysis, NMR spectroscopy, and theoretical calculations (Wiederstein and Sippl 2007).

Conformational B-Cell Epitope Prediction and Disulfide Engineering

The ElliPro tool of the IEDB server (http://tools.iedb.org/ellipro/) was used to examine the conformational B-cell epitopes present in the final vaccine construct utilizing default parameters (minimum score: 0.5; maximum distance: 6 Å) [57]. Disulfide by Design 2 v2.13 tool was used for the disulfide engineering of the vaccine candidates, which predicts residues that can form a novel bond after mutating to cysteines [39]. The Chi3 (χ3 ) angle was kept between + 97° ± 30° and − 87° ± 30°, and residue pairs with energy less than 2.2 kcal/mol were selected for mutation. Further, the Cα-Cβ-Sγ angle was set as the default value of 114.6° ± 10.

Protein-Protein Docking and Molecular Dynamics Simulation

Molecular protein-protein docking was performed to assess the binding affinity between vaccine candidates and human immune receptors, e.g., TLR using ClusPro v2.0 (https://cluspro.bu.edu/) (Kozakov et al. 2017) and PatchDock (https://bioinfo3d.cs.tau.ac.il/PatchDock/php.php) (Schneidman-Duhovny et al. 2005) servers. The refined vaccine construct was used as the ligand, while TLR2 (PDB id: 3A7B), TLR3 (PDB id: 1ZIW), TLR4 (PDB id: 3FXI) were used as receptors for docking analysis. FireDock tool was used to refine the complexes obtained from PatchDock and generate a global energy score of the best solution (Andrusier et al. 2007). Further, the iMODS server (http://imods.chaconlab.org/) was used to conduct normal mode analysis (NMA) in internal coordinates to explain the collactive motion of proteins (López-Blanco et al. 2014). The server predicted deformability, B-factors, eigenvalue, and covariance.

In Silico Immune Simulation

In silico immune simulation was conducted using the online simulation server C-ImmSim to evaluate the immunogenicity of the designed vaccine (https://www.iac.cnr.it/~filippo/c-immsim/index.html) (Rapin et al. 2010). Three dosages of the target CCHFV vaccine were administered at a time interval of 4 weeks. The periods were set at 1, 84, and 168 where each time step is equal to 8 h. The simulation step was set at 1050 and other parameters were fixed at default.

Codon Adaptation and In Silico Cloning of the Designed Vaccine

The codon for the finally selected vaccine protein was optimized using Java Codon Adaptation Tool (JCAT, http://www.jcat.de/) for Escherichia coli strain K12 prokaryotic expression system (Grote et al. 2005). Due to the difference between the native and foreign gene expression, host codon usages, rho-independent transcription terminators and prokaryotic ribosome binding sites were avoided. Codon adaptation index (CAI) value and GC content of the adapted sequence was used to determine the protein expression level (Sharp and Li 1987). Sites of HindIII and BamHI restriction enzymes were added at the N-terminal and C-terminal sites, respectively, to clone the finally adapted nucleotide sequence in the pET-28(+) expression vector. A graphical presentation of the workflow followed in the study is presented in Fig. 1.

Fig. 1

 A schematic presentation of the working procedure and tools used to develop multiepitope-based subunit vaccine against CCHFV

Results

Antigenic Protein Selection and Physicochemical Characterization

The CCHFV strain (Nigeria/IbAr10200/1970) was identified and selected from Uniprot database. Three protein sequences such as EPP, RDRP, and NP were retrieved (Additional File 1). The conservancy of these selected proteins was analyzed in different serotypes of the virus, where RDRP was found to be the most conserved among all of them (Table 1). These sequences were further analyzed using the VaxiJen server to determine the antigenicity at a threshold of 0.4. The score showed that EPP and RDRP had an antigenicity score greater than 0.4 and thus, selected for further analysis. Subcellular localization of EPP, RDRP and NP was predicted to be extracellular, nuclear, and cytoplasmic, respectively (Table 1).

Table 1

Retrieved CCHFV protein sequences from UniProt database along with their UniProt ID, length, antigenicity scores, subcellular localization, and protein BLAST results

Sl no.	Protein name	Uniprot ID	Length	Antigenicity score	MinimumIdentity %	MaximumIdentity %	Subcellular localization	Selection status
1.	EPP	Q8JSZ3	1684	0.5145	85.2	99.82	Extracellular	Selected
2.	RDRP	Q6TQR6	3945	0.4306	96.68	99.92	Nuclear	Selected
3.	NP	P27317	482	0.3210	96.06	100.00	Cytoplasmic	Not selected

The selected antigenic proteins (EPP and RDRP) were analyzed using ExPASy-ProtParam to determine their physical and chemical parameters. EPP had less molecular weight than RDRP, and both had an estimated in vitro half-life of 30 h in mammalian reticulocytes. The GRAVY value represents the peptides’ hydrophobicity, and in this case, both proteins showed a negative value which means their hydrophobic nature (Table S1).

T-Lymphocyte Epitope Prediction and Selection from Target Proteins

CTL are vital players of cellular immune response and can recognize immunogenic antigens on the virus-infected cell surface. The antigen-specific T-cell receptor binds CTL epitope and MHC class-I molecules that form complex on virus-infected cell surfaces (Andersen et al. 2006). Therefore, 9-mer CTL epitopes were predicted for 12 different class-I MHC supertypes using NetCTL 1.2 server at a default threshold score. A total of 524 potential CTL epitopes were selected, including 271 from EPP and 253 from RDRP. Based on the antigenicity, immunogenicity, and toxicity, 184 epitopes were selected. Among them, 83 epitopes were found to be non-allergenic (Table S2) and used for further analysis.

The HTL play a central role in adaptive immune response, including stimulation of B-cell to secret antibody and other T-cell activation (Alberts et al. 2002). Here, a total of 155 HTL epitopes (74 epitopes from EPP and 81 epitopes from RDRP) each of 15 aa long were predicted from the IEDB server using the NN-align method against the complete set of 27 human alleles. These predicted epitopes were analyzed for antigenicity, allergenicity, and toxicity. However, only 60 epitopes were found to be antigenic, non-allergenic, and non-toxic. Among these 60 epitopes, only 7 epitopes were showing the cytokine-inducing ability by the induction of IFN-γ, IL-4, and IL-10 (Table S3).

Prediction of B-Cell Epitopes from Target Proteins

Primarily a total of 160 potential linear B-cell epitopes (EPP- 90, RDRP- 90) were identified from the ABCpred server at a default threshold of 0.51. Among them, 55 epitopes were found to be antigenic, non-allergenic and non-toxic (Table S4). Parker hydrophilicity prediction algorithms, Kolaskar and Tongaonkar Antigenicity scale, Karplus & Schulz flexibility prediction tool, Emini surface accessibility prediction method, and Chou and Fasman beta-turn prediction algorithm was used to measure some essential parameters, including hydrophilicity, antigenic propensity, flexibility, surface accessibility, and beta-turn, respectively (Figs. S1 and S2). Then, DiscoTope 2.0 was used to identify conformational B-cell epitopes from both A and B chains of the targeted EPP and RDRP at -3.700 thresholds and 90% specificity. A total of 52 discontinuous epitopes (25 from EPP protein and 27 from RDRP protein) were selected, (Table S5).

Analysis of Epitope Conservancy and Top-Ranked Epitope Selection

Epitope conservancy across other strains is important as it determines the efficacy and confers broad-spectrum immunity (Bui et al. 2007). The previously selected CTL (83) and HTL (60) epitopes were subjected to conservancy analysis. For vaccine construction, 8 CTL epitopes (Table 2) and 4 HTL epitopes (Table 3) were selected based on their high conservancy and antigenicity score. From previously screened 52 epitopes, we have selected 10 top-ranked B-cell epitopes for final vaccine construction (Table 4). The cross-reactivity of the selected CTL and HTL epitopes with human proteomes and interaction with human MHC class alleles was evaluated (Tables 2 and 3). To predict MHC-I and MHC-II binding alleles of the selected CTL and HTL epitopes, the MHC-I and MHC-II binding prediction tools from the IEDB webserver were used, respectively. Among all the CTL epitopes, YLYIVITLY can bind to a large number of alleles including HLA-B*15:01, HLA-A*30:02, HLA-A*03:01, HLA-A*26:01, HLA-A*01:01, HLA-B*35:01, HLA-A*32:01, HLA-B*57:01, HLA-A*68:01, HLA-B*58:01, HLA-A*11:01, HLA-B*53:01, HLA-A*33:01, HLA-A*02:01, HLA-A*23:01, HLA-B*44:02, HLA-A*24:02, and HLA-B*44:03. The results of MHC-I and MHC-II binding alleles were then used for population coverage analysis.

Table 2

Final selected CTL epitopes for multi-epitope-based vaccine construction with their antigenicity, immunogenicity, toxicity, allergenicity, conservancy, cross-reactivity with human proteome and interacting MHC class-I allele prediction

Protein	Epitope	Position	Super Type	CombinedScore	Antigenicity	Immunogenicity	Toxicity	Allergenicity	Epitope Conservancy (%)	Interaction MHC class I alleles	Cross-reactivity
EPP	KKIEISGLK	1438–1446	B27	0.9605	0.9491	0.13994	Non-toxin	NA	100.00	HLA-A*30:01, HLA-A*03:01, HLA-A*11:01	Zero
	HREVEINVL	411–419	B39	2.5097	0.7096	0.31548	Non-toxin	NA	100.00	HLA-B*08:01, HLA-B*40:01	Zero
	TEAIVCVEL	1231–1239	B44	1.4217	1.1768	0.22941	Non-toxin	NA	100.00	HLA-B*40:01, HLA-B*44:03, HLA-B*44:02	Zero
RDRP	YLYIVITLY	1411–1419	A1	1.7622	0.9588	0.32448	Non-toxin	NA	100.00	HLA-B*15:01, HLA-A*30:02, HLA-A*03:01, HLA-A*26:01, HLA-A*01:01, HLA-B*35:01, HLA-A*32:01, HLA-B*57:01, HLA-A*68:01, HLA-B*58:01, HLA-A*11:01, HLA-B*53:01, HLA-A*33:01, HLA-A*02:01, HLA-A*23:01, HLA-B*44:02, HLA-A*24:02, HLA-B*44:03	Zero
	VTDIVVGAI	902–910	A1	1.4576	1.0628	0.27164	Non-toxin	NA	100.00	HLA-A*01:01, HLA-B*51:01, HLA-A*68:02, HLA-A*02:06, HLA-A*32:01	Zero
	LMSFLNWRV	3263–3271	A2	1.3555	1.1247	0.26439	Non-toxin	NA	100.00	HLA-A*02:01, HLA-A*02:03, HLA-A*02:06, HLA-A*68:02, HLA-A*32:01	Zero
	YHSIAELTM	42–50	B39	2.0951	0.6095	0.22589	Non-toxin	NA	100.00	HLA-B*35:01, HLA-B*08:01, HLA-B*51:01, HLA-A*24:02, HLA-B*53:01, HLA-A*23:01, HLA-B*40:01	Zero
	IQNRITGLY	2775–2783	B62	1.4651	0.6585	0.23824	Non-toxin	NA	100.00	HLA-B*15:01, HLA-A*30:02, HLA-A*01:01, HLA-A*32:01, HLA-A*03:01, HLA-B*57:01, HLA-A*26:01, HLA-B*58:01, HLA-A*30:01, HLA-A*11:01, HLA-B*35:01, HLA-B*44:03, HLA-B*44:02, HLA-B*53:01	Zero

NA non-allergen

Table 3

Predicted final HTL epitopes with their antigenicity, allergenicity, toxicity, cytokine inducing properties, conservancy, cross-reactivity, and interacting MHC class-II alleles

Protein	Epitopes	Position	IC50(nM)	Antigenicity	Allergenicity	Toxicity	IL4pred	IL10pred	IFN-epitope	Conservancy (%)	Interacting MHC class II alleles	Cross-reactivity
EPP	KAFSAMPKTSLCFYI	1508–1522	10	0.5604	NA	Non-toxin	Inducer	Inducer	Positive	100.00	HLA-DRB1*09:01, HLA-DRB1*11:01	Zero
	ILFFMFGWRILFCFK	1610–1624	6.3	1.1537	NA	Non-toxin	Inducer	Inducer	Positive	86.67	HLA-DQA1*01:01/DQB1*05:01, HLA-DRB3*01:01, HLA-DPA1*01:03/DPB1*02:01, HLA-DPA1*03:01/DPB1*04:02	Zero
RDRP	EVLINIRNSLKARSE	1240–1254	13.3	1.1352	NA	Non-toxin	Inducer	Inducer	Positive	100.00	HLA-DRB1*08:02, HLA-DRB3*02:02, HLA-DRB1*03:01, HLA-DRB5*01:01, HLA-DRB1*13:02, HLA-DRB1*15:01, HLA-DRB1*11:01	Zero
	VHFKRNQIPSVENVQ	700–714	5	1.2275	NA	Non-toxin	Inducer	Inducer	Positive	80.00	HLA-DRB3*02:02, HLA-DRB1*09:01	Zero

NA non-allergen

Table 4

Results of B-cell binding linear epitopes prediction with their antigenicity, allergenicity, and toxicity analysis to design the multi-epitope based subunit vaccine construct

Protein	Epitopes	Start Position	Score	Antigenicity	Allergenicity	Toxicity
EPP	CLHKEWPHSRNWRCNP	1180	0.97	0.5677	NA	Non-toxin
	SPPITTGSLPLSETTP	79	0.91	0.9648	NA	Non-toxin
	KIVIDKKNKLNDRCTL	466	0.91	1.0197	NA	Non-toxin
	GVTQHNHASFVNLLNI	1372	0.91	1.0254	NA	Non-toxin
	TSLSIEAPWGAINVQS	1060	0.90	1.1190	NA	Non-toxin
RDRP	EVTIPCFTVYGTFVNS	1588	0.93	0.5888	NA	Non-toxin
	NKLIDHCVDMEKKREA	3861	0.91	0.7054	NA	Non-toxin
	TGQLITHGRVSAKHND	1148	0.91	1.0537	NA	Non-toxin
	LHPEFRDLTPDFSLTQ	683	0.90	1.7810	NA	Non-toxin
	KELINSTGLSDLELES	3715	0.90	0.8086	NA	Non-toxin

NA non-allergen

Population Coverage of the Epitopes

The immunogenic response of epitope-based vaccine may fluctuate due to the high polymorphic nature of human HLA in different ethnicities and regions (Bui et al. 2006). Therefore, the percentage of population coverage of the selected 8 CTL and 4 HTL epitopes and their MHC binding alleles were calculated, both individually and combinedly, applying IEDB provided population coverage analysis tool. The most significant coverage was found in Europe (99.68%) and the least population coverage was obtained for central America (42.81%) (Table 5). The results of the analysis revealed that the predicted epitopes combinedly could cover 98.94% of the world population, thus, expected to work effectively throughout the world. Furthermore, MHC class-I and class-II restricted epitopes could cover 97.69% and 54.13% of the worldwide population, respectively (Table S6).

Table 5

Analysis of population coverage of the selected T-cell epitopes based on their respective HLA alleles

Population/area	Combined (class I and class II)
	Coverage (%)a	Average hitb	PC90c
Central Africa	90.84	3.91	1.14
Central America	42.81	0.64	0.17
East Africa	94.78	4.52	2.14
East Asia	98.35	5.43	2.72
Europe	99.68	6.67	3.97
North Africa	97.60	5.12	2.49
North America	99.33	6.18	3.38
Northeast Asia	96.40	4.79	2.33
Oceania	96.39	4.60	2.32
South Africa	94.10	4.21	2.01
South America	90.10	3.35	1.01
South Asia	96.16	4.69	2.20
Southeast Asia	96.55	4.86	2.35
Southwest Asia	93.03	4.17	1.91
West Africa	97.18	5.17	2.44
West Indies	99.07	6.17	3.29
World	98.94	5.86	3.12

Here MHC classes I and II were combined to calculate the percentage of population coverage in different geographical locations

aProjected population coverage

bAverage number of epitope hits / HLA combinations recognized by the population

cMinimum number of epitope hits / HLA combinations recognized by 90% of the population

Tertiary Structure Prediction and Epitope-HLA Docking Analysis

The 3D structure of the finally selected T-cell epitopes (8 CTL and 4 HTL) was predicted from the PEP-FOLD server and visualized using Discovery Studio (Fig. S3). The PEP-FOLD online peptide structure prediction tool created five 3D models for each epitope and the best model was identified as a ligand for docking analysis. The docking was performed to evaluate the binding capacity of the selected epitopes with class-I and class-II MHC molecules. Here, structure-based molecular docking was conducted to calculate the binding affinity using AutoDock Vina in PyRx. The result showed that all the epitopes bound efficiently with the binding groove of the respective alleles (Fig. 2). For MHC class-I, the epitope ‘TEAIVCVEL’ exhibited the highest binding affinity of − 8.3 kcal/mol with HLA A*01:01 and the epitope ‘VHFKRNQIPSVENVQ’ strongly interacted with MHC class-II molecule HLA- DRB1*15:01, with a binding affinity of − 6.9 kcal/mol (Table S7). Moreover, the binding affinity of the control peptides was − 8.4 kcal/mol and − 7.1 kcal/mol for MHC class-I and MHC class-II molecules, respectively. BIOVIA Discovery Studio Visualizer was used to identify the interacting residues and binding pattern visualization.

Fig. 2

Molecular docking of finally selected CTL and HTL epitopes. Molecular docking of finally selected CTL and HTL epitopes with HLA A*01:01 and HLA- DRB1*15:01 allele, respectively. Here, A control peptide (NP44), B KKIEISGLK, C HREVEINVL, D TEAIVCVEL, E YLYIVITLY, F VTDIVVGAI, G LMSFLNWRV, H YHSIAELTM, I IQNRITGLY represent docking of CTL epitopes. Similarly, J control peptide (human myelin basic protein) and K KAFSAMPKTSLCFYI, L ILFFMFGWRILFCFK, M EVLINIRNSLKARSE, N VHFKRNQIPSVENVQ represent docking of HTL epitopes. Hydrogen bond interactions are highlighted as green dotted lines

Vaccine Construction Considers Various Parameters Including Antigenicity, Allergenicity, Toxicity, and Solubility Prediction

A total of 8 CTL, 4 HTL, and 10 B-cell epitopes were selected from the EPP and RDRP of the target strain to construct a multi-epitope subunit vaccine against CCHFV. The epitopes were linked together using GGGS, GPGPG and KK linkers, respectively. In this study, we primarily constructed four vaccines using different adjuvants and screened them based on antigenicity, allergenicity, toxicity, and solubility to obtain an effective vaccine (Table S8). Based on these criteria, V1 was selected as a potential vaccine for further analysis, and it contained 45 aa long adjuvant β-defensin at the N terminal end, followed by a PADRE sequence. β-defensin enhances vaccine immunogenicity by stimulating the activation of toll-like receptors and the PADRE sequence improve vaccine potency (Funderburg et al. 2007). The results revealed that our vaccine construct is highly antigenic (0.7724) and non-allergenic in nature. The scaled solubility value was found to be 0.530. This result showed that the vaccine had a solubility score greater than the average threshold value of 0.45, which indicates the high solubility potential of the vaccine. Additionally, the vaccine was also found to be non-toxic. The finally designed multi-epitope vaccine consisted of 427 aa (Fig. 3A) and the position of each epitope was listed in Table S9. ProtParam tool was used to determine several physiochemical properties of the designed vaccine. The molecular weight of the vaccine was 46.425 kDa. The theoretical pI value was calculated as 9.83. The estimated half-life was 30 h in mammalian reticulocytes in vitro, > 20 h in yeast in vivo, and > 10 h in E. coli in vivo. The instability index of the vaccine was found to be 30.30 (< 40), which classifies the protein as stable. The value of the aliphatic index and GRAVY was calculated as 81.26 and − 0.316, respectively. The GRAVY value suggests that the vaccine is hydrophilic and the aliphatic index indicates the thermostability of the vaccine. The results of these above analyses showed our vaccine construct as a potential candidate.

Fig. 3

Multi-epitope based subunit vaccine against CCHFV.  A Graphical presentation of the vaccine construct. The position of each epitope was listed in Table S9. B Secondary structure prediction of the vaccine protein where the alpha-helix, beta-sheet, and coil structure are indicated by pink, yellow, and ash colour, respectively. C The vaccine consists of 427 aa, where an adjuvant is linked at the N-terminal of sequence with the aid of EAAAK linker followed by PADRE sequence. GGGS linker joined CTL epitopes, GPGPG linker connected HTL epitopes and KK linker joined B-cell epitopes. (D) The tertiary structure of the vaccine after refinement showing α-helix (sailor blue), beta strands (crimson), and loops (aqua)

PRISPRED online tool was used to predict the secondary structure of the vaccine construct. The results showed that 89 (20.85%) residues were involved in α-helix creation, and 113 (26.46%) residues were involved in β-strands formation. For coil structure formation, 225 aa were involved, consisting of 52.69% of the vaccine sequence (Fig. 3B). The overall peptide sequence of the vaccine construct is presented in Fig. 3C. To predict the tertiary structure, RaptorX server was used. The structure was constructed using distance-based protein folding method as no templet resembles our designed vaccine. Further, the structure was refined using GalaxyRefine server to improve the quality of the designed vaccine protein which generates five models. Among them, model 1 was selected as the best model (Fig. 3D) as it showed a GDT-HA score of 0.9906, RMSD 0.282, MolProbity 1.910, Clash score 14.9, and Poor rotamers 0.0 (Table S10). Further, the Ramachondon plot of the improved structure revealed that 89% of the residues are in favoured regions (represented by red colour), and 11% of the residues are in the additional allowed regions (represented by yellow colour) (Fig. S4A). ProSA-web and ERRAT were also used to validate the improved vaccine protein model. The results of the quality analysis revealed that the ERRAT score was 73.538 and the Z score of the model quality was − 7.14 (Fig. S4B), indicating a good quality of the vaccine construct.

Screening of Conformational B-Cell Epitopes and Vaccine Disulfide Engineering

B-cells play a critical role in the adaptive humoral immune system and produce antigen-specific immunoglobulin when B cell receptors bind specific pathogen-derived antigens. We found five conformational B-cell epitopes with scores varying from 0.546 to 0.808 in the vaccine, consisting of 247 residues (Fig. 4A–E and Table S11).

Fig. 4

Conformational B-cell epitopes of the vaccine construct. A–E The predicted five epitopes are indicated by faux chrome yellow colour and the rest of the residues are in aqua colour. F Disulfide engineering of the designed vaccine protein. Here only two mutated pairs are compatible for disulfide bond formation, shown in blue and crimson colour

The structural protein model was subjected to Disulfide by Design online tool to identify the pairs that can be mutated to cysteines to form disulfide bonds. In our study, a total of 34 pairs of residues have been identified with the potential to form a disulfide bond in disulfide engineering (Table S12). However, only 2 pairs of residues i.e., 136 PHE-149 ILE and 278 LEU-281 THR were selected for disulfide bond formation considering the energy score and χ3 angle (Fig. 4F).

Molecular Docking Analysis Between the Vaccine and Immune Receptors

Molecular docking of the constructed vaccine with human immune receptors TLR2, TLR3, and TLR4 was conducted using two online tools, i.e., ClusPro 2.0 and PatchDock. FireDock server was used to analyze results obtained from PatchDock and to calculate the global energy of the docked complexes. A total of thirty docked complexes denoted as clusters were generated by ClusPro for each vaccine-TLR docking. Among them, the cluster with the lowest binding energy was selected. The lowest binding energy was observed − 1242.7 for TLR3 and the lowest global energy was − 16.34 for TLR4 (Table S13). Moreover, the binding energy were − 1138.7 and − 1140.0 for TLR2 and TLR4, respectively. Molecular protein-protein docking of the vaccine constructs and TLR3 immune receptor was selected for further analysis (Fig. 5A). Normal mode analysis is a popular approach for demonstrating large scale protein mobility and stability. The iMODS server was used for NMA for the internal coordinates of the constructed vaccine protein and TLR3 complexes. The eigenvalue of the complex was found 8.246566e−06 (Fig. 5B). This value is associated with motion stiffness of each normal mode and directly related to require energy for structure deformation. The chain hinges indicate the region of deformability of the complex (Fig. 5C). The B-factor was obtained from NMA and docked complex, which is resembling the RMS (Fig. 5D). The covariance matrix of the coupling pair residues where correlated, uncorrelated, and anti-correlated motions are represented by red, white, and blue colour (Fig. 5E). The eigenvalue is also inversely correlated to the variance, which is associated with each normal mode (Fig. 5F). An elastic network model is also obtained from the dynamic analysis, which indicates the pair of atoms are connected by spring. One spring is represented by each dot and coloured according to the degree of stiffness (Fig. 5G).

Fig. 5

Molecular dynamics simulation of the designed vaccine-TLR3 complex. A Docked complex of vaccine construct (aqua) and TLR3 receptor (green) in cartoon format, B eigenvalue, C deformability of the complex (chain hinges indicates the region of deformability), D B-factor, E covariance map (here, correlated, uncorrelated, and anti-correlated motions are represented by red, white, and blue colour), F variance (individual and cumulative variances are indicated by red and green colour, respectively), G elastic network (darker grey colour represents mone stiffer regions)

Immune Response Simulation

The simulation study revealed that the designed vaccine could induce higher secondary and tertiary responses than the primary response. The higher concentration of IgM was characterized as a primary immune response. The secondary and tertiary immune response produced elevated immunoglobulins such as IgG1 + IgG2, IgM, and IgG + IgM with concurrent antigen reduction (Fig. 6A). The increased memory B cell was evident (Fig. 6B) that ensure antigen clearance upon subsequent exposures and the activated B cell counts supported long lasting effectivity (Fig. 6C). In addition, a clear increase in both HTL and CTL were observed with respective memory development (Fig. 6D–F). Besides, elevated macrophage populations (Fig. 6G) as well as higher antigen-presenting dendritic cells (Fig. 6H) were observed after vaccine administration. Further, higher IFN-γ and IL-2 levels (Fig. 6I) also support the antiviral potency of the target vaccine after immunization. Overall, the results suggest that the currently designed vaccine might effectively induce immune response and immunity against the virus in question.

Fig. 6

In silico immune simulation of the final vaccine construct. A Immunoglobulin production after vaccine administration, B B cell population, C B cell population per entity-state, D T-helper cell population after vaccination, E T-helper cell population per entity-state, F cytotoxic T-cell population per state, G macrophages population after vaccine administration, H dendritic cell population per state, I concentration of cytokines and interleukins. The insert plot shows IL-2 level with the Simpson index, D (a measure of diversity). The smaller the D value, the lower the diversity. An increase in D over time indicates the emergence of different epitope-specific dominant clones of T-cells. TC cytotoxic T-cell, TH Helper T-cell, DC dendritic cell, and MA macrophage

Codon adaptation increases the translation efficacy of any gene in a heterogenous host. The codon of the constructed vaccine was optimized for the E. coli K12 expression system.

. The reverse transcription of the vaccine protein sequence generated an optimized DNA sequence (Additional File 2) of 1281 nucleotides where the CAI value was found 0.992 and the GC content was calculated as 48.712%. The optimized range of CAI value and GC content lies between 0.8 and 1.0 and 30–70%, respectively (Shey et al. 2019). Therefore, our constructed vaccine reveals a high possibility of protein expression in the E. coli host. In the following step, restriction sites were attached to both ends of the adapted codon sequences for restriction enzymes HindIII (173) and BamHI (198) (Fig. 7A). Finally, the restriction cloning module of the SnapGene software was utilized to insert the DNA sequence into pET28a (+) vector (Fig. 7B). The length of the newly cloned plasmid was 6635 bp.

Fig. 7

In silico restriction cloning of the multi-epitope vaccine sequence into pET28a(+) expression vector. A Restriction strategy for the cloning of newly designed subunit vaccine into the vector pET28a(+) using HindIII and BamHI restriction enzyme sites. B The adapted vaccine sequence is indicated by red colour while keeping the plasmid backbone in black colour

Discussion

Vaccine is one the most effective medical intervention that remarkably prevents infectious diseases and improves human health. The conventionally designed vaccine has an indubitable contribution to global health but has some limitations. The development process is based mainly on inactivated or live-attenuated pathogen, inactivated toxins, and the recombinant subunits that are costly and take years to accomplish (Finco and Rappuoli 2014). Moreover, traditional vaccine development is complex against pathogens with antigenic hypervariability due to its inadequate immune response. The improper attenuation of the pathogens may also result in adverse immune reactions. Besides, a few conventional vaccines include whole organisms or large macromolecules that may contain unnecessary antigenic load and induce unwanted protective immune responses or allergic reactions (Seib et al. 2012). However, in modern times, the advancement in genome sequencing and bioinformatic techniques have created the opportunity for researchers to overcome such difficulties by following different approaches, including fusion of immunoinformatics and reverse vaccinomics for vaccine development (María et al. 2017). In recent years, many alternative strategies have been applied to design competent vaccines to overcome the drawbacks of traditional vaccines. Peptide-based vaccine development is one of the critical applications of immunoinformatics that relies on incorporating different antigenic parts of pathogens, as a result, induce specific immune responses. However, successful multi-epitope vaccine construction is still challenging because of immunodominant epitope selection, specific delivery system requirements, and adjuvant (Li et al. 2014). Therefore, the selection of appropriate epitopes of target antigen is crucial for multi-epitope driven vaccine designing. Researchers around the world have used immunoinformatics and reverse vaccinology approaches to develop multi-epitope based vaccines against the Hepatitis C virus (Ikram et al. 2018), Dengue virus (Ali et al. 2017), Ebola virus (Dash et al. 2017), Oropouche virus (Adhikari et al. 2018), and Lassa virus (Sayed et al. 2020).

This study was designed to construct an effective vaccine against CCHFV. With the wide-spreading tick vector of the virus and the high fatality rate, the virus is marked as one of the most dangerous viruses in the world (Mertens et al. 2013). A few studies were conducted previously using structural Gc, Gn, RDRP, and NP of CCHFV individually and combinedly to design epitope-based vaccines incorporating both B cell and T cell epitopes (Oany et al. 2015; Nosrati et al. 2019; Tahir Ul Qamar et al. 2021; Khan et al. 2021b). A multi-peptide vaccine candidate was designed against CCHFV by combining several CD8+, CD4+, and B-cell epitopes with appropriate linker and Mycobacterium tuberculosis lipoprotein LprG as an adjuvant that were predicted to cover worldwide 98% population (Khan et al. 2021b). Similarly, another novel multi-epitope CCHFV vaccine was developed by taking 9 immunodominant epitopes from four proteins to raise the worldwide coverage to 99.74% (Tahir Ul Qamar et al. 2021). EPP is present on the viral surface and can attach to the host cell receptor (Bergeron et al. 2007). The other protein, RDRP is responsible for the replication and transcription (Kinsella et al. 2004). On the other hand, NP failed to elicit adequate antigenicity, therefore, excluded from further analysis. However, a recombinant candidate vaccine called Modified Vaccinia virus Ankara, expressing the CCHFV NP generated humoral immune responses in susceptible mice model but failed to protect the host from the virus (Dowall et al. 2016). A functional multiple-subunit vaccine should include epitopes of both B and T cells, therefore inducing adequate immune responses. CTL, also known as CD8+ cells, kill infected cells that express surface antigen as a marker by inducing apoptosis (Lanzavecchia 1998). Like CTLs, HTLs (CD4+ cells) also play a key role in various immunologic processes. HTLs can activate cytotoxic T-cells, B-cells, and cells of the innate immune system by secreting different types of cytokines such as IFN-γ, IL-4, and IL-10 (Bevington et al. 2017). Cytokine production is crucial in defining the effector functions and development of T cell types. For instance, IFN-γ is secreted by TH1 cell, promotes further differentiation and activates macrophages as well as induces B cells isotype switching to opsonizing and complement-binding IgG subclasses (Ahlers et al. 2001). TH2 cell produces IL-4 and IL-10. IL-4 enhances further differentiation and proliferation of TH2 as well as stimulates B cells to produce IgE antibodies that bind to mast cells (Walker and McKenzie 2018). In contrast, IL-10 downregulates host innate and cell-mediated immune responses, particularly macrophage activation and TH1 mediated reactions (Fiorentino et al. 1991). Therefore, the HTL epitopes were evaluated for their specific cytokine inducing ability (Table 3). B-cells play a central role in the humoral adaptive immune system by producing antibodies.

The finally selected vaccine was constructed by linking the top-ranked epitopes sequentially with GGGS, GPGPG, and KK linkers (Fig. 3A). However, multi-epitope based vaccines are poorly immunogenic when administered without any additional immunogenic molecules. Therefore, a 45 aa long adjuvant β-defensin was added at the N-terminal end, followed by a pan HLA-DR epitope (PADRE) sequence joined with EAAAK linker. The adjuvant β-defensin have immunomodulatory and antimicrobial properties and the PADRE sequence has antigen-specific CD4+ T cell enhancing properties. Besides, the linkers were also used by many researchers to enhance the expression, stability, and folding of the vaccine protein (Sayed et al. 2020; Hasan et al. 2020).

Another important property of the multi-epitope based vaccine is the population coverage. The highly polymorphic nature of MHC class I and II retain the response of CTL and HTL epitopes in different ethnic communities. MHC molecules are also known as HLA which are cell surface receptor protein that binds and displays pathogen-derived peptide fragments for recognition by T-cells. The population coverage of a vaccine depends on the number of HLA alleles they bind. Therefore, to predict the worldwide distribution of the alleles, we have selected the CTL and HTL epitopes with their respective HLA alleles. Our selected epitopes combinedly could cover 98.94% of the world population (Table 5), with maximum coverage in Europe (99.68%) and North America (99.33%). CCHFV cases most significantly happened in Sub-Saharan Africa, Southern Europe, and Central and Southeast Asia. Therefore, vaccine candidates are vital to protect individuals against CCHFV infection in these geographical regions. The potential vaccine could cover over 95% population in the African region and 96.16% in South Asia (Table 5).

Antigenicity, allergenicity, and toxicity analysis revealed the vaccine construct as a potent and effective vaccine. The result of solubility analysis indicated the high solubility potential of the vaccine protein (Table S8). Immunoinformatics tools also showed the physiochemical properties of the vaccine protein. The overall properties are consistent with already reported multi-epitope subunit vaccines (Tahir ul Qamar et al. 2020; Khan et al. 2021a). The refined tertiary structure was validated with Ramachondon plot analysis and the results showed that most residues were presented in favoured and allowed regions. The Z score (− 7.14) and ERRAT score (73.538) revealed the good quality of the overall structure of the refined vaccine (Fig. S4).

The immune simulation study was conducted to identify the immunogenic potency of the vaccine construct against CCHFV. The obtained results showed consistency with typical immune responses. The antigen enhances immune responses following administration at a different time interval. Development of memory B cells and T cells were observed along with several months lasting memory B cells. In a previous study, in vivo administration of a CTL defined antigen containing lentivector resulted in strong CD8+ T cell responses which is relevant to our simulation results (Esslinger et al. 2003). After administering the vaccine, the higher concentration of IFN-γ and IL-2 indicated elevated levels of TH cells and efficient immunoglobulin production supporting humoral immune response. Similarly, an elevated level of antibody titers was found in both healthy young and elderly subjects four weeks after influenza vaccination, where cytokines such as IFN-γ, IL-10, and IL-6 were found significantly higher in young than elderly individuals (Bernstein et al. 1998). In addition, sufficient activity of macrophage and dendritic cells were observed in the present study (Fig. 6). The Simpson index also depicts a possibility of diverse immune responses based on the clonal specificity.

The codon optimization was carried out to attain a high level of vaccine protein expression in the E. coli K12 expression system (Additional File 2). Generally, foreign gene expression may vary inside the host due to inconsistency of the mRNA codon. Therefore, codon optimization is crucial to enhance translational efficacy. The GC content (48.712%) and CAI value (0.992) were fortunately within the optimum limit, promising for the high-level expression of a protein in the E. coli K12 strain. Moreover, we have engineered the vaccine disulfide bond to improve protein thermostability. In this study, novel immunoinformatics tools were combined to design a novel multi-epitope based subunit vaccine against CCHFV. The designed vaccine consists of B-cell and T-cell epitopes could efficiently induce both innate and adaptive immune systems of the host.

Conclusions

CCHFV is one of the deadliest RNA viruses that are responsible for causing hemorrhagic manifestation in vertebrates, including human. In this study, a possible multi-epitope based subunit vaccine was developed through reverse vaccinology and immunoinformatic techniques against CCHFV. Different analytical experiments identified both potential B-cell and T-cell epitopes from the selected protein sequences to produce effective immune responses that could be used to construct the vaccine candidates. The designed vaccine has possessed all the criteria of an ideal vaccine in terms of its ability to generate cell-mediated and humoral immune responses, and structural orientation. Molecular docking analysis confirmed a higher binding strength between human immune receptors and vaccine candidate. Furthermore, molecular dynamic and immune simulation study demonstrated the stability, large scale mobility, and immune responses in real event. The proposed vaccine will help to develop a novel vaccine against CCHFV as therapeutics and preventive measures, ultimately reducing the mortality and morbidity caused by it. Therefore, further in vitro and in vivo investigations using appropriate model organisms are needed to ensure the true potentiality of the designed vaccine construct.

Below is the link to the electronic supplementary material.

Supplementary material 1 (PDF 6023.0 kb)