Computational design and characterization of a multiepitope vaccine against carbapenemase-producing Klebsiella pneumoniae strains, derived from antigens identified through reverse vaccinology.

Klebsiella pneumoniae is a Gram-negative pathogen of clinical relevance, which can provoke serious urinary and blood infections and pneumonia. This bacterium is a major public health threat due to its resistance to several antibiotic classes. Using a reverse vaccinology approach, 7 potential antigens were identified, of which 4 were present in most of the sequences of Italian carbapenem-resistant K. pneumoniae clinical isolates. Bioinformatics tools demonstrated the antigenic potential of these bacterial proteins and allowed for the identification of T and B cell epitopes. This led to a rational design and in silico characterization of a multiepitope vaccine against carbapenem-resistant K. pneumoniae strains. As adjuvant, the mycobacterial heparin-binding hemagglutinin adhesin (HBHA), which is a Toll-like receptor 4 (TLR-4) agonist, was included, to increase the immunogenicity of the construct. The multiepitope vaccine candidate was analyzed by bioinformatics tools to assess its antigenicity, solubility, allergenicity, toxicity, physical and chemical parameters, and secondary and tertiary structures. Molecular docking binding energies to TLR-2 and TLR-4, two important innate immunity receptors involved in the immune response against K. pneumoniae infections, and molecular dynamics simulations of such complexes supported active interactions. A codon optimized multiepitope sequence cloning strategy is proposed, for production of recombinant vaccine in classical bacterial vectors. Finally, a 3 dose-immunization simulation with the multiepitope construct induced both cellular and humoral immune responses. These results suggest that this multiepitope construct has potential as a vaccination strategy against carbapenem-resistant K. pneumoniae and deserves further validation.

Introduction

Antimicrobial resistance is a current public health major threat and World Health Organization (WHO) foresees that, by 2050, 10 million deaths will occur, due to the constant rise of multidrug resistant pathogens, especially in clinical settings [1]. Among them, Klebsiella pneumoniae possesses an extremely developed plastic capability of acquiring resistance to many antibiotics, even to those considered last resort, such as carbapenems [2]. Carbapenem-resistant K. pneumoniae strains have been reported worldwide and, in 2014, they were responsible for>2 million critical infections, mostly in clinical settings [3].

While a vaccination strategy against this bacterial pathogen and its multidrug resistant strains is highly desired, in particular for patients that require medium to long hospitalization stays, there is no vaccine available against K. pneumoniae to date [4]. Currently, the only vaccine against K. pneumoniae in clinical trials is a tetravalent strategy based on LPS O-antigens and adjuvanted with the squalene-based AS03, named Kleb4V, which entered Phase I/II in 2021 [5]. In the past, formulations using inactivated K. pneumoniae have shown some efficacy, but required daily administration during several months, rendering these immunization strategies impractical from a logistical point of view [elegantly reviewed in 4]. Additionally, these preparations may lead to undesired side effects, due to the presence of toxins such as lipopolysaccharides (LPS), capsular polysaccharides (CPS), among other proteins, which may cause exacerbated inflammation [2]. Furthermore, multidrug resistant K. pneumoniae clinical isolates belong to many, diverse sequence types (ST), especially ST258 and ST512, among others [2], [6], [7], and it is unlikely that whole, inactivated bacterial preparations based on a single strain induce a potent immunological response against the many relevant clinical ST isolates currently circulating in endemic areas. Other immunization strategies reported were based in purified CPS and LPS O-antigens, in some cases conjugated to proteins to increase their immunogenicity, in outer membrane vesicles (OMV) and protein-based approaches [4]. However, an efficient immunization approach which can protect against multidrug resistant K. pneumoniae strains, particularly carbapenem-resistant clinical isolates is still lacking and is urgently needed, especially for inpatients who require hospital stays and are under immunosuppression drugs regimen, such as transplantation recipients.

In the last decades, technological developments, following the availability of sequenced bacterial genomes, led to an innovative approach to identify novel, potential antigens, called reverse vaccinology. This approach is based on in silico analysis to screen which proteins within a pathogen’s translated genome have the potential of possessing antigenic properties and be included in vaccination strategies. Reverse vaccinology was a crucial technique employed in the development of the Neisseria meningitidis serogroup B vaccine Bexsero [8], and is being currently applied to several other pathogens [9], [10], [11], [12], [13], [14], [15]. Consequently, many tools were developed to help characterize antigens identified by reverse vaccinology, either in terms of properties or at structural level [8], [9], [14], [16]. Although these bioinformatics techniques allow screening of thousands of potential antigens and selection of the few most promising for further costly animal studies, reverse vaccinology and associated tools also have limitations. Bacterial protective antigens containing polysaccharides, which often constitute protective antigens, are often not predicted [17]. When considering pathogens that require a T cell response, such as malaria, tuberculosis and hand, foot and mouth disease (HFMD), epitopes prediction software can be combined with reverse vaccinology tools, but will require expensive wet lab validation and are time-consuming [13], [15], [17]. In addition, there are many layers of antigenic complexity that need to be addressed such as accessibility of surface antigens, antigen masking, protein degradation, posttranslational modifications, exon–intron structures, alternative splicing, especially for eukaryotic antigens [8]. In the case of viral antigens, reverse vaccinology still needs to address how to predict accurately protein complexes and expression levels [8].

In this work, novel, potential K. pneumoniae antigens with a high predominance in carbapenem-resistant clinical isolates were identified and characterized in silico, using a reverse vaccinology approach. Using several available bioinformatics tools (Table S1) [8], [9], [10], [11], [12], [13], [14], [15], [16], [17], a multiepitope vaccine containing T and B cell epitopes of the antigens was designed and characterized for its antigenic potential. We also assessed the molecular docking between our multiepitope vaccine and TLR-2 and TLR-4, due to their importance in K. pneumoniae infection [2], [10], [14], [18], [19]. A cloning strategy to produce the suggested multiepitope approach is suggested, and an immunization simulation was also performed and analyzed.

Material and methods

Reverse vaccinology

