Immune system-based design and potential use of a universal cancer vaccine targeting RAS proteins

The entire details of pan-RAS mRNA cancer vaccine is given in the patent file number WO 2022/081764 A1 (PCT/US 2021/054859) (25). The FASTA sequence of SEQ ID NO 70 with following mutations G12D, G13D, L19F, A59T, G60D, Q61H, K117N and A146T is selected for the study. See Supplementary material. As described earlier (26), NetMHCpan v4.0 (27)1 is used to calculate the binding affinity of both mutated and matching wild-type peptides with HLA (MHC-I) alleles. For further analysis, a peptide was classified as a neoantigen only if the mutated form showed strong binding affinity (IC50 < 500 nM) while the corresponding wild-type peptide exhibited weak binding (IC50 > 500 nM), and such peptides were included in the analysis.

Cytotoxic T lymphocytes generally play a role in the containment of infections (28), and the prevalence of cytotoxic T cells usually correlates with the rate of pathogen clearance. Since these peptides are CD8 + T cell epitopes, they first need to be processed by the proteasome in antigen-presenting cells (APCs), and then transported by the Transporter Associated with Antigen Processing (TAP) system to HLA molecules. The resulting HLA–peptide complex is recognized by the T cell receptor (TCR), which activates CD8 + T cells. NetCTL v1.22 (20) is an immunoinformatics tool used to predict cytotoxic T lymphocyte (CTL) epitopes, which combines three factors: proteasomal cleavage, TAP transport efficiency, and MHC-I binding affinity. In NetCTL v1.2, MHC-I supertypes (A1, A2, A3, A24, A26, B7, B8, B27, B39, B44, B59, and B62) are included for predicting 9-mer CTL epitopes on neoantigens. The tool provides a combined score, and a cut-off value of 0.75 was used in this study. At this threshold, the sensitivity and specificity for CTL epitope prediction are 80 and 97%, respectively (29).

HTL epitopes are important for proper immune function because they release cytokines such as IFN-γ, IL-2, IL-4, and IL-10. Among them, IL-2 and IL-4 play key roles in activating CD8 + T cells. When HTL cells bind to the MHC II–peptide complex, they become activated and can either support CD8 + T cell–mediated immune responses (which kill tumor cells) or B cell–mediated antibody responses. Since CD8 + T cell–mediated immunity is central in cancer treatment, we focus our studies in identifying effective HTL epitopes. NetMHCpan v4.0 tool is used to predict the HTLs for the given vaccine sequence. The antigenicity, allergenicity and toxicity of the HTLs are again predicted by VaxiJen v2.0, AllerTOP v2.0 and ToxinPred respectively, though evaluating antigenicity is a prime step in vaccine design.

The Expasy ProtParam server7 is used to calculate several physicochemical properties of the multi-epitope vaccine, including amino acid composition, theoretical pI, molecular weight, instability index, half-life, aliphatic index, and the grand average of hydropathicity (GRAVY) (34).

The secondary structure of pan-RAS mRNA cancer vaccine is predicted by Prabi server.8 (35) Prabi server provides all the regions which produce (coiled, helical, and extended stranded) structures.

The Alphafold39 (36) server is used to predict the 3D structure of the pan-RAS mRNA cancer vaccine. Alphafold3 can generate highly accurate predictions of the structures of biomolecular complexes. To check the accuracy of the model, we performed validation using the ProSA-web server10 (37) and the SAVES v6.0 server.11 ProSA-web calculates a Z-score that reflects the overall quality of the protein model if the Z-score falls outside the normal range for native proteins; it suggests possible structural errors (37, 38). The PROCHECK tool within SAVES v6.0 further evaluates stereochemical quality by analysing residue geometry and the overall 3D structure of the protein (39, 40).

