Introduction
Background to mRNA vaccines
mRNA vaccines have emerged as a promising alternative solution to overcome the challenges of traditional vaccines, offering a faster and more flexible approach to vaccine development. Their ability to elicit strong immune responses, coupled with adaptability to different viruses, makes them particularly attractive for combating viral diseases [5]. Unlike traditional vaccines, mRNA vaccines (Fig. 1B) eliminate the need for live pathogens by using synthetic mRNA—a small fragment of genetic material from the virus—that encodes instructions for the host’s cells to produce a specific protein, usually a viral antigen, which triggers an immune response [6]. Once inside the cells, the mRNA is translated into the viral protein by the host’s cellular machinery, which the immune system recognizes as foreign, subsequently prompting antibody production, T-cell activation, and immune memory formation that can quickly respond to future infections if the actual virus enters the body [4]. Consequently, owing to their rapid development and adaptability, mRNA vaccines such as Pfizer-BioNTech’s BNT162b2 and Moderna’s mRNA-1273 were developed against COVID-19 and deployed with remarkable speed during the pandemic [7, 8].
Despite their promise, the design of mRNA vaccine faces several challenges that must be addressed to ensure their success, including mRNA stability, codon optimization, and untranslated regions (UTRs) incorporation [9]. As the field of mRNA vaccine development continues to evolve, existing protocols often tackle these issues in isolation—focusing, for example, mostly on codon optimization or structural stability—without providing a comprehensive and integrated guidelines [10, 11]. Such an approach overlooks critical interdependencies among various design elements, such as the interplay between codon usage, RNA structure, immunogenicity, and UTR selection [12, 13]. A unified protocol that systematically guides researchers from target selection through in silico optimization and validation remains an unmet need in the field [14].
This protocol aims to fill the existing gap by providing a comprehensive computational framework for mRNA vaccine design that integrates key immunoinformatics and structural principles. It defines the design space, establishes essential functional criteria, introduces concepts for optimizing vaccine design, and provides a detailed procedure incorporating codon optimization and structural stability analysis. Using the SARS-CoV-2 spike protein as a model, the protocol delineates critical design parameters, functional efficiency criteria, and provides a detailed step-by-step procedural guide for in silico vaccine sequence development. By unifying all key elements into a cohesive pipeline, it standardizes the design process and supports the development of stable, immunogenic, and translationally efficient mRNA vaccines across diverse pathogens.
Demonstrating a comparative figure of traditional vaccines and mRNA vaccines. () traditional vaccines introduce pre-formed antigens derived from inactivated or attenuated pathogens. These antigens are captured by antigen-presenting cells (APCs), including dendritic cells and macrophages, processed intracellularly, and presented as peptides on major histocompatibility complex (MHC) molecules. Peptide–MHC class II heterodimers activate CD4T helper cells, which differentiate into TH1 and TH2 subsets; TH1 cells promote cytotoxic T-lymphocyte (CTL) activation, whereas TH2 cells support B-cell activation and antibody production. Cross-presentation by dendritic cells enables peptide loading onto MHC class I molecules, leading to CD8T-cell priming. Activated CTLs mediate effector functions through perforin and granzyme release, while plasma cells derived from activated B cells secrete neutralizing antibodies that circulate systemically in a soluble form. () mRNA vaccines deliver synthetic mRNA encoding the target antigen into host cells via lipid nanoparticles. Following endocytosis and endosomal escape, the mRNA is translated in the cytoplasm without entering the nucleus. The encoded antigen often contains a signal peptide that directs it to the secretory pathway, allowing both intracellular processing and extracellular release. Intracellularly processed peptides are presented on MHC class I molecules to activate CD8T cells, whereas secreted or extracellular antigens can be taken up by APCs and presented on MHC class II molecules to activate CD4T helper cells or directly recognized by B cells. Both pathways induce coordinated cellular and humoral immune responses. Created with A B + + + + BioRender.com
Defining design space for mRNA vaccine efficiency
Design space for mRNA vaccine
Optimizing the design space is crucial for developing mRNA vaccines that encompass several key components including UTRs, coding sequence (CDS), mRNA secondary structure, sequence length and GC content. The 5‘UTR and 3’ UTR are essential in regulating translation efficiency, stability, and overall mRNA expression, with the 5‘UTR specifically contributing to ribosome recruitment and binding and 3’ UTR influencing mRNA stability and intracellular localization [15]. At the 5′ end, the mRNA includes a 5′ cap, which facilitates ribosome recognition, translation initiation, and protection from exonuclease degradation [16]. Following the 5‘cap is the 5’ UTR, within which the Kozak sequence is embedded—a conserved nucleotide motif surrounding the start codon, that guides the ribosome to the correct initiation site [17].
Following the 5’ UTR, the CDS encodes a protein that often contains structural domains with distinct functions governing protein folding, localization, and biological activity. Within these domains, conserved short patterns may occur that contribute to specific functional roles and are commonly referred to as motifs. In the context of vaccine design, particular attention is given to antigenic regions of proteins that may contain immunogenic epitopes. These immunogenic epitopes are often found within specific functional domains of the protein, particularly within surface-exposed or extracellular regions that are accessible to immune surveillance and represents the primary targets of vaccine-induced immune responses [5, 18]. Following translation, intracellularly expressed antigens are processed through proteasomal or endosomal pathways, generating peptide epitopes that are presented on the cell surface in association with MHC class I or class II molecules for T-cell recognition, while intact or secreted antigens may be directly recognized by B cells [4]. Toward the end of the CDS, the downstream 3′ UTR plays an important role in regulating mRNA stability, localization, recycling, and translation termination [19]. This region interacts with cellular machinery to protect the mRNA from degradation, thereby extending its half-life and supporting sustained protein production. Proper optimization of the 3′ UTR involves the incorporation of regulatory elements that enhance translation efficiency and prolong mRNA stability, with specific sequence features enabling productive interactions with protective cellular factors [20].
Secondary structure, length and GC content of mRNA vaccine
The secondary structure of mRNA is another key design element, as its folding can directly influence translation efficiency. Features such as hairpins and loops, particularly near the 5′ UTR or start codon, can impede translation efficiency [21]. UTRs are not merely linear sequences but often form these complex secondary structures that can create barriers to ribosome scanning and translation initiation, especially in the 5′ UTR [15]. Other elements in the design space are the length and nucleotide composition of the mRNA sequence that significantly impact mRNA stability and translation efficiency. Longer mRNA sequences are generally more prone to degradation and often require stabilizing modifications and length optimization [22]. Additionally, the GC content, or the proportion of guanine and cytosine nucleotides, affects the mRNA’s stability and secondary structure [23]. High GC content may form stable secondary structures that hinder translation, while low GC content can result in less stable mRNA [24]. Thus, efficient mRNA design must holistically consider 5′ and 3′ UTRs, CDS, secondary structure, sequence length, and GC content to ensure stability, strong protein expression, and minimal immunogenicity.
Essential criteria for mRNA vaccine functionality
The successful design of mRNA vaccine sequences hinges on critical criteria that ensure their functionality within the host, including the 5′ and 3′ UTRs, CDS, secondary structure stability, and immunogenicity. These factors determine how effectively the mRNA produces the target antigen while maintaining stability and safety. Optimizing the 5’ UTR selection must be carefully done to promote effective ribosome binding and avoid forming highly stable secondary structures that could hinder translation initiation or may prevent translation [4]. Selecting UTRs from host sources known to enhance mRNA stability and translation ensures compatibility with the host’s translational machinery, thereby promoting efficient protein synthesis [25]. Other criteria for UTR selection include length and structural complexity, with shorter 5’ UTRs generally preferred, as longer, highly structured regions may reduce translation efficiency [15]. Upstream start codons (AUGs) or upstream open reading frames (uORFs) within regulatory regions, especially in the 5′ UTR, should be avoided as they interfere with translation of the main open reading frame (ORF) through premature initiation or ribosome scanning disruption [26, 27]. Codon usage can also influence translation efficiency—mainly within the CDS—by modulating elongation rates, as will be further illustrated in this section. The combination of the MFE of both the 5′ UTR and the CDS should be calculated to assess overall stability without hindering translation, while still protecting the mRNA from degradation [28]. MFE should be evaluated comprehensively since UTR–CDS interactions affect both stability and translation. UTR optimization further requires balancing this structural stability, avoiding excessive folding that impairs translation while maintaining adequate protection [29, 30]. The 5′ UTR of human TMSB10 or Ces1d, for example, have been shown to enhance antigen-specific humoral and cellular immune responses when used in mRNA vaccines encoding the receptor-binding domain (RBD) of SARS-CoV-2, demonstrating superior expression in dendritic cells [31]. Likewise, the 5′ UTR of human β-globin is well-documented for its strong ability to promote cap-dependent translation initiation, thereby boosting protein synthesis in mRNA applications [32].
The 3′ UTR also requires careful consideration to avoid unintended immune responses and prolong its half-life, as it is vital for mRNA stability and cellular localization. Selecting 3′ UTR sequences that promote stability while avoiding immune-triggering motifs or premature degradation is critical. For instance, incorporating polyadenylation signals within the 3’ UTR can extend mRNA half-life and ensure sufficient protein production [33, 34]. The 3′ UTRs of human Apo A-II and AP3B1 have proven effective in extending mRNA half-life and supporting high protein output in vivo, especially when paired with strong 5′ UTRs [31]. In comparative screening studies assessing multiple UTR combinations, the pairing of the 5′ UTR from TMSB10 with the 3′ UTR of Apo A-II emerged as one of the most effective human UTR configurations when combined with viral coding sequences such as that of the SARS-CoV-2 spike protein, yielding significantly increased antigen expression and immunogenicity. These findings support the rationale for utilizing empirically validated UTRs in combination with viral CDS to optimize mRNA vaccine performance.
The CDS within mRNA contains crucial elements that significantly impact vaccine effectiveness. For instance, short motifs of 6–15 amino acids are often associated with receptor binding or immune activity, whereas longer motifs as 20–50 amino acids may correspond to structural elements such as fusion peptides or heptad repeats [35]. Additionally, certain regions of the encoded proteins, particularly those that are surface-exposed, may harbor the antigenic epitopes recognized by the immune system [36]. Epitopes are generally classified as linear and conformational peptide types. Linear epitopes are often favored, as they consist of continuous amino acid sequences meaning the residues are adjacent in the primary sequence, making them easier to identify and more reliably predicted for peptide-based vaccine design [18]. In contrast, conformational epitopes are discontinuous, formed by amino acids that are distant in sequence but brought together in the folded protein structure. Effective epitope selection should include both CD8+ and CD4+ T-cell epitopes, requiring compatibility with MHC class I and class II human leukocyte antigen (HLA) molecules to ensure a coordinated cellular immune response [37, 38]. HLA molecules are cell-surface proteins responsible for presenting peptide antigens to T cells, thereby enabling immune recognition and activation. The binding affinity between the epitope and HLA determines how efficiently the antigen is recognized and triggers immune response. Signal peptides optimization is another essential aspect for directing newly synthesized proteins into the secretory pathway, influencing protein folding, secretion, and overall expression levels [39]. These short, N-terminal sequences are found in both prokaryotes and eukaryotes and play a vital role in recombinant protein production [39, 40].
The immunogenicity must also be optimized to trigger an appropriate immune response without causing adverse effects such as autoimmunity [41]. The innate immune system recognizes exogenous mRNA via pattern recognition receptors, including endosomal Toll-like receptors (TLRs) and cytosolic RNA sensors such as RIG-I and MDA5, which can lead to inflammation or transcript degradation if excessively activated [42]. Overactivation of these pathways may amplify type I interferon signaling, promote local inflammatory responses and ultimately reduce antigen expression [43]. Consequently, candidate antigen sequences should be systematically screened for antigenicity, allergenicity, and toxicity as essential safety filters during vaccine design [44]. Population coverage analysis also represents a critical design criterion in epitope-based vaccines. It refers to the proportion of individuals within a given population who are predicted to mount an immune response based on HLA allele distribution. Because HLA polymorphisms vary across ethnicities and geographic regions, selecting epitopes that bind to multiple common HLA alleles increases the likelihood of broad vaccine effectiveness at the population level [45].
To enhance immunogenic breadth and ensure the induction of both cellular (T-cell-mediated) and humoral (antibody-mediated) immune responses, multiple immunodominant epitopes are integrated into a multi-epitope construct for MEVC formation, allowing simultaneous presentation of cytotoxic T lymphocytes (CTL), helper T lymphocytes (HTL) and B-cell epitopes within a single sequence [46]. The inclusion of an adjuvant is essential, as short peptide epitopes alone often exhibit limited intrinsic immunogenicity; therefore, molecular adjuvants are incorporated to enhance immune activation. Among these, β-defensin is commonly used due to its ability to stimulate innate immunity, promote dendritic cell maturation, enhance T-cell responses and its use has been consistently reported in multi-epitope vaccine studies [47–49]. Essential aspect here is the rational linker selection which is critical to prevent the formation of junctional epitopes and to facilitate correct epitope presentation. Once the multi-epitope sequence is assembled, it should be evaluated for immunogenicity and for physicochemical properties, including molecular weight, stability, hydrophobicity (GRAVY, the average hydropathy score, where positive values indicate hydrophobicity and negative values indicate hydrophilicity), secondary structure tendencies, and predicted in vivo and in vitro half-life [50].
Another critical aspect of mRNA design is codon optimization which is fundamental in the design space. Codons determine amino acid incorporation during protein synthesis. Organisms exhibit species-specific codon usage preferences, referred to as codon bias, whereby certain codons are more efficiently recognized by the translational machinery, resulting in enhanced protein expression [51]. Optimization of codon usage bias should be aligned with host preferences as human codon usage to enhance translation efficiency and reduce ribosomal stalling [52]. Chemical nucleotide modification is an important factor in mRNA functionality that complements sequence-level optimization in modern mRNA vaccine design [53]. In current platforms, standard ribonucleotides are often replaced with natural or synthetic analogues to enhance stability, translational efficiency, and tolerability. Uridine is commonly substituted with N1-methylpseudouridine (m1 Ψ) or pseudouridine (Ψ), as used in clinically approved SARS-CoV-2 mRNA vaccines, reducing innate immune activation while enhancing protein expression in vivo [53]. Additional modifications, including 5-methylcytidine (m5C), N6-methyladenosine (m6A), and 2-thiouridine (s2U), have also been explored for their roles in improving mRNA stability, translation, and immunogenicity. Comparative studies show that Ψ- and m1 Ψ-modified mRNA increase duplex stability and outperform unmodified constructs in protein production while maintaining improved safety profiles [54–56].
The secondary structure is another criterion that must be carefully considered as naturally folded mRNA can hinder ribosome access or promote degradation. Predicting folding patterns of MEVC/CDS vaccine and identifying stable regions allows adjustments to minimize unwanted folding [57]. Structural stability is commonly assessed by minimum free energy (MFE), a thermodynamic measure of RNA folding; lower MFE values indicate more stable structures that can protect mRNA from degradation and promote efficient translation while high MFE increase degradation risk [11, 58]. As MFE scales with RNA length and complexity, values must be interpreted relative to length and function. The length and sequence composition of UTRs and CDS also influence secondary structure formation, which must be minimized to facilitate ribosome binding [15]. Optimization of UTRs and the 5’ cap structure helps reduce unintended innate immune sensing and improves transcript stability [59, 60]. These steps should be performed alongside assessment of target accessibility to ensure mRNA availability to the translational machinery and additional two-dimensional (2D) structure analysis for the optimized MEVC to predict α-helices, β-sheets, and coil regions within the protein sequence of the vaccine construct [58].
Cloning efficiency is another important consideration, as suboptimal insertion can compromise construct amplification and necessitate further optimization. Accordingly, the MEVC or CDS constructs are introduced into a suitable expression system, most commonly an Escherichia coli plasmid such as pET-28a(+), to evaluate protein production efficiency [61]. In this context, the BamHI restriction enzyme recognition site is typically incorporated at the N-terminal 5′ end of the insert, while the XhoI site is introduced at the C-terminal 3′ end, ensuring correct orientation of the construct relative to the promoter and minimizing the likelihood of vector self-ligation [62–64].
Once the final construct is obtained, several aspects should guide the evaluation of vaccine-induced immune responses. Immune simulation serves as an essential validation step, where the immunization strategy—particularly dose, route, and scheduling, including prime–boost regimens—shapes the response [63]. Both humoral and cellular immunity should be assessed, including antibody production, neutralization capacity, and T-cell activation and phenotype. Cytokine and interleukin profiles provide insight into immune signaling and the balance between effective and excessive inflammation. The development of immunological memory through B- and T-cell responses should also be evaluated, together with antigen clearance to determine immune efficiency. Safety and tolerability remain critical to ensure a well-regulated response. The mRNA delivery system, most commonly lipid nanoparticles, further supports vaccine performance by protecting mRNA, promoting cellular uptake, and enabling endosomal escape for efficient translation [4, 9, 65]. Advances in lipid nanoparticle composition and formulation have demonstrated that delivery efficiency strongly influences antigen expression levels, immune activation, and overall vaccine performance [66, 67]. Although delivery optimization is beyond the scope of this protocol, its impact on vaccine performance is acknowledged as a critical factor for achieving effective and safe antigen expression in vivo.
Design mRNA vaccine for specificity
Epitope identification is subsequently performed using computational tools, which can be categorized into cytotoxic T lymphocyte (CTL; MHC class I-restricted), helper T lymphocyte (HTL; MHC class II-restricted), linear B-cell, and conformational B-cell epitope prediction [70]. For T cells, Binding affinity is evaluated using IC₅₀ values, which represent the half-maximal inhibitory concentration and reflect peptide–MHC binding strength; IC₅₀ thus serves as a key immunological criterion for prioritizing protein regions most likely to be presented by prevalent HLA molecules, guiding their inclusion in vaccine constructs [71]. Lower IC₅₀ values indicate stronger binding and higher likelihood of immune recognition; epitopes with IC₅₀ <50 nM are classified very strong binders, while those < 500 nM are considered regular binders [72]. Epitopes are further prioritized based on prediction-derived metrics, particularly the percentile rank and prediction score. Lower percentile ranks correspond to higher binding affinity, with values below 1% considered strong binders and values up to ~1.5–2% regarded as moderate- to high-affinity binders, while higher prediction scores indicate stronger peptide–MHC interactions [73]. For B-cell epitope prediction, higher scores are preferred, as they indicate a greater likelihood that the region represents a B-cell epitope and is recognized by antibodies. In vaccine design, epitope selection and filtering are based on a combination of percentile rank, high score and the frequency of experimental validation, such as the number of supporting assays reported in Immune Epitope Database (IEDB) [74]. Each assay represents an independent experimental validation performed ex vivo, in vitro, or in vivo and is linked to a reference identifier corresponding to the original study. Epitopes supported byhigh number of experimental assays are considered more reliable and prioritized for robust vaccine design [75]. The target’s specificity and immunogenicity are then validated by confirming that the selected epitopes elicit strong immune response, avoiding unintended immune interactions and allergenicity. Once the non-allergenic, non-toxic and antigenic epitopes are identified, their coding selection and structural stability are also evaluated by assessing translation efficiency, folding stability, and preservation of the original amino acid content. Following this, population coverage analysis is performed in epitope-based vaccine design, as it evaluates whether the selected T-cell epitopes MHC class I and MHC II and their associated strong binding HLA alleles can provide broad immune coverage across different ethnicities and geographical regions worldwide.
After pinpointing the immunologically relevant regions, the CDS corresponding to the selected protein or MEVC is retrieved. If starting from a consensus protein sequence, the equivalent CDS is derived through reverse translation. This CDS is then optimized in the next step through humanization and codon optimization to improve translation efficiency, aligning with host as human codon usage preferences which is normally assessed by a codon adaptation index (CAI) value. The next step in the design involves an iterative refinement process using the LinearDesign algorithm [11]. The LinearDesign algorithm focuses on optimizing the CAI and the MFE of the CDS region of the selected protein or MEVC [78]. CAI values near 1.0 indicating strong optimization and typically correlating with higher translation efficiency and protein expression, while GC content should typically fall within the range of 30–70% to ensure transcriptional stability [8, 78]. Following the optimization of the coding sequence based on CAI, calculate CAI and MFE values for the MEVC in addition to the full CDS of the spike protein. When using algorithms or tools that perform codon optimization without incorporating MFE calculations, then mRNA secondary structure and accessibility should be assessed after coding optimization. Additional 2D structure analysis is recommended for the optimized MEVC using tools based on position-specific scoring matrices to predict α-helices, β-sheets, and coil regions, while assigning confidence scores to each residue, where higher scores indicate greater prediction reliability [79].
Computational approaches, such as in silico screening, identify optimal signal peptides, ensuring the vaccine protein is efficiently expressed and presented to the immune system [80, 81]. Subsequently, in silico cloning is performed as a validation step to assess the feasibility of construct integration into an expression system and to ensure the success of subsequent experimental validation. In this context, the designed sequence, referred to as the insert, is virtually integrated into an expression vector, typically a plasmid. The pET-28a(+) vector is commonly used due to its strong T7 promoter, efficient transcriptional control, and compatibility with E. coli expression strains [63, 64]. During in silico cloning, restriction enzymes such as BamHI and XhoI are selected to ensure correct orientation and prevent self-ligation, as they are present in the multiple cloning site (MCS) and generate non-compatible cohesive ends for directional cloning. Specifically, the BamHI recognition site is introduced at the N-terminal (5′) end of the insert, while the XhoI site is incorporated at the C-terminal (3′) end, ensuring proper orientation of the construct within the vector. Following in silico cloning, immune simulation serves as the last in silico validation step that provides a preliminary assessment of the vaccine’s ability to induce an immune response prior to experimental testing [63]. A key consideration in this process is the simulation of a realistic vaccination regimen, typically following a single-dose or multi-dose (prime–boost) schedule, where repeated antigen exposures are used to mimic initial immunization and subsequent booster doses. This includes assessment of antibody production, T-cell activation, cytokine signaling, memory formation, and antigen clearance under simulated immunization conditions [82]. A progressive reduction in antigen levels after each simulated dose is generally indicative of an effective immune response.
The next optional step is UTR optimization, integrating 5′ and 3′ UTRs with the CDS/MEVC to balance structural stability and translation efficiency. In this optional step, Host-derived UTRs as human UTRs are combined with viral CDS, and multiple UTR combinations are systematically tested, with total combinations calculated as (number of 5′ UTRs × number of 3′ UTRs). As for mRNA vaccine applications (Fig. 3), the optimized coding sequence or MEVC is incorporated into a complete mRNA architecture that includes a 5′ cap structure to initiate translation, a 5′ untranslated region to enhance ribosome binding, the coding sequence encoding the multi-epitope construct, a 3′ untranslated region to improve transcript stability, and a poly(A) tail to further enhance stability and translation efficiency [83, 84].
Optimization of mRNA vaccine design. The crucial three steps in mRNA vaccine sequence optimization are presented in the figure: (1) sequence annotation & UTR optimization, where viral mRNA coding sequence is annotated into domains, motifs, surface-exposed regions, and epitopes. The 5′ and 3′ UTRs are optimized for stability and translation efficiency. (2) secondary structure optimization, and stability of different sequences calculated using MFE, a lower MFE indicating higher stability. Sequence E, MFE (−300), most stable. Codon adaptation index (CAI) is also computed, and sequence E has the maximum CAI (0.95), indicating greater translation efficiency. (3) target accessibility, where mRNA is optimized for ribosomal binding and circumvented for undesirable interactions. Accessible mRNA is too accessible and can engage immune sensors like TLR7 or lead to RNase-mediated degradation and compromise vaccine efficacy. Created with BioRender.com
Schematic representation of the multi-epitope mRNA vaccine construct. The design includes an N-terminal adjuvant followed by HTL, CTL, and B-cell epitopes linked using EAAAK, GPGPG, AAY, and KK linkers to facilitate proper processing and presentation. The coding sequence is flanked by regulatory elements, including the 5′ cap and 5′ UTR upstream and the 3′ UTR with a poly (A) tail downstream, ensuring efficient translation and stability. Created withand drawio.com MS word
Procedure: applying the protocol to the spike protein of SARS-CoV-2
Procedure workflow for mRNA vaccine design. The figure outlines a multi-step pipeline for designing the multi-epitope vaccine construct (MEVC) and the full coding sequence (CDS), starting with target sequence retrieval from UniProt, NCBI protein, and GenBank. Protein domains, motifs, and epitopes are annotated and predicted using pfam, MEME suite, NetMHCpan, and IEDB. Predicted epitopes and the full-length CDS are evaluated for antigenicity, allergenicity, and toxicity using VaxiJen, AllerTOP, and ToxinPred, followed by physicochemical analysis with ProtParam. The MEVC is then constructed and re-evaluated for immunogenicity using the same tools. The coding sequence of MEVC is generated by reverse translation and optimized using LinearDesign and NovoPro ExpOptimizer, followed by RNAfold analysis for MFE and structural accessibility assessment then PSIPRED for 2D secondary structure prediction. Finally, validation with in silico cloning is performed using SnapGene, and immune simulation validation using C-ImmSim. Created with drawio.com
Step 1: Target Gene Assignment and Sequence Retrieval
The first step is to obtain the target sequence, in this protocol, the spike protein of SARS-CoV-2 (accession number P0DTC2). To retrieve it from the UniProt database (https://www.uniprot.org), search using the accession number and download the sequence in FASTA format. The same sequence can be obtained from the NCBI Protein database (https://www.ncbi.nlm.nih.gov) by selecting “Protein” in the search menu, entering the accession number, and downloading the FASTA file. If only the gene name is available, use NCBI or GenBank (https://www.ncbi.nlm.nih.gov/genbank/) to locate the corresponding protein entry by filtering results for the correct species, then obtain its accession number and FASTA sequence. This FASTA file will serve as the reference for all downstream optimization steps.
⚠ CRITICAL STEP: If the target protein is unknown, as in unannotated viral genomes, the entire genome is scanned to identify potential targets. Conserved regions across multiple strains are located using Clustal Omega (https://www.ebi.ac.uk/jdispatcher/msa/clustalo). Genome annotation is performed with Prokka (https://github.com/tseemann/prokka) to extract ORFs and protein-coding regions from the FASTA genome sequence. Predicted ORFs are translated via ExPASy Translate (https://web.expasy.org/translate/) and compared with viral protein databases (https://www.ncbi.nlm.nih.gov/labs/virus/) to exclude non-functional sequences, after which the downstream workflow pipeline is applied to generate the final vaccine sequence for both the MEVC and full length CDS.
Step 2: Annotation of Protein Domains, Motifs and Epitopes
This step involves annotating the sequence to identify conserved domains, motifs, and potential epitopes using annotation and prediction tools, as outlined below:
Identify conserved protein domains and extracellular domains
Go to the Pfam website (https://pfam.xfam.org) and paste the P0DTC2 sequence then click on the “Search” tab. Pfam will then search its database of known protein domains and return a list of any matches found in the input sequence. For each match, Pfam provides detailed information about the domain, including its name, function, and location within the sequence. The output is a comprehensive annotation of protein domains within the spike protein, which can be used for further structural and functional analysis.
Discover conserved motifs
Open MEME Suite website (https://meme-suite.org) and select the “MEME” tool. Paste the P0DTC2 sequence and set the parameters based on known motif characteristics: 5–10 motifs, width 6–50 amino acids. This range captures key immunological and structural regions as short motifs (6–15 aa) often correspond to receptor binding domains or epitopes, while longer ones (20–50 aa) may represent structural domains like the fusion peptide or heptad repeats. Choose “zero or one occurrence per sequence” in protein mode to extract unique motifs, or “zero or more occurrences” to detect recurring motifs. Click “Submit”. MEME outputs statistically significant motifs with consensus sequences, highlighting conserved functional elements relevant to viral function and immune recognition.
MHC class I (CTL) and MHC class II (HTL) T cell epitope prediction
For MHC class I epitope prediction, open the IEDB (https://www.iedb.org/), select the T-cell epitope prediction tool (https://nextgen-tools.iedb.org/pipeline?tool=tc1) and input the sequence. Set all parameters to default and select a peptide length of 9 amino acids and choose the 27 human HLA class I alleles that provide approximately 97% population coverage. Run the prediction and export the results. For MHC class II epitope prediction, within the IEDB, open the tool (https://nextgen-tools.iedb.org/pipeline?tool=tc2), input the sequence and keep all parameters at default. Set the peptide length to 15 amino acids, and select the 27 human HLA class II alleles to ensure comparable population coverage. Run the prediction and export the results.
Linear/Confirmational B cell epitope prediction
For linear B-cell epitope prediction, within the IEDB open the tool (https://tools.iedb.org/bcell/) and input sequence selecting a prediction method such as BepiPred or Karplus–Schulz flexibility, which evaluates residue flexibility as an indicator of antigenicity. Run with default parameters. For conformational/discontinues B-cell epitope prediction, within IEDB, use ElliPro tool (https://tools.iedb.org/ellipro/) and input the validated 3D structure of the SARS-CoV-2 spike protein (obtain it from the Protein Data Bank; https://www.rcsb.org/; PDB ID: 6VXX). In ElliPro, input the PDB ID and run with default parameters, including minimum score (protrusion index/PI threshold) required for a residue to be considered part of an epitope and maximum distance (Å) to define the clustering radius for grouping residues into epitopes. Submit the job, select all chains (A, B, and C), and resubmit to obtain results.
Alternative method
Visualization of SARS-CoV-2 spike protein annotations using the UniProt feature viewer. The UniProt feature viewer displays the spike protein sequence (1–1273), highlighting curated annotations including domains, motifs, topology (surface-exposed or extracellular regions of transmembrane proteins), and epitopes. Features are color-coded and mapped to sequence positions. GFF-formatted annotations with start/end positions and descriptions can be downloaded. A categorized panel allows selection of annotation layers such as molecular processing, topology, and structural features. Created with UniProt
Step 3: Specificity and Immunogenicity Filtering, Validation, and MEVC Design
MHC class I (CTL) and MHC class II (HTL) T cell filtering
For MHC class I/II epitopes, filter the results based on two main thresholds criteria. First, select epitopes with the lowest percentile rank (NetMHCpan EL percentile), as lower values indicate stronger binding affinity; consider values < 1% as strong binders and values between 1 and 2% as moderate binders. Second, select epitopes with the highest prediction score (NetMHCpan EL score), where values closer to 1 indicate stronger predicted binding affinity. If IC₅₀ information is available through the IEDB then a third optional filtration criteria can be applied, by selecting those with IC₅₀ less than 50 nM for very strong binders or less than 500 nM for average binders. IC₅₀ information could also be access through NetMHCpan server (https://services.healthtech.dtu.dk/services/NetMHCpan-4.1/) if not available in the IEDB.
epitope filtering Linear/Confirmational B cell
For linear B-cell epitopes, the tool assigns scores to residues or peptide segments, where higher scores indicate greater likelihood of antibody recognition. A threshold is automatically defined and region above it is considered potential epitopes. Select continuous high-scoring regions above the threshold as the most probable linear B-cell epitopes. For conformational/discontinues B-cell epitopes, the tool provides predicted epitopes with scores reflecting surface exposure and protrusion, where higher scores indicate greater accessibility. Select top-scoring epitopes and visualize them using the “View 3D structure” option to assess spatial distribution and surface accessibility.
IEDB assays validation
After filtering, validate the selected epitopes using the IEDB search interface by setting the host to “Human,” selecting both MHC class I and class II, then select all assay types (B-cell, T-cell, MHC Ligand) and filter for SARS-CoV-2 spike protein epitopes. Select epitopes based on binding affinity and immune assays. Prioritize epitopes based on the number and type of associated experimental assays as T-cell activation and MHC binding assays, then review their ID to reference research articles to assess validation strength.
Validation of target specificity and immunogenicity
Assess sequences for antigenicity, toxicity, and allergenicity; if unsuitable, return to epitope selection. Use VaxiJen (threshold 0.4, “Virus” setting; http://www.ddg-pharmfac.net/vaxijen/VaxiJen/VaxiJen.html) to predict antigenicity scoring from FASTA sequences. Evaluate toxicity with ToxinPred (https://webs.iiitd.edu.in/raghava/toxinpred/) using “protein scanning” for full CDS, or the “Batch Submission” for the epitopes (excluding toxic peptides). Assess allergenicity using AllerTOP (https://www.ddg-pharmfac.net/AllerTOP/) by submitting FASTA sequences. Apply all evaluations to both selected epitopes and the full CDS.
Population coverage analysis
Set the number of epitopes according to the final selected list. Choose the query option area_country_ethnicity and select the targeted population (here we selected World to evaluate global population coverage). From MHC binding prediction results (NetMHC for class I and NetMHCII for class II), identify the most relevant HLA alleles for each epitope and select only the strongest binders based on binding affinity rank (Supplementary table S1), apply a strict cutoff of percentile rank ≤ 1.6. For each epitope, enter the sequence with its corresponding selected HLA alleles in the input table, ensuring correct matching. Select Class I and II combined to assess overall immune coverage, keep other parameters as default, and submit the analysis. After submission, review outputs including a graphical distribution of epitope recognition, population coverage percentage, average hit and PC90 that indicates the minimum number of epitope–HLA combinations recognized by 90% of the population.
vaccine Multi-epitope mRNA construct
Construct the multi-epitope vaccine by concatenating selected CTL, HTL, and B-cell epitopes into a continuous amino acid sequence using a biologically guided strategy, ensuring high binding affinity, antigenicity, non-allergenicity, and non-toxicity while avoiding redundancy. Begin with an immunostimulatory adjuvant at the N-terminus, such as human β-defensin 2 (HBD2), retrieve its sequence from (https://www.rcsb.org/structure/1FD3), followed by a rigid EAAAK linker to maintain structural separation. Arrange epitopes in a defined order: HTL epitopes first (linked by GPGPG) to support MHC class II presentation, followed by CTL epitopes (AAY) to enhance proteasomal processing and MHC class I presentation, and finally B-cell epitopes (KK) to preserve antibody accessibility. Assemble the construct as: adjuvant → EAAAK → HTL (GPGPG) → CTL (AAY) → B-cell (KK). Verify sequence continuity, correct reading frame, and absence of premature stop codons, and evaluate the construct as a single protein. Although simplified designs without linkers are possible, rational linker-based assembly is preferred; thus, the MEVC strategy was applied throughout this protocol.
To assess immunogenicity of the MEVC sequence, use the same tools as in the previous step: VaxiJen for antigenicity, AllerTOP for allergenicity, and ToxinPred for toxicity. In addition, use ProtParam (https://web.expasy.org/protparam/) to evaluate physicochemical properties of MEVC. Paste the MEVC sequence and submit. Check from the output, the molecular weight, stability, hydrophobicity, secondary structure tendencies, and in vivo and in vitro half-life, to confirm structural stability and suitability for expression. Verify that linker inclusion does not disrupt epitope integrity or reduce predicted binding affinity. In this protocol, five epitopes were selected from each category (MHC class I, MHC class II, linear B-cell, and conformational B-cell). Initial filtering was based on high binding affinity (percentile rank ≤ 1.6) and highest prediction scores, with antigenicity assessed only for this subset. Final selection considered immunogenicity, antigenicity, and validation in IEDB assays. Linear B-cell epitopes were further concatenated in pairs to enhance antigenic coverage and potentially strengthen humoral responses, consistent with validated IEDB regions. Although optional, this step was applied to improve epitope representation and immunogenic potential.
Step 4: Codon Optimization and MFE Assessment
First, retrieve the CDS for the selected protein or epitopes/MEVC. Download directly from databases when available; otherwise, generate it using a reverse translation tool. For SARS-CoV-2 spike, obtain accession P0DTC2 from UniProt and follow cross-references to the GenBank record (NC_045512.2), where the “surface glycoprotein” CDS (21563– 25,384) can be downloaded in FASTA format. Alternatively, use NCBI (Gene ID: 43,740,568 or NC_045512.2) to retrieve the CDS via genome annotation tools. Both full-length CDS and MEVC-based sequences are used in parallel.
The following step involves codon optimization of the CDS or MEVC to enhance translational efficiency. Use the LinearDesign algorithm (https://github.com/LinearDesignSoftware/LinearDesign, command-line version) to calculate CAI and MFE, running the CDS with λ = 3, which favors translation efficiency while maintaining structural stability (λ = 0 prioritizes structure; λ = 1–2 balances both; λ >3 biases codon usage). Optimize the codon usage to match the host organism’s preference, by selecting human as a host. After optimization, append UTRs and evaluate the full mRNA using RNAfold (http://rna.tbi.univie.ac.at//cgi-bin/RNAWebSuite/RNAfold.cgi) to assess MFE and secondary structure stability. In this protocol, LinearDesign was applied to the spike CDS (YP_009724390.1; Gene ID 43,740,568) and subsequently to the MEVC. Alternatively, using NovoPro ExpOptimizer tool (https://www.novoprolabs.com/tools/codon-optimization), input the MEVC amino acid sequence, select protein as input type, and choose Human and select restriction sites to be avoided (e.g., BamHI, XhoI) to prevent internal cleavage during cloning. Multiple runs may be performed to obtain the most optimized sequence. NovoPro provides CAI and GC content but not MFE; therefore, RNAfold was used for structural evaluation. Both CDS and MEVC were optimized, with emphasis on MEVC for cloning.
⚠ CRITICAL STEP: Target accessibility assessment ensures that key mRNA regions, especially around ribosome binding sites and epitopes, remain structurally accessible for translation and immune recognition. Evaluate optimized CDS and MEVC sequences using ViennaRNA tools (http://rna.tbi.univie.ac.at/): RNAfold for MFE secondary structure and accessibility, RNAeval for thermodynamic analysis, and RNAplfold for local base-pairing probabilities. As optional validation, analyze the MEVC at the protein level using PSIPRED (http://bioinf.cs.ucl.ac.uk/psipred/) for secondary structure prediction, perform molecular docking with immune receptors (e.g., TLR4) using ClusPro (https://cluspro.org/) with visualization in PyMOL (https://pymol.org/), and assess complex stability via molecular dynamics simulations using GROMACS (http://www.gromacs.org/). In this protocol, an additional step was performed for 2D structure analysis using PSIPRED to predict the coil structure, alpha-helix, and beta-sheet of the optimized MEVC. To perform that, open PSIPRED, choose to predict secondary structure, translate the optimized MEVC DNA sequence (resulted by NovoPro) into protein by Expasy translate tool (https://web.expasy.org/translate/) to input the optimized MEVC sequence and run the analysis.
Step 5: In silico Validation
To verify cloning and expression feasibility of the optimized mRNA construct perform in silico cloning, starting with selecting a suitable expression vector such as pET-28a(+), and download its sequence from Addgene (https://www.addgene.org), in SnapGene format. Import the vector into SnapGene (https://www.snapgene.com/) and review its annotated features. Open the optimized insert MEVC separately and perform restriction analysis on both vector and insert. Select BamHI (GGATCC) and XhoI (CTCGAG) based on their presence in the MCS, ability to enable directional cloning, and absence of internal sites in the insert. Add BamHI to the 5′ end and XhoI to the 3′ end of the CDS. Simulate digestion of both vector and insert with BamHI and XhoI to generate compatible ends. To perform this, open the vector file in SnapGene and select “Actions → Restriction and Insertion Cloning → Insert Fragment” to ligate the modified MEVC insert. Cut both the vector and insert with BamHI and XhoI then flip the fragment orientation if required to ensure correct directional insertion, as indicated by SnapGene orientation warnings. Validate the construct by confirming ORF continuity, correct start codon positioning, absence of internal stop codons, and the expected total plasmid size (~5.3 kb + insert).
For Immune simulation validation, evaluate the immunogenic potential of the designed MEVC using the C-IMMSIM server (https://kraken.iac.rm.cnr.it/C-IMMSIM/). Input the final amino acid sequence into the submission field then run the simulation using default parameters unless specific conditions are required. Define simulation steps to represent vaccine doses as single or multiple injections by defining time intervals corresponding to primary and booster immunizations. In this protocol, a single-dose injection model was used. Configure simulation volume and random seed if needed, then run the analysis. The platform uses PSSM and machine learning to model immune responses, including antigen processing and immune cell activation. Analyze outputs such as immunoglobulin levels (IgM, IgG), B- and T-cell populations, cytokine profiles (as IFN-γ, IL-2), and memory cell formation. A strong response is indicated by elevated antibody levels, robust T-cell activation, and sustained memory development.
⚠ CRITICAL STEP: Following the main workflow, optional UTR and CDS/MEVC optimization can be performed for experimental readiness. This involves testing combinations of 5′ and 3′ UTRs with the CDS or MEVC (total = number of 5′ UTRs × number of 3′ UTRs). UTRs are retrieved from UTRdb (https://utrdb.cloud.ba.infn.it/utrdb/index_107.html) or Ensembl (https://www.ensembl.org), combined with the CDS/MEVC (5′ UTR + CDS/MEVC + 3′ UTR), and analyzed using RNAfold (http://rna.tbi.univie.ac.at/cgi-bin/RNAWebSuite/RNAfold.cgi) to obtain MFE, secondary structure, and base-pair probabilities, where lower MFE indicates higher stability. For MEVC, reverse translation (https://www.bioinformatics.org/sms2/rev_trans.html) and mRNA transcription using Biomodel (https://biomodel.uah.es/en/lab/cybertory/analysis/trans.htm) are performed before UTR integration. For human expression, validated UTR pairs such as TMSB10 (5′) with APOA2 or AP3B1 (3′) are prioritized. An additional optional step evaluates nucleotide modifications, such as pseudouridine, using RNAfold energy corrections to assess their effects on structure and stability.
Results and discussion
Antigenicity assessment of the full-length CDS demonstrated a high antigenicity score of 0.47, while the MEVC exhibited a markedly higher antigenicity score of 1.70, supporting the effectiveness of multi-epitope vaccine design strategies in enhancing predicted immunogenic potential [104]. Toxicity screening of the CDS identified only three short peptide fragments as potentially toxic, namely VCGPKKSTNLVKNKC, CGPKKSTNLVKNKCV, and GPKKSTNLVKNKCVN; however, these fragments were already excluded from the final construct of MEVC, reflecting an improved safety profile. Physicochemical analysis of the MEVC (309 amino acids; 33.74 kDa) indicated favorable properties, including stability (instability index: 26.94 and >30 h half-life in human cells), thermostability (aliphatic index: 74.24), and hydrophilicity (GRAVY: −0.146, indicating a hydrophilic nature and enhanced solubility), along with suitable expression stability across systems.
Population coverage analysis demonstrated near-universal representation across diverse HLA alleles, with multiple epitope–HLA interactions per individual. Specifically, population coverage reached 99.82%, with an average hit of 5.18 and a PC90 value of 3.39, indicating that at least 90% of the population can recognize approximately 3–4 epitopes through HLA combinations, providing redundancy in immune recognition and reducing the likelihood of immune escape [105]. The distribution of 4–6 recognized epitopes per individual and a cumulative coverage approaching 100% further support robust and inclusive immune coverage (Fig. 6A). The findings support the suitability of the selected epitopes for multi-epitope vaccine design, as they maximize immune coverage across diverse HLA genotypes with redundancy in antigen recognition.
The time required to execute this procedure for a target protein such as the SARS-CoV-2 spike protein (P0DTC2), was approximately 2–3 hours, depending on server availability and computational load. In contrast, analysis of unknown targets typically requires additional computational steps and longer processing time. Therefore, the development of automated software solutions would be valuable to streamline and optimize mRNA vaccine design workflows, particularly for emerging pathogens where rapid response is critical.
Population coverage and conformational B-cell epitope mapping. () global world population coverage analysis of the selected epitopes based on HLA binding, showing high worldwide coverage, average hit, and PC90 values, indicating broad population inclusivity and the ability of most individuals to recognize multiple epitope–HLA combinations. Created with. () structural mapping of predicted conformational B-cell epitopes on the SARS-CoV-2 spike protein, where yellow highlights indicate epitope residues. The first panel corresponds to the top-ranked epitope spanning chains A and C (Supplementary table), whereas the remaining panels represent selected epitopes mapped on chains A, B, and C, respectively. Created with A B IEDB/Population IEDB/Ellipro S1
Design of the multi-epitope vaccine construct (MEVC). () schematic representation of the arrangement of adjuvant (HBD2), linkers, and selected HTL, CTL, and B-cell epitopes. () final amino acid sequence of the MEVC showing the concatenated epitopes and linker regions. Color coding: adjuvant (orange), HTL epitopes (blue), CTL epitopes (green), B-cell epitopes (yellow), and linkers (red). Created withand A B drawio.com MS word
RNA secondary structure analysis of optimized MEVC and full-length CDS. () predicted secondary structures of the optimized MEVC, where the minimum free energy (MFE) structure is shown on the left and the centroid structure on the right. () predicted secondary structures of the optimized full-length CDS, with the MFE structure on the left and the centroid structure on the right, illustrating increased structural complexity due to the longer sequence length. In both () and (), structures are colored according to base-pairing probabilities (0–1 scale), where red indicates highly stable paired regions and blue/green indicates lower pairing probability and higher accessibility. () mountain plot representation of the MEVC showing MFE (red), partition function (green), and centroid (blue) profiles, indicating stable folding behavior with a balance between structured and accessible regions. () mountain plot of the full-length CDS displaying similar trends with greater structural variation, reflecting its extended length. Created with A B A B C D RNAfold
Secondary structure (2D) prediction of the optimized MEVC. The figure illustrates the predicted secondary structural elements along the protein sequence, including α-helices (pink), β-sheets (yellow), and coil regions (gray). The upper panel represents the confidence score for each amino acid residue, where darker shades indicate higher prediction confidence. The distribution of structural elements shows the presence of both ordered regions (α-helices and β-sheets) and flexible coil regions, supporting the structural stability and functional flexibility of the optimized MEVC. Created with PSIPRED
In silico cloning of the MEVC into the pET-28a(+) vector. () in silico cloning of the adapted vaccine sequence into the pET-28a(+) expression vector, showing the selected cloning region highlighted in red between XhoI (158) at the N-terminal and BamHI (1091) at the C-terminal, while the vector backbone is represented in black. () cloning construction workflow illustrating the insertion of the optimized MEVC fragment (939 bp) into the pET-28a(+) vector (5369 bp), resulting in a recombinant plasmid of 6262 bp. Directional cloning using BamHI (C-terminal) and XhoI (N-terminal) ensures correct orientation of the insert within the vector. Created with A B SnapGene
Immune simulation results of the optimized designed MEVC. () B-cell and plasma cell dynamics, including total B-cell count, memory B cells, immunoglobulin isotypes (IgM, IgG1, IgG2), and B-cell population states (active, presenting, duplicating, and anergic), illustrating activation, proliferation, differentiation, and memory formation following antigen exposure. () cellular immune response and innate immune dynamics, including total T-cell populations (helper and cytotoxic T cells), T-cell memory subsets, natural killer (NK) cells, dendritic cells (DCs), macrophages, and epithelial cells, demonstrating coordinated immune activation, antigen presentation, and early innate responses. () antigen and antibody profile, showing rapid antigen (ag) clearance accompanied by a typical primary immune response with IgM production followed by class switching to IgG subclasses (IgG1 and IgG2), indicating effective humoral immunity. () cytokine and interleukin response, highlighting elevated levels of IFN-γ, IL-2, and other key cytokines (e.g., IL-4, IL-6, TNF-α), reflecting strong immune activation and regulation of cellular immune responses. Created with A B C D C-IMMSIM
Conclusion
This protocol establishes a comprehensive computational framework for mRNA vaccine design by integrating immunoinformatics, structural biology, and molecular optimization into a unified workflow. The pipeline ensures efficient antigen selection, stable mRNA design, and robust immune activation. The protocol defines the design space, key selection criteria, and stepwise procedures required for the development of both full-length CDS and MEVC, which represent a reproducible and adaptable strategy applicable to different pathogens, particularly in scenarios requiring rapid vaccine development. The application of this framework to the SARS-CoV-2 spike protein demonstrated its ability to generate constructs with high antigenicity, broad population coverage, optimized codon usage, and stable RNA secondary structures, supported by in silico cloning and immune simulation validation analyses of the MEVC. These results are consistent with experimentally validated findings and highlight the effectiveness of combining epitope-driven design with sequence-level optimization strategies. While computational predictions indicate strong immunogenic potential and translational feasibility, further in vitro and in vivo studies are necessary to validate the safety, efficacy, and protective performance of the designed construct.
Electronic supplementary material
Below is the link to the electronic supplementary material.
Supplementary material 1
Supplementary material 2