To identify potential K. pneumoniae antigens, 8 complete chromosomal genomes of strains of this bacterial pathogen available on the National Center for Biotechnology Center (NCBI) Genome website were downloaded. Genome HS11286 was chosen as reference, and genomes were selected from strains belonging to clusters not in the immediate vicinity of each other (Fig. S1 and Table S2). In parallel, bacterial proteins predicted subcellular localizations for the same strains were downloaded from PSORTdb 4.0 (https://db.psort.org/) [20], and the genomes were annotated using RAST 2.0 online software [21], using as reference the genome of the strain HS11286. Potential antigens were chosen based on their protein identity among the 8 K. pneumoniae genome strains selected (higher than 80 %), peptide chain length (higher than 250 amino acids), low protein identity with E. coli K-12 (lower than 50 %), bacteriophages’ sequences excluded, subcellular location (outer membrane or secreted proteins) and 2 or less internal helices (to less affect negatively proteins’ solubility). Short-listed antigens were blasted against all genomes available at NCBI [22] (https://blast.ncbi.nlm.nih.gov/Blast.cgi?PROGRAM = blastp&PAGE_TYPE = BlastSearch&LINK_LOC = blasthome), excluding klebsiellal genomes, to ensure no homologous sequences with human proteins or other species were present, to avoid generation of self-immune or cross-species, unspecific epitopes. The presence of potential antigens that met the referred criteria was assessed in 467 sequenced genomes of carbapenem-resistant K. pneumoniae clinical isolates from 2 health structures, IRCCS-ISMETT (Istituto Mediterraneo per i Trapianti e Terapie ad Alta Specializzazione, Palermo, Italy) and Azienda Ospedaliera Universitaria Careggi (Florence, Italy). All 467 strains genomes displayed the genetic sequence of the carbapenemase KPC-2 (accession number QBQ66786), and were fully described and characterized elsewhere [6], [7], [23], [24], [25], [26], [27], [28]. The antigens present in>95 % of clinical isolates were considered for further analysis. This procedure is schematized in Fig. S2 (Part I).

Prediction of physicochemical properties, antigenicity, allergenicity, toxicity and solubility of potential antigens

Following identification of potential K. pneumoniae antigens, their physicochemical characteristics were predicted by Expasy Protparam (https://web.expasy.org/protparam/), such as amino acid composition, theoretical isoelectric point, half-life (both in vitro and in vivo), instability and aliphatic indexes, molecular weight (in kDa) and grand average of hydropathicity index (GRAVY) [29], [30].

Before assessing the antigenic potential of the predicted proteins, their sequences were submitted to AntigenDB (http://crdd.osdd.net/raghava/antigendb/), an antigen database [31], to confirm if the predicted klebsiellal antigens were identified previously and submitted to this database. AntigenDB is linked to several important protein databases, and allows for the extraction, analysis, and validation of antigens. Following confirmation that these potential antigens were not included in AntigenDB, their antigenic potential was then evaluated using 2 available software online, ANTIGENpro (https://scratch.proteomics.ics.uci.edu/) [32] and VaxiJen 2.0 (https://www.ddg-pharmfac.net/vaxijen/VaxiJen/VaxiJen.html) [33]. Allergenic potential of the predicted antigens was determined using the online available tools AllerTOP 2.0 (https://www.ddg-pharmfac.net/AllerTOP/) [34] and AllergenFP 1.0 (https://ddg-pharmfac.net/AllergenFP/index.html) [35].

The hypothetical solubility of the predicted klebsiellal antigens was determined using Protein–Sol server (https://protein-sol.manchester.ac.uk) [36], and the virulence potential of the antigens was assessed through VirulentPred [37] (http://bioinfo.icgeb.res.in/virulent/index.html).

Following the design of a multiepitope strategy against carbapenem-resistant K. pneumoniae strains based on the antigens identified in this work, the above tools were also used to characterize the proposed multiepitope construct.

A flow diagram illustrating the process used in the antigen characterization can be found online on Fig. S2 (Part II).

Cloning and expression of the predicted klebsiellal antigens

The selected antigens were cloned in the expression vector pET-9d (Merck Life Science S.r.l., Milano (MI), Italy). Their sequences (Accession numbers in Table 1) were codon optimized for expression in E. coli K-12, using the online software Java Codon Adaptation Tool (JCat) server (https://www.jcat.de/) (Fig. S2, Part III) [38]. All codon optimizations showed a codon optimization index (CAI) superior to 0.85 and a GC content between 50 % and 55 %, which is considered favorable for foreign protein expression [38]. To allow for the purification and to improve the solubility of the selected antigens, a double tag, 6xhistidin and maltose-binding protein (6His-MBP), was introduced at the N-terminal, as described elsewhere [39]. The optimized sequences were synthetized by an external supplier (Eurofins Genomics Italy, Vimodrone (MI), Italy), using the restriction sites NcoI and BamHI for insertion in pET-9d. Cloning and expression of the recombinant proteins was done as previously detailed [40], [41]. All cloning and expression reagents were purchased from New England Biolabs Inc. (Ipswich, MA, USA) and Fisher Italia (Rodano (MI), Italy), respectively.

Table 1

Carbapenem-resistant K. pneumoniae antigens identified by reverse vaccinology.

Antigen ID	Location	Accession number	Protein name	Presence in clinical isolates (%)
Klebs#1	Outer membrane	YP_005225656	putative TonB-dependent siderophore receptor	96.44
Klebs#2	Extracellular/ secreted	YP_005228444	β-1,4-mannanase	99.55
Klebs#3	Extracellular	WP_002916123	type 3 fimbria adhesin subunit MrkD	99.78
Klebs#4	Outer membrane	YP_005229348	putative cellulose synthase/ BcsC protein	95.55

Epitope prediction

Human helper T lymphocytes (CD4+) 15-mer epitopes were predicted using the online available software NetMHCIIpan-3.2 (https://www.cbs.dtu.dk/services/NetMHCIIpan-3.2/), and chosen from the classical 7-allele reference set [42], [43]. The T cell helper epitopes present in the klebsiellal antigens identified by reverse vaccinology were also analyzed for their ability to induce IL-4 secretion, which skews the immune response to a Th2-like response, desired for antibody production against extracellular pathogens such as K. pneumoniae, through the online tool IL4pred (https://webs.iiitd.edu.in/raghava/il4pred/index.php) [44]. IL-17 inducing epitopes were also predicted using the online tool IL17eScan (https://metabiosys.iiserb.ac.in/IL17eScan/) [45].

B cell epitopes 16-mers were predicted with the online software ABCpred (https://webs.iiitd.edu.in/raghava/abcpred/index.html) [46]. Only epitopes with a score higher than 0.9 were considered for further evaluation.

The schematics of epitope prediction software used, based on the ME_Klebs amino acido sequence, can be observed on Fig. S2 (Part IV).

Multiepitope vaccine construct design and characterization

A multiepitope vaccine construct was designed, using the predicted CD4+ and B cell epitopes with higher scores, and named ME_Klebs. To increase the immunogenicity of the construct, the mycobacterial heparin-binding hemagglutinin protein (HBHA, accession number AAC26052), a known TLR-4 agonist [47], was added at the N terminal of the construct. Epitopes’ sequences were connected using the linker GPGPG among them and the linker EAAAK between HBHA and the first epitope, as suggested by previous works [10], [13], [14]. The resulting ME_Klebs was characterized as described above for the predicted antigens (Fig. S2, part II). Toxicity of the predicted epitopes was evaluated using the ToxinPred online tool (https://webs.iiitd.edu.in/raghava/toxinpred/multi_submit.php) [48].

ME_Klebs secondary structure prediction

To predict and characterize the secondary structure of ME_Klebs, 2 online tools were used, PSIPRED (https://bioinf.cs.ucl.ac.uk/psipred/) [49] and RaptorX Property (https://raptorx.uchicago.edu/StructurePropertyPred/predict/) [50] (Fig. S2, part V). PSIPRED allows the analysis of the secondary structure and provides information on the protein folding, transmembrane helix packing and topology, domain and function [49]. RaptorX, on the other hand, provides predicted secondary structure α-helices, β-sheets and coils, disorder regions and solvent accessibility [50].

ME_Klebs three-dimensional structure prediction, refinement, and validation

To determine the three-dimensional (3D) structure of ME_Klebs, the online software I-TASSER (Iterative Threading ASSEmbly Refinement) tool was used (https://zhanggroup.org/I-TASSER/). This tool is defined as an in silico hierarchical procedure which predicts 3D protein structures and, based in the latter, their function. Such assemblies are compared with similar structures previously deposited in the Protein Data Bank (PDB; https://www.rcsb.org/) and ranked in terms of template modeling score (TM-score), in which values higher than 0.5 are predicted as topologically accurate [51]. In parallel, ME_Klebs tertiary structure was also estimated with AlphaFold Colab, a simplified version of AlphaFold v2.1.0 (https://colab.research.google.com/github/deepmind/alphafold/blob/main/notebooks/AlphaFold.ipynb#scrollTo=2tTeTTsLKPjB) [52], an artificial intelligence tool that predicts with high accuracy the structure of proteins based on their genetic sequence and through comparison to similar proteins.

The obtained ME_Klebs 3D structures were subsequently refined, i.e., further analyzed with online tools to obtain a structure which would resemble more the native, real conformation of the protein, by predicting, identifying and correcting errors derived from 3D structure predictions bioinformatics tools [53]. Refinement was done using the GalaxyRefine online tool (https://galaxy.seoklab.org/cgi-bin/submit.cgi?type=REFINE), which initially restructures and rearranges amino acid side chains, leading to a general structure relaxation. GalaxyRefine uses molecular dynamics (MD) simulations which originate improved structure quality, at both local and global levels [54].

Following refinement, it is advised that 3D protein structures be submitted to further validation, allowing identification and correction of errors not detected in the earlier simulations [10], [13]. Validation of ME_Klebs 3D structure was done through ProSA (Protein Structural Analysis, https://prosa.services.came.sbg.ac.at/prosa.php) [55], [56] and trough ERRAT (https://saves.mbi.ucla.edu/) [57]. In addition, a Ramachandran plot for ME_Klebs was obtained through the VADAR (Volume, Area, Dihedral Angle Reporter) online tool (http://vadar.wishartlab.com/). This type of plot allows a visual depiction of how the amino acids conformations behave in the principal peptide core chain and whether they fit with the preferred, more stable conformations region or not [58], [59].

The bioinformatics tools used in the prediction, refinement and validation of the tertiary structure and respective scheme can be observed in Fig. S2 (Part V).

Conformational B cell epitopes prediction

Continuous B cell epitopes, as described above, are predicted based solely in the amino acid sequence of a protein. However, some epitopes can also be conformational, or discontinuous, meaning they depend on the 3D structure of the antigenic protein and not on its linear sequence [60]. Once ME_Klebs 3D structure was obtained, refined and validated, its conformational B cell epitopes were determined using the online tool ElliPro (based on the words Ellipsoid and Protrusion; http://tools.iedb.org/ellipro/) (Fig. S2, part V).

Molecular docking of ME_Klebs with TLR-2 and with TLR-4

As mentioned earlier, TLR-2 and TLR-4 are innate immunity receptors which play an important role in host protection during K. pneumoniae infection [2], [10], [14], [18], [19]. Thus, the potential interaction between ME_Klebs and these two receptors was also simulated using online tools. Both TLR-2 and TLR-4 structures were obtained from PDB (ID 2Z7X and ID 3FXI, respectively). First, ME_Klebs, TLR-2 and TLR-4 were submitted to CPORT (Consensus Prediction Of interface Residues in Transient complexes), to identify active residues in their structures (https://alcazar.science.uu.nl/services/CPORT/). This initial step precedes molecular docking experiments and it is advised in conditions where experimental information is not available [61], as in the specific case of ME_Klebs. The obtained structures were then submitted to the web tool HADDOCK 2.4 (High Ambiguity Driven protein–protein DOCKing), which predicts the docking between 2 biomolecular entities, based on previously available data from either predicted or identified protein interactions in ambiguous interaction restraints (AIRs) (https://bianca.science.uu.nl/haddock2.4/) [62]. In parallel, both dockings were also predicted with ClusPro (https://cluspro.bu.edu/home.php) [63] and with PatchDock (https://bioinfo3d.cs.tau.ac.il/PatchDock/php.php) [64], [65].

Dockings of ME_Klebs with either TLR-2 or TLR-4 predicted with either HADDOCK 2.4 or ClusPro were subsequently refined with the online tool FireDock (https://bioinfo3d.cs.tau.ac.il/FireDock/) [66], [67] and PRODIGY (PROtein binDIng enerGY prediction; https://we nmr.science.uu.nl/prodigy/) [68].

Finally, both refined ME_Klebs-TLR-2 and ME_Klebs-TLR-4 dockings were further improved by submitting them to HawkDock web tool (http://cadd.zju.edu.cn/hawkdock/), a multifunctional platform which unifies several software, including the molecular mechanics/generalized born surface area (MM/GBSA) free energy decomposition score, for prediction of protein–protein interactions [69]. The docking prediction with the lowest MM/GBSA score is considered the one which resembles most the native protein–protein interaction [69].

The process used for the molecular docking of ME_Klebs with TLR-2 and TLR-4 is schematized in Fig. S2 (part VI).

Molecular dynamics prediction of ME_Klebs dockings with either TLR-2 or TLR-4

All MD simulations of the TLR-2-ME_Klebs and TLR-4-ME_Klebs complexes were performed using GROMACS 2021.5 [70], [71], on an NVIDIA TESLA Volta 100 graphical processor unit (GPU) installed on the high performance computing (HPC) server at Bezmialem Vakif University (Istanbul, Turkey). Simulations were done as described elsewhere [72], [73], [74], [75], [76], [77], [78]. Three independent MD runs were performed for 50 ns, and the coordinates were recorded every 10.0 ps. The trajectories of the 2 simulated complexes were analyzed for both root mean square deviation (RMSD) of the structure over the course of the simulation and root mean squared fluctuation (RMSF) of residues over time, using embedded GROMACS packages [79]. The average RMSD and RMSF of the 3 independent runs and respective graphs’ plotting were done using Microsoft 365 Excel (Redmond, WA, USA). The Proteins Interaction Calculator webserver was used to generate the interaction report between ME_Klebs and the 2 innate immune receptors, TLR-2 and TLR-4 [80]. The free binding energy of the complexes, molecular mechanics/Poisson–Boltzmann surface area (MM/PBSA), was calculated by use of the g_mmpbsa GROMACS package, as described elsewhere [81], and the reported values are the average of 3 independent runs.

The molecular dynamics processes follow the results obtained previously for the molecular docking predictions and can be observed in Fig. S2 (Part VI).

In silico codon optimization and cloning of ME_Klebs

ME_Klebs genetic sequence was codon optimized with JCat, for expression in E. coli K12 and without the restriction sites NcoI and BamHI, as described above. In silico cloning was done by inserting the ME_Klebs codon optimized sequence in the expression vector pET-9d (Merck) sequence, using the aforementioned restriction sites sequences in the 5′ and 3′ ends of the multiepitope sequence, respectively, using the SnapGene® software (from Insightful Science; available at https://www.snapgene.com). The codon optimization and cloning schemes can be observed on Fig. S2 (Part VII).

In silico immunization simulation with ME_Klebs

To simulate an immune response induced by ME_Klebs, the C-ImmSim web tool was used (https://kraken.iac.rm.cnr.it/C-IMMSIM/). It allows the prediction of both cellular and humoral immune responses against an antigen, simulating dynamics in agreement with experimental observed basic immunological phenomena in mammalian hosts [82]. Simulated vaccination with ME_Klebs consisted in three doses given one month apart each (days 0, 28 and 56), and volume, steps of simulation and number of antigens to inject were 50, 1,000 and 1,000, respectively, with no adjuvant (lipopolysaccharides, LPS). Host HLA types used were A0101 and B0702 (MHC I) and DRB1_0101 (MHC II). The in silico immunization scheme is included in the overall flow diagram together with the cloning process, in Fig. S2 (Part VII).

Results and discussion

Identification of carbapenem-resistant K. Pneumoniae antigens by reverse vaccinology and subsequent characterization

Through reverse vaccinology, 7 potential K. pneumoniae antigens were identified, of which 4 were present in the vast majority of sequenced genomes of carbapenem-resistant K. pneumoniae Italian clinical isolates (Table 1). Putative TonB-dependent siderophore receptor, henceforth called Klebs#1, is an outer membrane protein which has an auxiliary role in iron sequestration, a process required for K. pneumoniae survival in the mammal host [2], and has been identified as a potential antigen by previous in silico analyses studies [14], [83], [84]. β-1,4-mannanase, hereafter referred as Klebs#2, is a glycoside hydrolase which catalyzes hydrolysis of mannan (linear, galacto-, gluco- and galactogluco-), in the extracellular milieu of several microorganisms’ cell wall [85]. Type 3 fimbria adhesin subunit MrkD has been widely studied by its virulence potential, due to its role in biofilm formation in antibiotic-resistant K. pneumoniae strains [7], [23], [86], [87] and it is widely disseminated in K. pneumoniae clinical isolates [7], [88], [89], [90]. Not surprisingly, this virulence factor, from now on designated as Klebs#3, has been used in immunization strategies in mouse models and showed partial protection both in active and passive immunizations against intranasal challenge with a pathogenic K. pneumoniae strain [91]. Several immunoinformatic works have analyzed and characterized Klebs#3 as a potential antigen [84], [92], [93]. Cellulose biosynthesis protein BcsC, henceforward described as Klebs#4, is a protein crucial to maximize bacterial cellulose synthesis, and downregulation of the bcsC gene severely impacted the ability of K. pneumoniae to form viable biofilms [94]. To the best of our knowledge, neither Klebs#2 nor Klebs#4 have been identified by in silico approaches or considered potential antigens by previous works.

Physicochemical properties of the four selected antigens can be observed in Table 2, such as molecular weight (MW, in kDa), theoretical isoelectric point (pI) and half-live in mammalian, yeast and bacterial cells. All proteins, except for Klebs#4, show an instability index lower than 40, suggesting they are stable [29]. The aliphatic index was higher than 64 for all antigens, indicating they are thermostable between 20 °C and 45 °C [29], [30]. Hydrophilicity of the identified antigens was determined by their GRAVY index, which were all negative, showing that all antigens are, in theory, hydrophilic [29].

Table 2

Physicochemical properties of reverse vaccinology selected klebsiellal antigens and multiepitope vaccine construct.

	theoretical isoelectric point (pI)	in vitro half-life	in vivo half-life	Instability index	aliphatic index	molecular weight (MW)	grand average of hydropathicity (GRAVY)
Klebs#1	5.54	>20 h (yeast, in vivo)>10 h (E. coli, in vivo)	30 h (mammalian reticulocytes, in vitro)	26.38	73.61	79.25	−0.403
Klebs#2	5.71	>20 h (yeast, in vivo)>10 h (E. coli, in vivo)	30 h (mammalian reticulocytes, in vitro)	29.47	71.14	81.19	−0.450
Klebs#3	9.28	>20 h (yeast, in vivo)>10 h (E. coli, in vivo)	30 h (mammalian reticulocytes, in vitro)	32.58	81.87	35.06	−0.009
Klebs#4	6.16	>20 h (yeast, in vivo)>10 h (E. coli, in vivo)	30 h (mammalian reticulocytes, in vitro)	41.82	73.19	145.10	−0.477
ME_Klebs	8.77	>20 h (yeast, in vivo)>10 h (E. coli, in vivo)	30 h (mammalian reticulocytes, in vitro)	31.02	65.46	55.49	−0.636

The Klebs#1–4 antigens were also submitted to prediction software, to understand their virulence and antigenicity potential, solubility and allergenicity (Table 3). Although none of these antigens were present in the bacterial antigens database AntigenDB (data not shown), they were all considered potential antigens by both ANTIGENpro and VaxiJen 2.0 online tools. All selected klebsiellal proteins have an antigenicity score higher than 0.4, indicating they are indeed potentially antigens. In terms of solubility, only Klebs#3 possesses a solubility index superior to 0.45, which makes it the only identified antigen with a higher solubility possibility. On the other hand, Klebs#3 is also the only of these antigens to be classified as an allergen, hindering its potential immunogenicity. In terms of virulence, all tested antigens but Klebs#2 are considered virulent.

Table 3

Prediction of vaccinology desired properties for the four antigens identified by reverse vaccinology and the multiepitope construct.

	Virulence	Antigenicity		Solubility	Allergenicity	Toxicitya
		ANTIGENpro	VaxiJen
Klebs#1	Yes	0.85	0.67	0.35	Non-allergen	ND
Klebs#2	No	0.97	0.42	0.32	Non-allergen	ND
Klebs#3	Yes	0.90	0.83	0.50	Allergen	ND
Klebs#4	Yes	0.97	0.56	0.35	Non-allergen	ND
ME_Klebs	Yes	0.88	0.76	0.49	Non-allergen	Non-toxin

, Toxicity only assessed for the epitopes present in ME_Klebs; ND, not determined.

Preliminary wet-lab expression experiments were performed with the Klebs#1–4 potential antigens (Fig. S2, part III). It was observed that, even with the MBP tag, all 4 proteins were found mostly in the insoluble fraction, including Klebs#3. Interestingly, although this protein has been predicted as soluble (Table 3), other works fused it with glutathione S-transferase (GST) [91], [95], which is used, as MBP, to improve protein solubility [96], suggesting the referred authors dealt with the same issue of Klebs#3 solubility and that the in silico prediction for this particular antigen differs from the empiric evidence obtained.

Although the 4 predicted klebsiellal antigens show a high antigenicity potential index as estimated by both ANTIGENpro and VaxiJen 2.0, other factors render them less appealing to be used in a vaccine against carbapenem-resistant K. pneumoniae strains. In terms of virulence, it is known that many antigens are pathogenic virulence factors, due to their external localization and to the fact that a vaccine against these factors can prevent their virulent phenotype and, consequently, the pathological manifestations of the disease caused by the pathogens [37]. In agreement, all studied antigens in this work showed virulence potential, except for Klebs#2. Ideally, recombinant proteins with biotherapeutical potential applications are soluble, to avoid colloidal instability and complex and costly purification protocols with aggressive agents, such as urea, and protein correct refolding [36], [97]. Our experimental approach confirmed the predicted insolubility for most of the antigens identified by reverse vaccinology and further confirmed previous results indicating the lack of solubility of Klebs#3 [91], [95]. Moreover, Klebs#3 was predicted to be an allergen, which should not be included in vaccine preparations for safety issues. Overall, none of the 4 antigens under study possesses all characteristics which would classify them at theoretical level as ideal antigens.

Design and characterization of a multiepitope vaccine against carbapenem-resistant K. Pneumoniae strains

Prediction of human helper CD4+ T cell, B cell, IL-4 inducing and IL-17 inducing epitopes

Since none of the klebsiellal proteins identified by reverse vaccinology fulfilled all theoretical, predicted attributes to be fully considered a potential antigen, a multiepitope approach was subsequently pursued, using epitopes present in the proteins from Table 1. Similar approaches have been done for K. pneumoniae and other pathogens, such as Mycobacterium tuberculosis and HFMD virus [10], [11], [12], [13], [14], [15]. K. pneumoniae is an extracellular pathogen, for which a humoral response and Th-17 cells, but not a cytotoxic T cell action, are required to prevent further dissemination of this bacterium within the host [2]. For this rational multiepitope vaccine design, CD8+ T cell epitopes were not considered since these immune cells have no major action against this bacterium. K. pneumoniae is an extracellular pathogen, and an immunological response against this type of pathogens relies more on a humoral response and not on cytotoxic T cells, a subset of CD8+ T cells which kill infected cells [98]. Instead, immunity against K. pneumoniae requires B cells and helper CD4+ T cells. Also, it has been seen that both interferon-γ (INF-γ) and cytotoxic T cells only have a limited action on K. pneumoniae initial invasion, but not during subsequent infection [2]. Moreover, volunteers immunized with Dentavax, a therapeutical and immunoprophylactic preparation consisting of several killed bacterial species, including K. pneumoniae, increased the levels of CD57+CD8+ T cells [99], a cellular subpopulation associated with short-term cardiovascular mortality in acute myocardial infarction patients due to their pro-inflammatory phenotype and high cytotoxicity potential [100]. Thus, not only CD8+ T cells seem to contribute scarcely against K. pneumoniae infection, but they could also, in certain conditions, induce inflammation-related pathologies, which render them undesirable as a component of an immune response induced by a vaccination strategy against this bacterial pathogen. For these reasons, no CD8+ T cell or INF-γ epitopes were selected as components of ME_Klebs. It is assumed that the selected T cell and B cell epitopes, if applied together in a synergistic way, can mount a Th-2 response, which will induce the production of opsonophagocytic antibodies against carbapenem-resistant K. pneumoniae strains and clear its infection and dissemination or, at least, prevent the pathologies associated with this bacterium. Thus, CD4+ T cell, IL-4 inducing, IL-17 inducing, and B cell epitopes present in the sequences of the 4 K. pneumoniae antigens derived from reverse vaccinology were predicted (Table 4, Table 5). For each protein, the 2 helper 15mer T cell epitopes which showed both low percentile scores (which means the highest binding affinity probability) [101], and support vector machine (SVM) score for predicted IL-4 inducers, were chosen, covering 5 different HLA-DR loci: HLA-DRB1*03:01, HLA-DRB3*01:01, HLA-DRB3*02:02, HLA-DRB4*01:01 and HLA-DRB5*01:01 (Table 4). Interestingly, no IL-17 inducing epitopes were predicted in any of the Klebs#1–4 antigens’ sequences. Similarly, 2 16mer B cell epitopes from each klebsiellal antigen were selected, with a binding score higher than 0.9 (Table 5), corresponding to high binding affinity [102].

Table 4

Human helper T-cell epitopes selected for multiepitope construct.

Antigen ID	Epitope sequence	Allele	Percentile score	IL4 inducer	SVMscore
Klebs#1	HQKVVVRNYSNATGL	HLA-DRB3*02:02	0.19	Yes	0.34
	TVGRLDTNRIIDAIS	HLA-DRB3*01:01	1.40	Yes	0.30
Klebs#2	DGTYKKINPFALTDP	HLA-DRB5*01:01	1.30	Yes	0.11
	PWSGQFLLNVTFNGK	HLA-DRB3*02:02	0.38	Yes	0.14
Klebs#3	SPTVMLDMVVGRVVV	HLA-DRB1*03:01	0.44	Yes	0.29
	ETYLSANAITVVSPS	HLA-DRB3*02:02	0.41	Yes	0.66
Klebs#4	SAKLFRADVLRHNKD	HLA-DRB3*01:01	0.62	Yes	0.11
	KLIRAMQSDPQNTDL	HLA-DRB4*01:01	0.58	Yes	0.15

Table 5

B-cell epitopes selected for multiepitope construct.

Antigen ID	Epitope Sequence	Allele	Percentile Score
Klebs#1	WSGIDDTGTYLFEGDP	HLA-DRB3*02:02	0.19
	HREGEGAIDNDKRRTT	HLA-DRB3*02:02	0.32
Klebs#2	TVVGRDISDGVNEWNP	HLA-DRB5*01:01	1.30
	LHRIEDDARGPCRADP	HLA-DRB3*02:02	0.38
Klebs#3	TVNIVYPDVFSSRVYN	HLA-DRB1*03:01	0.44
	DWTMSAPGGASYRCTS	HLA-DRB3*02:02	0.41
Klebs#4	AGTASGDAWRRYGANP	HLA-DRB3*01:01	0.62
	RRSGQARGQGYANLNS	HLA-DRB4*01:01	0.58

Design and characterization of a multiepitope vaccine against carbapenem-resistant K. Pneumoniae isolates

Once the T cell and B cell epitopes were chosen, the multiepitope vaccine, ME_Klebs, was designed. As a molecular adjuvant, the M. tuberculosis HBHA protein was included in the chimeric construct, at its N-terminal. This protein is a known TLR-4 agonist and has been used as adjuvant in cancer vaccines [47]. Moreover, HBHA was also chosen considering that TLR-4 can prevent mortality and K. pneumoniae growth in a mouse model [2]. To link the epitopes, the linker GPGPG was used. This linker has been widely used for induction of helper T cells both in multiepitope and DNA vaccines [103] and in several other works to connect B cell epitopes [10], [13], [14]. The linker EAAAK was used to join HBHA and the first epitope. EAAAK is a rigid linker which promotes protein expression and pH and thermal stability [104], thus being chosen in other multiepitope vaccine design works as well [10], [13], [14]. A 6His tag was added at the end, to facilitate purification of the recombinant chimeric construct [39]. The schematic representation of the proposed construct ME_Klebs and its sequence can be observed in Fig. 1 and S3, respectively.

Fig. 1

Schematic representation of ME_Klebs. Epitopes (ep) order follow the cardinal order of the Klebs#1–4 antigens by their position in Table 4, Table 5, and each color represents a particular antigen (Klebs#1, blue; Klebs#2, purple; Klebs#3, green; and Klebs#4, red). B cell epitopes were placed before CD4+ T cell epitopes. Linkers and 6His tag are placed as described in the main text and indicated in the figure. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Once ME_Klebs sequence was defined, it was analyzed like the original Klebs#1–4 antigens in terms of physicochemical properties, virulence, antigenicity, solubility and allergenicity (Table 2, Table 3). Also, the epitopes that were incorporated into ME_Klebs were also analyzed in terms of toxicity prediction (Table 3). This multiepitope construct has a MW of 55.49 kDa and a theoretical pI of 8.77 and has comparable half-lives in vitro and in vivo to the Klebs#1–4 antigens. Its instability index is lower than 40 and the aliphatic index higher than 64, thus being predicted as a stable protein within a temperature range of 20 °C and 45 °C, which indicates that, within the human host following immunization, ME_Klebs would be thermostable. Its GRAVY was negative, suggesting this protein is indeed hydrophilic (Table 2). As expected, due to its chimeric nature, ME_Klebs was also not present in AntigenDB (data not shown). Both antigenicity prediction software show a high antigenic potential for ME_Klebs. Furthermore, it displays a solubility index of 0.49, strengthening that, following expression, ME_Klebs is a soluble protein. Finally, ME_Klebs is not predicted to be neither an allergen nor possessing toxic-prone epitopes in its sequence (Table 3). Overall, all prediction tools results indicate that ME_Klebs fulfils all requirements to be an ideal antigen and, consequently, a good immunization strategy for carbapenem-resistant K. pneumoniae isolates.

Secondary ME_Klebs structure prediction

In terms of predicted secondary structure, ME_Klebs possesses 30 % α-helix, 11 % β-strand and 58 % coil (Fig. 2). In terms of exposition of amino acids to solvent, 61 % were predicted to be fully exposed, 19 % partially exposed and 19 % fully hidden from solvent. This prediction suggests that ME_Klebs shows potential as a structural antigen, due to its high number of α-helices and coil-coiled domains. In fact, the presence of α-helical coil-coiled domains in epitope vaccines is important, since it allows correct protein folding, mimicking the native proteins structures, and thus confer an adequate and functional humoral immunity against the target pathogen. Also, those antibodies would not only target epitopes present in solvent-exposed surfaces of single α-helices, but also epitopes in other α-helices in the vicinity [105].

Fig. 2

ME_Klebs secondary structure PSIPRED prediction, indicating it possesses 30% α-helix, 11% β-strand and 58% coil-coiled regions.

3D ME_Klebs structure prediction, refinement and validation

Five potential 3D ME_Klebs structures were predicted by I-TASSER online tool. The second option was chosen, since it had the highest Confidence score (C-score), equal to −2.34 (Fig. 3A). Usually, C-scores fall with the range (−5, 2) and the highest value is correlated with prediction confidence. Template modelling score (TM-score) and RMSD were 0.44 ± 0.14 and 13.2 ± 4.1 Å, respectively. In bioinformatics, TM-score measures how similar 2 protein structures are, and values around 0.5 or higher are indicators of correct protein topology, while the backbone RMSD measures the average distance between the corresponding backbone atoms of superimposed proteins [51]. The chosen ME_Klebs structure TM-score value suggests a prediction with topology of high accuracy.

Fig. 3

ME_Klebs 3D structure prediction, refinement and validation. A) ME_Klebs 3D structure obtained by I-TASSER and B) refined with GalaxyRefine. Superimposed refined structure (in colors) over the pre-refined structure (in gray). C) Ramachandran plot, in which red areas correspond to protein core, yellow to favorable areas and green to acceptable zones, while grey areas are disallowed areas. Squares represent residues in core and allowed regions, crosses amino acids in disallowed regions and triangles are glycine residues. Of the ME_Klebs residues, 79% are in the protein core, 14% in the allowed areas and 4% in the disallowed regions. D) Location of ME_Klebs (black circle) in a ProSA analysis graphic, indicating that the Z-score of this multiepitope construct is −2.33. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Interestingly, prediction of the ME_Klebs tertiary structure with AlphaFold Colab gave different results (Fig. S4A). The structure obtained is different from the one obtained with I-TASSER (Fig. 3A), and the predicted local-distance difference test (pLDDT) values are mostly under 70 (Fig. S4B), especially for the multiepitope sequence region, which classifies this structure model as having low confidence [52]. On the other hand, the predicted aligned error (PAE), as indicated by the residues’ pair (x, y) in the output graph produced by AlphaFold (Fig. S4C), is low for the ME_Klebs tertiary structure, which indicates the residues positions are well defined within the structure [52].

The ME_Klebs structures were subsequently refined by the GalaxyRefine online tool (Fig. 3B and Fig. S4D). Since both tertiary structure prediction tools (I-TASSER and AlphaFold Colab) refined structures showed similar values regarding the refined parameters (data not shown), only the structures derived from I-TASSER were further considered. Overall, I-TASSER predicted tertiary structures seem to be more reliable than those inferred by AlphaFold Colab, despite the high potential of the deep learning algorithms and high performance of the latter. However, the AlphaFold method has proven to have limitations in predicting 3D structures for proteins for which experimental data of related proteins is not available [106], [107]. In agreement, AlphaFold Colab multiple sequence alignment (MSA) did not predict many unique sequences for the epitope polypeptide chain of ME_Klebs (Fig. S4E), which is also the protein sequence with low confidence model prediction (Fig. S4B). Thus, this prediction tool does not seem to be ideal to characterize multiepitope vaccine constructs, due to their chimeric nature, unless precedent structure data about them or related constructs is available in online databases.

From the initial I-TASSER model 2 prediction, 5 different refined model structures were provided by GalaxyRefine, and model number 4 was chosen for further studies. This model displayed the best characteristics among all 5 models, with a global distance test - high accuracy (GDT-HA) of 0.8861, an RMSD of 0.668, MolProbity of 2.599 Å, a clashscore of 24.5, poor rotamers of 1.0 and Ramachandran favored value of 81.5. A GDT-HA score varies between 0 and 1, respectively the least and the most accurate possible. As default, the original structure GDT-HA value is 1, meaning the refined ME_Klebs structure does not deviate greatly from the initial model. MolProbity and associated criteria (clashscore, rotamers and Ramachandran) perform an contact analysis of all atoms, adding and optimizing all hydrogen atoms and then estimate their H-bond, steric clash, and favorable van der Waals forces [108]. Ramachandran value of 81.5 was an improvement from the original structure, which was of 62.6, meaning that at least 81.5 of the residues are in the most conformational stables zones. In agreement, the Ramachandran plot for ME_Klebs showed 79 % of the residues in the preferred zones, 14 % in the allowed areas and 4 % of the amino acids were outliers (Fig. 3C). ProSA analyses showed a Z-score of −2.33, visibly distant of the average Z-score range for proteins containing the same residues number (Fig. 3D). This could indicate that the ME_Klebs structure contains errors [55], [56]. However, when submitted to ERRAT, which compares the target structure to other protein structures derived from experimental crystallographic studies [57], ME_Klebs displayed an overall quality factor of 72.51. Above 50, the model is considered of high quality [57]. These results confirmed the high quality of the refined ME_Klebs 3D structure.

Conformational B cell epitopes present in ME_Klebs structure prediction

In terms of predicted conformational B cell epitopes, when submitting the validated ME_Klebs 3D structure to ElliPro, a total of 7 discontinuous epitopes were estimated, with a score higher than 0.5 (score range between 0.508 and 0.848) (Table 6). The predicted conformational B cell epitopes A-G varied in length (75, 92, 52, 4, 21, 9 and 4 residues) and in residues’ composition (Fig. 4 and Table 6). Individual scores for each predicted discontinuous B cell epitope, from A-G, can be observed in Fig. 5. These results corroborate the presence of 7 conformational B cell epitopes in the ME_Klebs structure with a high prediction score.

Table 6

Predicted ElliPro conformational B cell epitopes residues of ME_Klebs.

	Residues	Number of Residues	Score
A	A:M1, A:A2, A:E3, A:N4, A:S5, A:N6, A:I7, A:D8, A:D9, A:I10, A:K11, A:A12, A:P13, A:L14, A:L15, A:A16, A:A17, A:L18, A:G19, A:A20, A:A21, A:D22, A:L23, A:A24, A:L25, A:A26, A:T27, A:V28, A:N29, A:E30, A:L31, A:I32, A:T33, A:N34, A:L35, A:R36, A:R38, A:A39, A:E40, A:E41, A:T42, A:R43, A:T44, A:D45, A:T46, A:R47, A:S48, A:R49, A:V50, A:E51, A:E52, A:S53, A:R54, A:A55, A:R56, A:L57, A:T58, A:K59, A:L60, A:Q61, A:E62, A:D63, A:L64, A:P65, A:E66, A:Q67, A:L68, A:T69, A:E70, A:L71, A:K74, A:A77, A:E78, A:E79, A:K82	75	0.848
B	A:L397, A:L439, A:L440, A:N441, A:V442, A:T443, A:F444, A:N445, A:G446, A:K447, A:G448, A:P449, A:G450, A:P451, A:G452, A:S453, A:P454, A:T455, A:V456, A:M457, A:L458, A:D459, A:M460, A:V461, A:V462, A:G463, A:R464, A:V465, A:V466, A:V467, A:G468, A:P469, A:G470, A:P471, A:G472, A:E473, A:T474, A:Y475, A:L476, A:S477, A:A478, A:N479, A:A480, A:I481, A:T482, A:V483, A:V484, A:S485, A:P486, A:S487, A:G488, A:P489, A:G490, A:P491, A:G492, A:S493, A:A494, A:K495, A:L496, A:F497, A:R498, A:A499, A:D500, A:V501, A:L502, A:R503, A:H504, A:N505, A:K506, A:D507, A:G508, A:P509, A:G510, A:P511, A:G512, A:K513, A:L514, A:I515, A:R516, A:A517, A:M518, A:Q519, A:S520, A:D521, A:P522, A:N524, A:H528, A:H529, A:H530, A:H531, A:H532, A:H533	92	0.761
C	A:E89, A:A90, A:T92, A:S93, A:R94, A:Y95, A:N96, A:E97, A:L98, A:V99, A:E100, A:R101, A:G102, A:E103, A:A104, A:A105, A:L106, A:E107, A:L109, A:R110, A:A152, A:K153, A:L154, A:V155, A:G156, A:I157, A:E158, A:L159, A:P160, A:A166, A:K167, A:K168, A:A169, A:A170, A:P171, A:A172, A:K173, A:K174, A:A175, A:A176, A:P177, A:A178, A:K179, A:K180, A:A181, A:A182, A:A183, A:K184, A:A186, A:P187, A:A188, A:K189	52	0.67
D	A:R72, A:E73, A:F75, A:P121	4	0.559
E	A:A331, A:G332, A:T333, A:A334, A:S335, A:G336, A:D337, A:A338, A:W339, A:R340, A:Y342, A:A344, A:D398, A:T399, A:P431, A:G432, A:P433, A:S435, A:G436, A:Q437, A:F438	21	0.525
F	A:V302, A:G305, A:P306, A:G307, A:P308, A:D310, A:M313, A:S314, A:A315	9	0.523
G	A:A122, A:E123, A:G124, A:Y125	4	0.508

Fig. 4

Conformational B cell epitopes predicted in the ME_Klebs 3D structure. Each field corresponds to one discontinuous B cell epitope (A-G), in agreement with the results presented in Table 6.

Fig. 5

2-dimensional, individual ElliPro score chart, residue-by-residue, for each conformational B-cell epitope present in the ME_Klebs 3D structure.

Molecular docking of ME_Klebs to TLR-2 and TLR-4

In a pneumonia mouse model, both TLR-2 and TLR-4 are related to reduction of mortality and K. pneumoniae dissemination within the host [2], [19]. At an initial stage of infection, TLR-4 controls and impairs bacterial spreading while TLR-2 reduces the levels of inflammation caused by the infection but, at a later stage, only TLR-2 has also a role in preventing bacterial expansion [19]. Thus, as referred above, a vaccine against this pathogenic bacterium should activate both TLR-2 and TLR-4. Docking of ME_Klebs with either TLR-2 or TLR-4 was performed with different online tools, such as ClusPro, PatchDock and HawkDock, to ensure the most precise prediction was selected for further studies, as suggested elsewhere [13]. Docking complexes 3D structures calculated by ClusPro were submitted to PRODIGY, to estimate their binding affinity score. Similarly, the global energy of the docked ME_Klebs complexes with either TLR-2 or TLR-4 predicted by PatchDock was assessed by FireDock (Fig. 6). HawkDock ranks its docking complexes prediction also based in binding free energy, concomitantly with the lowest MM/GBSA scores obtained. Following ClusPro and PRODIGY analyses, the best ME_Klebs-TLR-2 and ME_Klebs-TLR-4 predicted 3D complex structures had a binding affinity energy of −44.23 and −32.85 kcal/mol, respectively. The global binding energy, provided by the PatchDock/FireDock best predictions, was of −41.14 and −37.54 kcal/mol, for ME_Klebs-TLR-2 and ME_Klebs-TLR-4, respectively. Finally, the relative binding free energy, as estimated by HawkDock, for the ME_Klebs-TLR-2 and ME_Klebs-TLR-4 complexes, was of −52.40 and −44.83 kcal/mol, respectively, suggesting these structures are the most thermodynamical favorable at the conditions used during prediction, usually 25 °C and 1 atm. These predictions indicate that ME_Klebs could activate both innate and adaptive immunological responses through stimulation of both TLR-2 and TLR-4, although only empirical data would confirm these assumptions.

Fig. 6

Docking complexes of ME_Klebs with TLR-2 (A,B) and TLR-4 (C,D), as predicted by FireDock. On green, ME_Klebs, on light blue TLR-2 and TLR-4. The residues that mediate interaction between the two proteins are highlighted on red. A general view of the docking is seen for ME_Klebs-TLR-2 (A) and ME_Klebs-TLR-4 (C). An enlarged view of the interaction site between the two proteins can also be observed both for ME_Klebs-TLR-2 (B) and ME_Klebs-TLR-4 (D). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Molecular dynamics simulation of ME_Klebs-TLR-2 and ME_Klebs-TLR-4 complexes

The structural changes between the docked and the simulated complexes at 50 ns can be observed in Fig. 7A, B. The ME_Klebs-TLR-4 complex exhibited a stable structure starting at 11 ns, while ME_Klebs-TLR-2 achieved stability after 33 ns. The RMSD of the simulation is usually indicative of the stability of the systems over the course of the simulation. Variations of less than 1 nm are highly stable, while variation values between 1 and 2 nm are considered stable. The RMSD of ME_Klebs-TLR-2 reached its initial convergence at 33 ns, with an RMSD of approximately 2.57 nm (Fig. 7C). In contrast, the ME_Klebs-TLR-4 complex reached a convergence towards a stable structure much earlier, at 11 ns, with an RMSD smaller than 1.3 nm (Fig. 7D).

Fig. 7

Molecular dynamics of the complexes of ME_Klebs with TLR-2 and TLR-4. Structural changes observed between the molecular docked (green) and the simulated complexes (red), at 50 ns, of ME_Klebs with TLR-2 (A) and TLR-4 (B). The average RMSD of complexes and ligands after 3 independent runs is shown, for ME_Klebs-TLR-2 (C) and ME_Klebs-TLR-4 (D). The comparison of average residue flexibility, RSMF, during the simulations, for both complexes, is also indicated (E). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Amino acid residues 1–200 of the ME_Klebs were not directly interacting with the TLR2 receptor, and displayed more fluctuation, as they were not restricted by any interactions (Fig. 7A, B, E). The residues D255, G256, W260, I271, E272, D273, R353, W434, S435, V465, V467, G468, T474, A478, and V483 that interacted with the TLR2 receptor exhibited restricted flexibility as compared to the first 200 amino acids of the ME_Klebs sequence. In contrast, several of the first 1 – 200 amino acids of the ME_Klebs sequence, such as R36, E70, K74, F75, E78, E79, E103, E107, R108, L109, R110, S111, Q134, T139, V140, A141, S142, T144, K194, K195, T197, Q198, K204, W205 and T211, interacted directly with TLR4 receptor and, consequently, exhibited limited flexibility (Fig. 7E). Notably, the ME_Klebs residues interacting with TLR-2 were more flexible than those of TLR-4 (Fig. 7E), most probably due to the high affinity of the chosen molecular mycobacterial adjuvant, HBHA, to TLR-4 [47]. The high flexibility of the residues would also explain the later convergence of the ME_Klebs-TLR-2 complex (Fig. 7C).

Analysis of the interactions between the vaccine candidate and the receptors showed notable hydrogen bonds and hydrophobic interactions between the corresponding side and main chains of each complex pair. Interactions between receptors and ligands are typically the most important indicator of the stability of a complex. The residues involved in hydrophobic and hydrogen bonding between the ME_Klebs and the innate immune receptors can be observed in Fig. 8 and are extensively listed in Table 7. Several of the residues of both the receptors and ME_Klebs exhibited hydrogen bond interactions to more than one residue, and the same was estimated for the hydrophobic interactions observed between TLR-2 and ME_Klebs residues.

Fig. 8

Prominent hydrogen bonds and hydrophobic interactions between the vaccine candidate ME_Klebs (blue) and the receptors (green) TLR-2 (A) and TLR-4 (B). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Table 7

ME_Klebs-receptor interactions over the course of the molecular dynamics simulation.

Hydrophobic interactions		All Hydrogen bonds
TLR-2 Residues	ME_Klebs Residues	Donor Residues TLR-2	Acceptor Residues ME_Klebs	Distance (Å)
L11	Y475	R321	D274	1.78–2.83
Y111	L440P491	Donor Residues ME_Klebs	AcceptorResiduesTLR-2	Distance (Å)
		R352*	D263E264	1.7–1.91.8–2.4
		R353*	D233D235	1.80–2.41.7–1.9
TLR-4 Residues	ME_Klebs Residues	Donor Residues TLR-4	Acceptor Residues ME_Klebs	Distance (Å)
F408	W205	Q616	S111	2.9–4.4
V411	F75W205F216	K362*	E78E79	2.732.71–3.42
		Q505	T197	2.18–3.31
		Donor Residue (ME_Klebs)	Acceptor (TLR4)
		S142	C20	2.2
		W205	R382	2.75–3.1
		Q198	Q505	2.68–3.77
		R110*	D580E608	1.59–3.051.68–4.0
		Q113	D580	2.88–3.11

* Multiple occupancy observed over the course of simulation with range of distance shown for the same residue where applicable.

The MM/PBSA free binding energy calculations were done as described elsewhere for the 2 complexes [109], [110], and indicate a relatively strongly bound ME_Klebs-TLR-4 complex (ΔG = − 3,786.2 ± 1,501.5 kJ/mol) complex as compared to TLR2-ME_Klebs (ΔG = − 1,390.5 ± 858.1 kJ/mol) (Table 8). The energy changes over the course of the simulation, indicating that the ME_Klebs-TLR-4 complex achieved a stable complex with lower higher binding energy earlier in the simulation, while the ME_Klebs-TLR-2 required comparatively higher binding energy to reach stability. This is in agreement with the RMSD simulation’ results (Fig. 7C, D). Among the several individual MM/PBSA binding energy values (Table 8), the main contributor was the electrostatic energy between the receptors TLR-4 and TLR-2 and the vaccine candidate ME_Klebs, with −7,114.1 ± 1,953.3 kJ/mol and −2,354.9 ± 853.5 kJ/mol, respectively. These results are in agreement with the multiple hydrogen bond interactions predicted (Table 7), which are part of the electrostatic interactions between the receptors and the ligand observed during the simulations. Conversely, this is in disagreement with the docking results, which predicted that the ME_Klebs-TLR2 complex is less stable than the ME_Klebs-TLR-4 complex. However, the difference in the relative binding free energies between the docked complexes is small (−7.57 kcal/mol), and a possible explanation could be because of the rigidity employed by the docking approach for the receptor, which means only the binding pocket is flexible, in contrast to an MD approach, where the entire complex is flexible.

Table 8

MM/PBSA binding energy terms for ME_Klebs-TLR-2 and ME_Klebs-TLR-4 complexes.

Energy Termsa	ME_Klebs-TLR-2	ME_Klebs-TLR-4
van der Waals	−1,092.6 ± 536.4	−1,142.2 ± 318.9
Electrostatic energy	−2,354.9 ± 853.5	−7,114.1 ± 1,953.3
Polar energy	2,163.5 ± 737.8	3,503.0 ± 1,427.2
Solvent-accessible surface area (SASA)	−35.2 ± 345.3	−167.9 ± 37.3
Binding Energy	−1,390.5 ± 858.1	−3,786.2 ± 1,501.5

Values are average of 3 independent runs.

The MD simulations’ data for both complexes further support the molecular docking results described above, strengthening the interaction potential of ME_Klebs with the 2 innate immune receptors.

Virtual ME_Klebs cloning and immunization

The precedent simulations strengthened the potential of ME_Klebs as a potent, stable antigen able to generate both innate and adaptive immunological responses. This antigen would need to be produced in adequate expression vectors. As described above for the Klebs#1–4 antigens, the expression vector pET-9d showed qualitative high levels of expression in the E. coli BL21-Gold (DE3) strain (Fig. S2). Consequently, ME_Klebs sequence was codon optimized for better expression in E. coli, showing a CAI of 0.87 and a GC content of 57 %, indicative of good expression in this bacterial host [38]. ME_Klebs optimized sequence was then inserted into pET-9d sequence at the NcoI and BamHI restriction sites (Fig. 9), and the resulting cloning product could be used to transform the referred E. coli strain to produce the antigen that, following purification as described elsewhere [40], [41], would be used in immunization protocols.

Fig. 9

Schematic pET-9d plasmid map containing the codon optimized ME_Klebs for expression in E. coli strains.

ME_Klebs was used as an immunogen in a simulated human immunization protocol, 3-doses one month apart each, without adjuvant, using the C-ImmSim online tool. As expected, after each immunization, an increase in B cells (Fig. 10A) and T helper cells (Fig. 10B) was observed, both effector and memory, and the former tend to decrease 4 months after initial immunization. Also, IgM+ B cells peaked before IgG+ B cells and T helper cells, immediately after each immunization and a few days before the latter cells (Fig. 10A, B). Regarding the humoral response, immunoglobulin levels (both IgG and IgM) tend to peak around 10 days after each immunization, until all antigen (ME_Klebs) was gone (Fig. 10C), due to cellular uptake and subsequent degradation, mounting an immune response against it. This is in agreement with experimental immunization studies against bacterial pathogens done in mice using either a recombinant protein or a multiepitope chimeric antigen, in which an initial IgM, T cell independent response is observed, followed by production of IgG, in a T cell dependent manner [111], [112]. Also, both memory B cells and memory T helper cells lasted several months in this simulated immunization (Fig. 10A, B), indicating that ME_Klebs is able to induce a sustained and long-lasting humoral response.

Fig. 10

Simulated human immune response induced by three immunizations with ME_Klebs at days 0, 28 and 56, with the online tool C-ImmSim. Both B cell (A) and T helper cell (B) immunological responses were evaluated up to 350 days after first immunization. Levels of different immunoglobulin classes (IgM, IgG1 and IgG2) and their relationship with the ME_Klebs antigen levels, after each immunization, were also assessed (C).

Conclusions

Identification of potential vaccine antigens is based on a long, laborious, and costly experimental procedures, such as identification of potential immunogens, their production and purification and immunogenicity studies in adequate animal models. In the last decades, several in silico approaches and tools became available, allowing a computational analysis of potential antigens and their characteristics, screening and filtering thousands of proteins and only pursuing at experimental level the most promising immunogens, reducing greatly the costs and time associated with their validation. In this work, potential antigens against carbapenem-resistant K. pneumoniae strains were identified by reverse vaccinology and genome mining. However, both in silico characterization and experimental procedures suggested that these antigens were not ideal, due to their insolubility and, in the case of Klebs#3/MrkD, potential allergenicity. Thus, a multiepitope approach was pursued instead and a vaccine construct was designed, ME_Klebs, containing B cell and helper T cell epitopes of the 4 identified antigens. Due to the importance of TLR-4 in protection against K. pneumoniae infection, a molecular adjuvant, which activates this innate immunity receptor, was added. In silico online tools prediction results indicated that ME_Klebs was indeed antigenic, soluble, non-allergen and non-toxic. Predictions of secondary and 3D structures produced models of high quality, and, using the latter model, several putative conformational B cell epitopes were identified. Predicted docking complexes structures and subsequent molecular dynamics allowed estimation of energies of ME_Klebs complexed with the relevant innate immunity receptors TLR-2 and TLR-4, which control K. pneumoniae infection and bacterial dissemination, indicated that this multiepitope construct stably binds and activates both receptors and, consequently, potentially triggers the innate and adaptive immunological responses initiated by both TLRs. Cloning of ME_Klebs through conventional cloning techniques appeared feasible, and the resulting purified protein is expected to be able to induce a Th-2-type immune response, eliciting a strong humoral response. In conclusion, the ME_Klebs extensive characterization in silico strongly suggests this multiepitope vaccine could be an efficient immunization approach against carbapenem-resistant K. pneumoniae clinical isolates.

CRediT authorship contribution statement

Nicola Cuscino: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Validation. Ayesha Fatima: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Validation, Writing – original draft. Vincenzo Di Pilato: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Validation. Matteo Bulati: Formal analysis, Writing – original draft. Caterina Alfano: Conceptualization, Funding acquisition, Project administration, Supervision. Elisa Monaca: Data curation, Investigation, Methodology, Resources. Giuseppina Di Mento: Data curation, Investigation, Methodology, Resources. Daniele Di Carlo: Conceptualization, Funding acquisition, Project administration, Supervision. Francesca Cardinale: Data curation, Investigation, Methodology, Resources. Francesco Monaco: Data curation, Investigation, Methodology, Resources. Gian Maria Rossolini: Conceptualization, Funding acquisition, Project administration, Supervision. Asif M. Khan: Conceptualization, Funding acquisition, Project administration, Supervision, Writing – review & editing. Pier Giulio Conaldi: Conceptualization, Funding acquisition, Project administration, Supervision. Bruno Douradinha: Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.