B lymphocytes are key players in the immune system because they produce antibodies, which help provide long-term immunity (41). B-cell epitopes are predicted to stimulate strong production of antigen-specific antibodies in the serum. Under natural conditions, antibody generation begins when immunogenic epitopes interact with B-cell receptors, triggering B-cell differentiation into plasma cells and memory cells. Plasma cells are primarily responsible for antibody secretion during the initial immune response, while memory cells contribute to rapid antibody production upon subsequent exposures. Therefore, in epitope-based vaccine design, a key goal is to identify and predict B-cell epitopes within the target protein. To identify B-cell epitopes, we predicted linear epitopes using the BCPREDS server.12 This tool uses a subsequence kernel-based SVM classifier with about 74.57% accuracy to predict linear B-cell epitopes (41–44). For discontinuous (conformational) epitopes, we used the ElliPro server.13 ElliPro applies residue clustering and Thornton’s method to predict these epitopes. Each predicted epitope is given a score called the PI (protrusion index) value, which shows how much the epitope stands out on the protein’s surface (45).

The solvent-accessible surface area (SASA) of the mRNA pan-RAS cancer vaccine is predicted using GETAREA, a tool from the FANTOM program (46).14 This tool measures different SASA parameters, such as the total surface area, solvation energy, and contributions from backbone and side chains, to check which residues are buried or exposed. Residues with high SASA values (>20) are more exposed to the solvent, which often suggests they may be more immunogenic.

In pancreatic cancer, Toll Like Receptor (TLR7/8) exhibits a complex dual role. TLR7 is frequently expressed in the tumor stroma and on pancreatic cancer cells, where it contributes to chronic inflammation and therapy resistance. Thus, therapeutic targeting of TLR7/8 requires careful contextual consideration to maximize its anti-tumor benefits while minimizing pro-tumor risks. The TLR7/8 agonists such as R848 can enhance anti-tumor immune responses and reduce cachexia in murine models, particularly when combined with radiotherapy. On the other hand, they may also drive pro-proliferative and pro-inflammatory effects, potentially promoting tumor progression in immunocompetent hosts. Strategies may include combination with modalities like Stereotactic Body Radiotherapy (SBRT) or the concurrent inhibition of tumor-promoting macrophages and exosomal microRNAs (22, 47).

The docking of TLR7/8 receptors with the mRNA panras cancer vaccine is carried out using Alphafold3 sever and then HDOCK server (45)15 is used for the docking of (TLR7/8-pan-RAS mRNA cancer vaccine-Resiquimod) complex. The AlphaFold3 model,16 with its advanced diffusion-based architecture can accurately predict the joint structure of complexes involving proteins, amino acids, and nucleic acids. The mRNA vaccine sequence docked with the relevant TLRs using the AlphaFold3 server was evaluated based on their pTM and ipTM scores. Subsequently, the docked (TLRs–mRNA vaccine complex) and resiquimod subjected to additional docking using the HDOCK server, employed a hybrid approach combining template-based modeling with free docking. This strategy allows cases affected by misleading templates to be corrected through the free docking protocol. The top-ranked docking model (TLR7/8-pan-RAS mRNA cancer vaccine-resiquimod) was then selected as an input for molecular dynamics (MD) simulations.

The TLR7 (5GMH) Organism(s): (Macaca mulatta), and TLR8 (5AWB) is taken from the protein data bank17 and then after removing the water molecules, ligands and ions, the FASTA sequence is again modelled with Alphafold3 server. TLR7 is a single-stranded RNA (ssRNA) sensor in innate immunity and also responds to guanosine and chemical ligands, such as imidazoquinoline compounds. However, TLR7 activation mechanism by these ligands remains largely unknown. TLR8, an intracellular immune receptor found in endosomes that recognizes ssRNA from viruses and certain bacteria, as well as synthetic compounds activates the immune system by promoting the production of cytokines and activating immune cells like monocytes and myeloid dendritic cells. TLR8 receptors play a key role in innate and adaptive immunity against viral infections and influence immune responses in bacterial infections and autoimmunity. In this study, TLR7/8 were considered as innate immune sensors involved in endosomal recognition of mRNA-derived RNA fragments, rather than as classical binding receptors for intact vaccine constructs. Consequently, docking of large protein or peptide constructs to the ectodomains of TLR7/8 does not recapitulate their established ligand recognition or activation mechanisms. Therefore, the presented docking results should be interpreted as an exploratory structural interaction analysis that provides qualitative insights into potential spatial compatibility or contact regions, rather than as mechanistic evidence of TLR7/8 activation.

All-atom classical MD simulations were performed for the (TLR7–pan-RAS mRNA cancer vaccine–R848) and (TLR8–pan-RAS mRNA cancer vaccine–R848) complexes using the AMBER20 (48) software package. The AMBER ff14SB force field was applied to parameterize the TLR7 and TLR8 protein structures, while the OL3 RNA force field was employed for the pan-RAS mRNA cancer vaccine to accurately describe RNA conformational dynamics. The small-molecule agonist R848 was parameterized using the General AMBER Force Field (GAFF), with ligand atom types and partial atomic charges assigned via the Antechamber module (49) using the AM1-BCC charge model. Each system was solvated in an explicit TIP3P water model within a octahedral periodic box, maintaining a minimum distance of 10 Å between the protein and the box boundary. System neutrality was achieved by adding 12 Cl⁻ ions to the TLR7 complex and 15 Na⁺ ions to the TLR8 complex, corresponding to the net charge of each system. Protonation states of ionizable residues were assigned to reflect physiological pH (7.4), with histidine protonation states manually inspected and assigned as appropriate.

After setup and parameterization of both systems, energy minimization was carried out in two stages to remove unfavorable steric interactions. Initially, 5,000 steps of minimization were performed on solvent molecules using steepest descent method while restraining solute heavy atoms using a harmonic potential. This was followed by 5,000 steps of unrestrained minimization of the entire system using conjugate gradient algorithms. The minimized systems were then gradually heated from 10 K to 310 K under the NVT ensemble over 50 ps, applying weak positional restraints of 4 kcal mol⁻¹ Å⁻² to complexes. Subsequently, density equilibration was conducted under the NPT ensemble at 310 K and 1 atm for 1 ns using the Langevin thermostat (50) and Berendsen barostat with a collision frequency of 2 ps⁻¹, followed by an additional 3 ns of unrestrained NPT equilibration to ensure proper stabilization. Production MD simulations were then performed for 500 ns under the NPT ensemble with a 2 fs time step. Long-range electrostatic interactions were treated using the Particle Mesh Ewald (PME) method with a real-space cutoff of 10 Å, while all covalent bonds involving hydrogen atoms were constrained using the SHAKE algorithm. Periodic boundary conditions were applied in all spatial directions throughout the simulations. Trajectory analyses, including assessments of structural stability, conformational dynamics, and intermolecular interactions, were performed using the CPPTRAJ module of AMBER20. All MD simulations were performed using the GPU code of PMEMD engine. MD production run has two replicas each system to ensure reproducibility and consistency. The open-source PyMOL software (51) was used for molecular visualization and preparation of publication-quality figures. ONIOM calculations were performed using Gaussian 16, Revision B.01 (52) with the AMBER force field. The optimization has been performed at B3LYP/6-311G*//AmberFF using Grimme Dispersion correction D3. The electronic embedding scheme was used while boundary was treated by employing hydrogen link atoms. Atomistic-level diagrams were generated using Chemcraft software (53).

Codon optimization for mRNA vaccine makes use of these preferred codons to improve protein production, so fewer copies of mRNA are needed. We used the Java Codon Adaptation Tool19 (55) to optimize the codons in our vaccine sequence.

To test how the vaccine might trigger an immune response, we used the C-ImmSim server20 (56), which has an agent-based algorithm incorporating a position-specific scoring matrix. The machine learning (ML) techniques is employed to estimate antigens and foreign particles affecting immune activity. Simulations are conducted using default parameters, except for the time steps, which were set at 1, 84, and 168 (with time step 1 representing the injection at time = 0, and each subsequent time step corresponding to 8 h). Accordingly, three injections are anticipated at 4-week intervals, consistent with the typical dosing schedule of most commercial vaccines. The Simpson’s diversity index (D) is calculated based on the resulting plot. This tool predicts how a mammalian immune system could react to a vaccine, including both antibody (humoral) and T-cell (cellular) responses (57–60). A three-dose vaccination schedule with 4 weeks between each dose is simulated. Most settings were kept at their default values: the number of adjuvants and antigen injections are set to 100 and 1,000, respectively (61). The simulation volume and steps were adjusted to 50 and 1,000, and the random seed was fixed at 12345 to ensure reproducibility. No lipopolysaccharides (LPS) were included in the setup.

The safety profile of vaccine with human proteins is performed with pBLAST21 (62) against human proteome (taxid: 9606) which is downloaded from the UniProt database. The similarity of vaccine with humans is estimated to avoid any autoimmune response (Homo Sapiens).

It has been observed that the TLR7/8 agonist R848 triggers anti-tumor activity and reduces cachexia in mouse models of PDAC. In vivo, tumors from two out of three tested cell lines responded to R848, leading to smaller tumor size, enhanced immune complexity, greater CD8 + T-cell infiltration and activity, and reduced Treg frequency. Mice treated with R848 also exhibited improvements in both behavioral and molecular features of cachexia, resulting in nearly twice the survival time. Knockout studies demonstrated that stromal but not tumor cell TLR7 is essential for R848’s effects. Analysis of patient samples further showed that TLR7 is consistently expressed in stromal cells throughout all stages of pancreatic neoplasia, whereas epithelial expression is relatively rare. Together, these findings suggest that immune-stimulating strategies such as R848 may hold promise for treating PDAC and cancer-associated cachexia. In combination therapies it is vital to predict the toxicity of the drug used in combination therapies for the cancer patients.

So the SMILES notation of resiquimod22 is taken from the drug bank and cardiotoxicity of drug is predicted with cardioToxCSM web server (61)23cardioToxCSM, which can predict six types of cardiac toxicity outcomes, including arrhythmia, cardiac failure, heart block, hERG toxicity, hypertension, and myocardial infarction, efficiently and accurately. Drug-induced liver injury (DILI) is a challenging endpoint in predictive toxicology because of the complex reactive metabolites that cause liver damage and the wide range of mechanisms involved in the development of the disease. The HepatoToxicity Portal24 (HTP) (63) is used to predict the liver toxicity of the drug. Embryotoxicity refers to the adverse and harmful effects of chemical, physical, or biological agents on an embryo’s development, leading to structural malformations, functional impairments, growth retardation, or even embryonic death. These effects can manifest in anatomical, physiological, or functional changes that are not always immediately apparent. Understanding embryotoxicity is crucial for identifying hazardous substances and implementing safety measures in various fields, including agriculture, pharmacology, and environmental protection. The embryoTox webserver25 (64) was trained and validated using in vitro bioactivity data of over 700 small molecules with teratogenicity effects. This simple and integrated platform provides a practical resource for optimising screening libraries to determine effective and safe molecules to use during pregnancy (65–70).

The FASTA sequence of SEQ ID NO 70 with following mutations G12D, G13D, L19F, A59T, G60D, Q61H, K117N and A146T is given in the. Supplementary material

Total 40 CTL epitopes are predicted, out of which only CVKIKKCII is toxin. One epitope VEDAFYTLY is immunogen while 39 epitopes are non-immunogen. 11 epitopes are non-allergenic whereas 29 are allergenic in nature. See Table 1.

Total 133 HTL epitopes are predicted and all are non-toxic. 81 HTL epitopes are non-allergenic and 69 HTL epitopes are non-antigens. See Table 2.

The predicted HTL epitope of pan-RAS mRNA cancer vaccine (WO 2022/081764 A1) (PCT/US 2021/054859) with the corresponding Alleles as predicted by NetMHCpan v4.0 tools is given in. Supplementary Table 2

The FASTA sequence with the above-mentioned mutations (mutations G12D, G13D, L19F, A59T, G60D, Q61H, K117N and A146T) is used to predict the biophysical and physicochemical properties of the pan-RAS mRNA cancer vaccine. Results showed that the total number of amino acids in pan-RAS mRNA cancer vaccine: 189, Molecular weight: 21,973.98 Da, and theoretical pI: 5.31. The molecular weight of pan-RAS mRNA cancer vaccine is lower than the standard MW (100 kDa). The formula predicted for the vaccine is C₉₇₀H₁₅₃₁N₂₅₉O₃₀₂S₁₀. Total number of negatively charged residues (Asp+Glu)-32; Total number of positively charged residues (Arg + Lys)-26 and the Total number of atoms is 3,072. The instability index (II) is computed to be 35.05, below the standard threshold of 40, indicating the stable nature of the protein. See Table 1. The Aliphatic index: 80.37 indicate the thermostable nature of the vaccine. The Grand average of hydropathicity (GRAVY) is (−0.496). The negative value GRAVY indicates the hydrophilic nature of the vaccine. The antigenicity of the vaccine, as predicted by the Vaxijen v2.0 server is 0.897175. As the value is > 0.4, pan-RAS mRNA cancer vaccine is probably an antigen. The solubility score of pan-RAS mRNA cancer vaccine is predicted by the SOLpro server, which is 0.733657. The solubility score above 0.5 (green color) indicates the soluble expression of the vaccine in Escherichia coli. AllerTOP v. 2.0 results indicated the non-allergenicity of the vaccine. (See Table 3).

The 2D structure of pan-RAS mRNA cancer vaccine is predicted by prabi server26 with all the coiled, helical, and extended stranded regions. Results indicate Alpha helix-71; 310 helix-0; Pi helix-0; Beta bridge-0, Extended strand-37; Beta turn-0; Beta region-0; Bend region-0; Random coil-81; Ambiguous state-0 and other states-0. Overall, the secondary structure of pan-RAS mRNA cancer vaccine has high globular conformability, stability and flexibility. See Figure 2.

The refined 3D structure of pan-RAS mRNA cancer vaccine as predicted by Alphafold3 server is given in Figure 3a. The pTM = 0.85 score showed stable structure of pan-RAS mRNA cancer vaccine. Ramachandran plots from the PROCHECK server are used for quality assessment of the 3D structure of the vaccine. Ramachandran plots analyse the conformational flexibility of proteins and the protein regions that causes structural disturbances. It can predict the structural abnormalities and conformations of protein loops. The Ramachandran plot predicts 91.4% core, 8. 6% allowed, 0.0% general, and 0.0% disallowed regions for the pan-RAS vaccine structure. All Ramachandrans: 89 labelled residues (out of 3,387) and Chi1-chi2 plots: 20 labelled residues (out of 1983) with better side-chain params. See Figure 3b. A Z-score of −7.32 indicates good quality of the pan-RAS mRNA cancer vaccine. The local quality shows the plotted graph of the residue scores. Negative values describe an absence of erroneous parts in the vaccine. See Figure 3c.

The Ellipro server is used to predict and analyse the immune response of conformational B-cell epitopes. The default parameters are – (minimum score- 0.5) and (maximum distance- 6 Å). The server predicted six new conformational B-cell epitopes with peptides IQNHFVDEYDPTIEDSYR, DHEEYSAM, FE, HHYREQIKRVKD, LPSRTVDTKQ and YFIKKISKEEKTPGCVKIKK with 18, 8, 2, 12, 10, 20 length. All the B-cell epitopes are non-toxin, non-immunogen with one peptide (HHYREQIKRVKD) as non-allergen. See Table 4.

Solvent Accessible Surface Area (SASA) analysis using GETAREA showed that the number of surface atoms is 911, with 630 as buried atoms and number of atoms with ASP = 0. Total area/energy is 10,911.86 Ų per kcal/mol. (Supplementary Table 1)

Docking of pan-RAS mRNA cancer vaccine with TLR7/8 is carried out with Alphafold3 sever. The ipTM = 0.21pTM = 0.77 is for (pan-RAS mRNA cancer vaccine -TLR7) complex and ipTM = 0.21pTM = 0.76 is observed for (pan-RAS mRNA cancer vaccine -TLR8) complex. The predicted template modeling (pTM) score and the interface predicted template modeling (ipTM) score are both based on the template modeling (TM) score, which measures how accurate a protein structure is (71, 72). A pTM score above 0.5 suggests that the overall fold of the predicted complex is likely close to the real structure. ipTM focuses on how well the relative positions of different subunits in the complex are predicted with Scores above 0.8 mean the prediction is very reliable; Scores below 0.6 usually mean the prediction failed and Scores between 0.6 and 0.8 fall into a gray area where the prediction could go either way. Since TM scoring is very strict for small proteins or short chains, pTM scores can be very low (<0.05) when fewer than 20 residues are involved. See Figure 4.

The docking structure of (TLR7/8–resiquimod–pan-RAS mRNA vaccine) as predicted by HDOCK server is given in Figure 5.

To elucidate the conformational stability and dynamic behaviour of TLR7 and TLR8 in complex with Resiquimod and the pan-RAS mRNA cancer vaccine, 500 ns all-atom molecular dynamics (MD) simulations were performed under explicit solvent conditions. The trajectories were analysed for structural stability, compactness, and residue-level flexibility to evaluate receptor adaptability and ligand-induced conformational changes. The root mean square deviation (RMSD) revealed stabilization throughout the production phase. The TLR7 complex displayed an average RMSD of approximately 1.6 Å, while the TLR8 complex showed a slightly higher deviation of ~1.9 Å, suggesting enhanced conformational flexibility upon ligand binding. The radius of gyration (RoG) demonstrated consistent compactness, with TLR7 maintaining an average RoG of ~20.1 Å and TLR8 exhibiting ~17.9 Å, indicative of stable folding and minimal expansion of the tertiary structure. The slightly lower RoG value for TLR8 implies a more compact conformation, correlating with its tighter ligand accommodation and stronger noncovalent interactions, consistent with ONIOM-derived electronic polarization results. The root mean square fluctuation (RMSF) analysis quantified the flexibility of individual residues within the receptor structures. Most residues showed minimal fluctuation (<2 Å), confirming stable secondary structure maintenance. Localized peaks in the RMSF profile of TLR8 correspond to surface loops near the ligand binding region, which are typically associated with allosteric adjustments facilitating ligand recognition. Similarly, minor fluctuations observed in TLR7 loops reflect adaptive flexibility for ligand-induced signalling activation. See Figure 6.

To further characterize the conformational dynamics, normal mode analysis (NMA) was performed to investigate the collective motions within the complexes. The deformability showed coordinate deviations within 0–1 Å, indicating limited flexibility confined to loop regions without perturbing the overall structural core. B-factor mobility profiles predicted atomic displacements within 0–1 Å, consistent with the MD-derived RMSF trends, confirming the robustness of both TLR7 and TLR8 during ligand engagement. See Figures 7 and 8 respectively. Additionally, the elastic network models depicted dense inter-residue connections, particularly around the binding interface, signifying strong cooperativity and structural coupling between the ligand and receptor residues. Collectively, these MD and NMA results indicate that the (TLR7/8–resiquimod–pan-RAS mRNA vaccine) complexes attain a stable and energetically favorable conformation, with restricted residue flexibility and robust intramolecular connectivity. The observed differences in RMSD and RoG between TLR7 and TLR8 highlight receptor-specific adaptations, wherein TLR8 exhibits greater compactness and electronic polarization, suggesting a stronger and more stable binding interaction with resiquimod and the pan-RAS mRNA vaccine components.

Binding Site Analysis of TLR7 and TLR8 with Resiquimod: The most populated (after Clustering of 500 ns simulations) snapshot were taken for structural visualization of the (TLR7–resiquimod) and (TLR8–resiquimod) complexes revealed distinct but complementary binding patterns within their respective active sites (Figure 9). In the TLR7–resiquimod complex, the ligand is deeply embedded within the hydrophobic cleft located between the leucine-rich repeat (LRR) domains. Key hydrogen bonding interactions were observed between the imidazoquinoline ring of resiquimod and polar residues such as K406, Y338, and Q328, which stabilize the ligand orientation. Additionally, hydrophobic contacts with residues V329, P409, and F325 contribute to the structural rigidity of the complex. The hydrogen bond network formed between the ligand’s hydroxyl and amine functional groups and the side chains of TLR7 residues enhances electrostatic complementarity and charge delocalization within the active site. The overall configuration reflects a stable and well-defined binding pocket conducive to receptor activation. The interaction geometry supports the observed electronic polarization and the smaller HOMO–LUMO gap (4.574 eV), confirming strong ligand–receptor coupling and structural integrity at the energy-minimized state.

In the TLR8–resiquimod complex, the ligand is situated within a broader yet deeper binding groove formed at the junction of the LRR domains. Here, the amine and carbonyl groups of Resiquimod engage in strong hydrogen bonds with F453, H315, and H319, while π–π stacking interactions with F453 and Y611 further stabilize the aromatic core of the ligand. Several nonpolar residues, including F28, and V317, form van der Waals contacts that enhance the overall hydrophobic stabilization of the complex. The compactness of the TLR8 binding cavity, as confirmed by the lower average RoG value (~17.9 Å), facilitates denser packing of residues and stronger intermolecular forces. This corresponds to a smaller electronic gap (4.023 eV), reflecting increased charge transfer and enhanced electronic communication between the ligand and the receptor. The optimized interactions observed in TLR8 suggest a higher binding affinity and stronger activation potential compared to TLR7, consistent with molecular dynamics findings and normal mode analysis showing greater structural rigidity and cooperative domain motions upon ligand binding.

Hydrogen Bond Interaction Analysis: To investigate the nature and persistence of weak interactions between R848 and TLR7/8, hydrogen bond analyses were performed over the 500 ns MD trajectories. The hydrogen bonding frequency, average donor–acceptor distances, and angles were calculated to assess the strength and stability of these noncovalent interactions (Figure 10). In the TLR7–resiquimod complex, several stable hydrogen bonds were observed throughout the simulation, primarily involving K406, Q328, Y338, P409, and V329. The strongest and most persistent hydrogen bonds were formed between the nitrogen atoms of resiquimod (R848_1545@N1) and the terminal amine hydrogens of K406 (HZ1–HZ3), with a cumulative occupancy of approximately 0.06 and average bond distances around 2.89 Å, and bond angles near 158°, indicating robust electrostatic stabilization. Additional hydrogen bonds between R848_1545@N2 and Q328@NE2 (occupancy 0.011, distance 2.93 Å) and Y338@OH (occupancy 0.0079, distance 2.85 Å) contribute to maintaining the ligand’s orientation within the hydrophobic cavity. Weak but consistent interactions were also detected with P409 and V329, reflecting the dynamic nature of secondary contacts that enhance overall binding affinity. See Figure 10.

In the TLR8–resiquimod complex, hydrogen bonding was more widely distributed across the receptor interface, suggesting a greater degree of adaptive interaction. The dominant contacts included those between R848_609@O1 and Y196@OH (distance 2.81 Å, angle 162°), as well as interactions between R848_609@N1 and R41@NH1 (distance 2.90 Å, angle 158°) and Q311@NE2 (distance 2.90 Å, angle 161°). These hydrogen bonds remained stable for a significant portion of the simulation, facilitating optimal positioning of the ligand within the binding cleft. Additional weaker contacts were observed with S287, S39, D30, and R474, characterized by average donor–acceptor distances between 2.78–2.93 Å, reinforcing structural stability through transient hydrogen networks. The slightly lower HOMO–LUMO gap (4.023 eV) and compactness of TLR8 correlate with this extensive hydrogen bonding network, highlighting enhanced charge delocalization and electronic communication between ligand and receptor residues. The average distance (Å) and average bond angle (°) for the donor (resiquimod) and receptors (TLR7/8–pan-RAS mRNA vaccine) complexes are given in, respectively. Supplementary Tables 3 and 4

Following the molecular dynamics simulations, ONIOM-based quantum mechanical/molecular mechanical (QM/MM) calculations were performed to gain deeper insights into the electronic structure and noncovalent interactions governing drug–receptor and vaccine–protein binding. These calculations were executed using Gaussian16 Rev. B.01, employing the B3LYP functional (73) with the LANL2DZ (74) basis set, augmented by Grimme’s D3 dispersion (75) correction to account for long-range van der Waals interactions. The QM region comprised the bound ligand and key residues within the binding pocket, while the surrounding protein matrix was treated at the molecular mechanics level using the AMBER force field parameters extracted from the equilibrated MD trajectory. The primary objective of these hybrid ONIOM calculations (76) was to elucidate the electronic structure, charge distribution, and weak interaction networks (including hydrogen bonding, π–π stacking, and electrostatic stabilization) that stabilize the drug within the active site. The HOMO–LUMO energy gap was computed to assess electronic stability and reactivity of the ligand in both TLR7 and TLR8 complexes. See Figure 11.

The results revealed a HOMO–LUMO gap of 4.574 eV for the TLR7–R848 complex and 4.023 eV for the TLR8–R848 complex, indicating enhanced electronic delocalization and stronger polarization interactions in the TLR8 binding pocket. The smaller energy gap for TLR8 suggests greater electronic flexibility, which may contribute to more effective receptor activation and drug binding affinity. Frequency calculations confirmed the absence of imaginary modes, validating the optimized geometries as true energy minima. The optimized structures exhibited distinct conformational stabilization driven by intramolecular hydrogen bonds between the drug’s heteroatoms and key amino acid residues, alongside π–cation interactions involving aromatic side chains. The inclusion of dispersion corrections was essential to accurately capture stacking interactions within the hydrophobic regions of the binding site, particularly near the imidazole and purine moieties of the ligand. The frontier molecular orbitals (HOMO and LUMO) demonstrated that electron density in the TLR7-bound R848 was primarily localized on the aromatic framework, whereas in the TLR8-bound complex, it extended toward the polar functional groups interacting with charged residues, enhancing polarization and binding strength. These observations are consistent with the electrostatic field variations observed during the MD simulations, further supporting the notion that electronic polarization and charge redistribution play a pivotal role in stabilizing the drug–receptor interface. The optimized geometries, frontier orbital visualizations, and key bond parameters are presented in Figure 11, illustrating the spatial orientation of the drug within both TLR7 and TLR8 active sites. These findings provide quantum-level insights into the structural basis of ligand recognition and stabilization, complementing the molecular dynamics results and supporting the observed immunoinformatics predictions.

The original GC content: 40.89% for the species: Homo Sapiens. The CAI-Value of the improved sequence: 0.95; GC-Content of the improved sequence: 61.72% for the species: Homo Sapiens. Similarly, the original GC content: 41.75% for the mus-musculus. The CAI-Value of the improved sequence: 0.72 and GC-Content of the improved sequence: 59.78% for mus musculus.

Immunoinformatics analysis using the CImmSim server predicted strong immune responses to the pan-Ras cancer vaccine, characterized by elevated levels of IgM, IgG1, and IgG3 antibodies, along with increased populations of B cells, helper T cells, and cytotoxic T cells. Immune response graphs also revealed heightened levels of cytokines, interleukins, and natural killer (NK) cells, indicating sustained adaptive immunity (Figure 12). Immune simulations were conducted to understand the dynamics of antibody production and cell-mediated immune responses following vaccination, mapping the initiation and progression of immune activity within the host. Simulation results from the C-ImmSim server indicated a robust induction of immune responses following administration of the vaccine. The study predicted that the primary immune response to the antigenic fragments increased progressively after each of the three vaccine doses, as evidenced by rising concentrations of multiple immunoglobulins. Secondary immune responses were also enhanced, demonstrating effective priming from the initial exposure. Elevated levels of IgM and IgG reflected a strong humoral response (Figure 12a). The B cell population exhibited an increase in memory B cells and IgM isotypes, accompanied by a reduction in antigen concentration (Figures 12b,c), indicating improved antigen clearance during subsequent exposures (Figures 12d,e).

Additionally, the vaccine elicited substantial activation of TH1 and TC cell populations, along with memory formation (Figure 12f). Increased numbers of dendritic cells and macrophages suggested efficient antigen presentation by these APCs (Figure 12g,h). The construct also stimulated production of multiple cytokines, including IFN-γ, IFN-β, IL-10, and IL-23, which are critical for initiating and sustaining immune responses (Figure 12i). These findings were consistent with the observed IFN-γ induction following vaccination. Furthermore, the low Simpson diversity index (D) indicated a broad and varied immune response, highlighting the vaccine’s capacity to generate diverse immunological activity (Figure 12i).

Blastp server showed alignment score with more than 100 domains, which suggest that the mRNA-encoded protein may have cross-reactivity potential, highlighting the need for careful epitope selection, immunogenicity testing, and safety evaluation before clinical application. If predicted epitopes are shared across multiple unrelated proteins, the immune system may not selectively target the intended tumor antigen, reducing therapeutic efficacy and increasing unintended reactivity. BLASTp hits against proteins from pathogens, toxins, or allergens may indicate the potential for allergic or toxic responses, depending on the nature of the aligned sequences (See). Supplementary Table 5

The expected Hepato-toxicity is 0.881633 and the cardiotoxicity predicts that it is safe for Cardiac failure, Heart Block, herG toxicity and Myocardial Infarction while it is unsafe for Arrhythmia and Hypertension. The embryotoxicity results predict its unsafe use during pregnancy. See Figure 13.

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