一种基于mRNA的新型狂犬病病毒多表位疫苗：通过反向疫苗学和免疫信息学进行计算机设计

序列检索和共有序列生成

狂犬病病毒蛋白（G、P、NP、M和L）的序列从美国国家生物技术信息中心数据库（NCBI）下载，然后通过Clustal omega V1.2.4（https://www.ebi.ac.uk/jdispatcher/msa/clustalo）与参考株序列进行多序列比对。随后，通过EMBOSS Cons CLI V6.6.0工具（https://www.ebi.ac.uk/Tools/）为每种狂犬病病毒蛋白生成共有序列。这些目标序列被用作预测和表位图谱分析的参考，并通过VaxiJen v2.0服务器（http://www.ddg-pharmfac.net/vaxijen/VaxiJen/VaxiJen.html，访问时间为2024年2月17日）评估其抗原性。该服务器系统通过不依赖比对的预测方法来预测蛋白质的抗原性，临界值大于0.45。VaxiJen服务器根据氨基酸的理化性质来显示蛋白质的抗原性质。该服务器的抗原预测方法基于将蛋白质序列通过自交叉协方差（ACC）转换为具有氨基酸特征的统一向量40。

线性B细胞表位的预测与分析

为提高表位预测的准确性，通过多个组合服务器分别预测了亚单位疫苗候选物的线性B细胞表位。

BCPREDS服务器（http://ailab.cs.iastate.edu/bcpreds/）采用了三种已开发的方法，即AAP、BCPred和FBCPred，来预测B细胞表位。BCpred使用后续的基于核函数的支持向量机（SVM）分类器，其准确率AUC值为0.758，截断分数为0.8，用于预测线性B细胞表位41。而ABCpred（http://crdd.osdd.net/raghava/abcpred/）则应用递归人工神经网络（ANN）方法对表位和非表位进行分类42。其阈值设定为0.8，相应的窗口长度为16聚体。

BepiPred2.0（http://services.healthtech.dtu.dk/service.php?BepiPred-2.0）采用随机森林算法，该算法基于从抗体-抗原蛋白质结构中注释的表位进行训练。BepiPred-2.0利用基于序列的表位预测方法，既对来自3D结构的表位数据进行预测，也对从IEDB数据库下载的大量线性表位集合进行预测43。BcePred（http://crdd.osdd.net/raghava/bcepred/）以蛋白质的一级序列作为输入，根据氨基酸的物理化学性质（如亲水性＞1.9、灵活性＞2.0、可及性极性＞1.8）以及暴露表面来预测B细胞表位44。

Epitope 1D（http://biosig.lab.uq.edu.au/epitope1d/）采用随机森林算法，同时利用基于图的蛋白质序列特征表示，并整合了生物本体论信息45。

所有工具均设置了默认参数，选择至少三个软件程序预测的顶级表位中部分或完全重叠的序列作为B细胞优势表位。使用计算工具评估B细胞表位的抗原性、过敏原性、毒性和免疫原性。使用VaxiJen v2.0预测抗原性，而AllerTop v2.0（http://www.ddg-pharmfac.net/AllerTOP）通过将蛋白质序列转换为固定长度的统一氨基酸组成（ACC）向量来评估过敏原性，并采用-0.4作为特异性阈值46。

ToxinPred（http://webs.iiitd.edu.in/raghava/toxinpred）通过针对1805种毒性肽的数据集分析理化性质来确定毒性47。最后，IEDB的免疫原性工具（http://tools.iedb.org/immunogenicity/）基于理化特征（包括侧链组成和氨基酸定位）预测免疫原性潜力48。

细胞毒性CD8+T淋巴细胞（CTL）表位预测

NetCTL1.2服务器（http://www.cbs.dtu.dk/services/NetCTL/）49用于预测G、P、N和M RV蛋白的9聚体长度CTL表位。针对1类12种MHC超型（A2、A3、A24、A26、B7、B8、B27、B39、B44、B58和B62）设定的阈值为0.75。该服务器采用默认值预测这些CTL表位，涉及三个特征的组合：TAP转运效率=0.15、C端切割=0.05以及MHC-I结合肽。

人工神经网络（ANNs）被用于预测主要组织相容性复合体（MHC）I类分子的结合以及蛋白酶体C端切割，而权重矩阵则被用于计算TAP转运蛋白的效能。神经网络预测MHC/肽结合的输出是一个与纳摩尔（nM）单位的IC50值相关的对数转换值。经确认与≥4种人类白细胞抗原（HLA）I类超型存在相互作用的表位会进行下游研究。

通过免疫表位共识（IEDB）（http://tools.iedb.org/mhci/）50方法进一步检测已识别的表位。在完整的HLA参考集中，百分位排名＜0.2且IC50值＜50 nM的肽段被视为强结合剂。

然后对预测出的细胞毒性T淋巴细胞表位进行了抗原性、致敏性和免疫原性筛选。只有那些在免疫原性方面显示出正值的表位被保留下来，用于下一阶段的评估。

辅助性T淋巴细胞表位预测

NetMHCIIpan–4.0（http://services.healthtech.dtu.dk/service.php?NetMHCIIpan-4.0）用于通过人工神经网络预测15聚体长度的肽与任何已知序列的MHC-II分子的结合情况51。在所有人群中最常见的HLA参考组分子DRB1∗01:01、DRB1∗03:01、DRB1∗04:01、DRB1∗07:01、DRB1∗08:01、DRB1∗11:01、DRB1∗13:01和DRB1∗15:01被选用于表位筛选。使用IEDB SMM方法对预测的表位进行排名，将IC50的临界值设为500 nM。根据抗原性、非过敏性、无毒性和免疫原性，对与≥3种超型等位基因相互作用的入围表位进行评估。

CTL和HTL表位与MHC等位基因之间的肽-蛋白对接

GalaxyWEB（https://galaxy.seoklab.org/）上提供的PepDock（https://galaxy.seoklab.org/cgi-bin/submit.cgi?type=PEPDOCK）52被用于CTL和HTL表位与其相关MHC等位基因的分子对接，以确定它们在细胞表面呈递以供T细胞识别的可能性。HLA等位基因的PDB结构，如HLA-A*01:01（PDB: 6MPP）、HLA-A*02:01（PDB: 2GIT）、HLA-A*02:06（PDB: 3OXR）、HLA-A*08:01（PDB: 7NVI）、HLA-B*15:01（PDB: 8ELG）、HLA-B*35:01（PDB: 4LNR）、HLA-DRB1*01:01（PDB: 5V4N）、HLA-DRB1*04:01（PDB: 5LAX）、HLA-DRB1*11:01（PDB: 6CPL）和HLA-DRB1*15:01（PDB: 5ELH），均从RCSB蛋白质数据库服务器（https://www.rcsb.org/）检索获得。使用ChimeraX v1.8去除了所有额外的肽、离子、水和配体。最后，利用HADDOCK服务器的蛋白质结合能预测工具PRODIGY（https://wenmr.science.uu.nl/prodigy/）53计算表位与等位基因之间的对接分数。该服务器包含以下计算参数：室温下复合物的结合亲和力（ΔG）和解离常数（Kd）。

人群覆盖率分析

将具有各自HLA等位基因的T细胞表位提交至IEDB人群覆盖率分析工具（http://tools.iedb.org/population/），以计算其在全球人群中的覆盖率。该工具用于根据HLA基因型频率以及与MHC结合或T细胞限制相关的数据，估算可能对特定表位集产生反应的个体比例54。

非人类同源表位筛选

为避免自身免疫的发生并筛选非人类同源表位，使用了BLASTp（https://blast.ncbi.nlm.nih.gov/Blast.cgiPAGE=Proteins）。选择在E值纳入阈值0.005以下无匹配项且一致性低于35%的表位作为非宿主/非同源蛋白。

多表位mRNA疫苗构建体设计

最终选定的LBL、HTL、CTL和表位通过KK、GPGPG和AAY连接子连接，以设计疫苗构建体。为了制造高效疫苗并增强免疫反应，具有强效TLR4激动剂作用的50S核糖体蛋白L7/L12佐剂序列被依次添加到入围的表位中。此外，在N端添加了组织型纤溶酶原激活剂（tPA）分泌信号序列（tPA），以帮助将疫苗转运出细胞；在C端添加了MHC I靶向结构域（MITD），用于将CTL表位引导至内质网的MHC-I区室。

mRNA疫苗构建体中还添加了5′端帽、Kozak序列、120-150个碱基的多聚（A）尾，以及3′和5′非翻译区（UTRs），这些对翻译和稳定性至关重要。

疫苗构建体的抗原性、免疫原性、变应原性和毒性预测

如前所述，通过在线服务器评估了所设计的mRNA候选疫苗的抗原性、过敏原性、毒性和免疫原性。

疫苗构建体的理化特性

通过ExPASy ProtParam服务器（http://us.expasy.org/tools/protparam.htm）55计算了该疫苗构建体的理化参数，包括理论等电点、分子量、正负残基总数、消光系数、不稳定性指数、脂肪族指数和总平均亲水性。

蛋白质溶解度预测

SOLpro服务器（http://scratch.proteomics.ics.uci.edu）用于预测疫苗构建体的溶解度。SOLpro通过基于一级序列多种表征的两阶段支持向量机（SVM）架构，预测蛋白质在大肠杆菌（Escherichia coli）中过表达时的可溶性倾向56。

二级结构预测

通过PSIPRED（4.0版本）（http://bioinf.cs.ucl.ac.uk/psipred/）预测了整个疫苗构建体的二级结构。PSIPRED是一种简单且准确的二级结构预测方法，它整合了两个前馈神经网络，对从PSI-BLAST（位置特异性迭代BLAST）获得的输出进行分析57。溶剂可及性（ACC）和无序区域（DISO）由RaptorX Property服务器预测58。

3D结构建模与优化

疫苗亚单位的三级结构通过transform restrained trRosetta服务器（https://yanglab.nankai.edu.cn/trRosetta/）59预测。这个基于网络的平台结合了深度学习和Rosetta，以实现快速且准确的蛋白质结构预测。通过采用受约束的Rosetta和直接能量最小化，以及深度神经网络的应用，成功构建了蛋白质结构，同时确定了残基的位置和方向。

为了提高从Rosetta服务器获得的3D模型的质量，使用了3Drefine网络服务器（http://sysbio.rnet.missouri.edu/3Drefine/）60。该优化协议基于两步流程：优化氢键网络和原子级能量最小化。因此，通过GDT-HA、GTD-TS、MolProbity、冲突和RMSD分数评估了改进后模型的质量。

疫苗三维结构验证

通过拉氏图以及从PROCHECK（http://saves.mbi.ucla.edu/）和ProSA（http://prosa.services.came.sbg.ac.at/prosa.php）服务器获得的Z分数，对优化后的3D结构模型可能存在的结构错误进行了验证。ERRAT用于检查非键合原子相互作用。

不连续（构象型）B细胞表位筛选

肽疫苗模型的精细化3D结构已提交至IEDB Ellipro工具（http://tools.iedb.org/ellipro/），以识别构象B细胞表位。ElliPro会基于相邻残基、形状（椭球体）估计和突出指数参数，为每个输出的表位打分64。

多表位亚单位疫苗与 toll 样受体（TLR4）的分子对接

使用ClusPro 2.0蛋白质-蛋白质对接在线服务器（http://cluspro.bu.edu/）65执行对接，以确定疫苗构建体与目标免疫细胞受体的结合亲和力。TLR4（PDB 3FXI）的三维结构从RCSB蛋白质数据库（RCSB PDB）（https://www.rcsb.org/）获取。该服务器通过运行PIPER对接程序进行刚体对接，能够基于不同的评分方案生成不同的模型。肽的结构通过PyMOL进行可视化66。具有最大簇大小和高结合亲和力（最低能量结构）的对接复合物进一步进行分子动力学模拟。

最后，使用PDBsum67对复合物中的蛋白质-蛋白质相互作用进行了详细分析。通过PRODIGY探究了对接复合物的结合强度。

疫苗-TLR4复合物分子动力学（MD）模拟

为了评估对接的TLR4-疫苗复合物的动态行为和稳定性，使用GROMACS 2023.1进行了明确的全原子分子动力学模拟（MDS）68，该软件编译时支持CUDA以利用GPU加速。模拟在突尼斯巴斯德研究所的高性能计算集群（称为omics）上运行，该集群配备了NVIDIA GPU和200个CPU核心，能够实现高效的并行计算69。

复合物的初始3D结构首先使用pdb2gmx模块结合CHARMM36力场进行处理，选择适合生理条件的质子化状态，并采用SPC/E水模型进行溶剂化。处理后的结构被置于立方模拟盒中心，与盒壁保持1.0 nm的间距，随后使用预定义的spc216水构型进行溶剂化。为了中和系统并模拟生理离子条件，使用genion工具通过替换溶剂分子添加了钠离子（Na+）和氯离子（Cl-）。采用最陡下降算法进行能量最小化，以消除空间位阻冲突并松弛系统。最小化后的结构随后经历两步平衡过程：首先在恒定体积和温度（NVT系综）下，然后在恒定压力和温度（NPT系综）下，每步均持续100 ps，同时约束重原子。生产分子动力学模拟在NPT系综中进行，总时长为200 ns，时间步长为2 fs，使用50个MPI线程，每个MPI进程配备4个OpenMP线程，并通过CUDA利用GPU加速。静电作用采用粒子网格埃瓦尔德（PME）方法处理，所有键长均使用LINCS算法约束以确保数值稳定性。轨迹文件每10 ps保存一次。模拟后分析用于评估系统稳定性和构象变化：计算主链原子（Cα）的均方根偏差（RMSD）以评估随时间的结构漂移；从能量曲线中提取平均温度以确保热稳定性；量化整个轨迹中主链原子间的氢键以评估分子内和分子间的内聚力；每10 ns提取代表性结构快照，以可视化复合物随时间变化的演变70。这一全面的方案确保了在明确溶剂条件下，对蛋白质-疫苗相互作用在较长时间尺度上的准确建模。

TLR4-疫苗复合物的正常模式分析（NMA）

通过iMODS服务器（http://imods.chaconlab.org/）71对TLR4疫苗结构进行了正常模式分析（NMA），以通过结构变形分析、相互作用残基的特征值（使结构变形所需的能量）、B因子值评分（迁移率分布）、单个残基间的协方差图以及弹性网络模型来研究疫苗-TLR4复合物的结构动力学和灵活性。基本界面与用于粗粒度模型原子表征的CA选项一起使用。

计算机免疫模拟

C-IMMSIM服务器（https://kraken.iac.rm.cnr.it/C-IMMSIM/）72（访问时间为2024年8月17日）被用于进行计算机免疫模拟，以阐明其触发宿主免疫系统对所设计肽疫苗产生反应的潜力。该服务器同时模拟哺乳动物免疫细胞的三个部分——骨髓、胸腺和淋巴结，并利用从机器学习技术中得出的位置特异性评分矩阵（PSSMs）来预测免疫相互作用。

每四周注射一次，共三剂，每剂含有1000个抗原蛋白且不含脂多糖（LPS），所有模拟参数均设为默认值，时间步长分别设为1、84和168个模拟步。一个时间步相当于日常生活中的8小时，第一剂在时间=0时注射。免疫模拟的输入参数包括10微升的标准体积、12345个随机种子和1050个模拟步。

mRNA疫苗的密码子优化与二级结构预测

通过GenScript（http://www.genscript.com/tools/gensmart-codon-optimization）（访问日期：2024年8月17日）对疫苗构建体的密码子进行优化，以优化该疫苗在人类宿主中的密码子使用情况。利用GenScript稀有密码子分析工具（http://www.genscript.com/tools/rare-codon-analysis）测定优化后DNA序列的GC含量、密码子适应指数（CAI）和密码子频率分布（CFD）。在宿主中实现充分表达的理想GC含量百分比范围为30%至70%。CAI是衡量目标序列同义密码子使用偏好的指标，其范围为0至1，数值高于临界值（>0.8）表明相关基因更有可能高效表达。密码子频率分布（CFD）用于表示异常串联密码子的出现情况。随后，通过翻译工具将优化后的DNA序列转换为RNA序列，并借助RNAFold网络服务器（http://rna.tbi.univie.ac.at/cgi-bin/RNAWebSuite/RNAfold.cgi）73进行二级结构预测，以开展热力学分析并获取最小自由能评分。

蛋白质组筛选和蛋白质优先级排序

从NCBI获取狂犬病病毒（RV）的G、P、NP、M和L蛋白的可用序列，并进行多序列比对。EMBOSS软件为每种蛋白生成一个共有序列，然后结合巴斯德株的相应参考序列对保守性进行分析。

结果显示，这五种蛋白质中有四种（G、P、N和M）具有抗原性（表1）。

这些抗原被选为疫苗靶点，用于针对狂犬病病毒（RV）的mRNA疫苗开发的计算机模拟分析。

线性B细胞表位预测

使用ABCPred、BCEpred、BCPREDS、Bepipred和epitope 1D服务器进行线性B细胞预测（图2，表S1）。保留至少在三个高分服务器中重叠且满足抗原性、毒性、过敏性和免疫原性标准的预测表位，用于最终疫苗构建。所有测试的蛋白质组分共预测出13个表位（表2）。

细胞毒性T淋巴细胞（CTL）表位预测

在选定的四种狂犬病病毒（RV）蛋白中，共有13个长度为9个氨基酸的细胞毒性T淋巴细胞（CTL）表位与至少四种主要组织相容性复合体（MHC）I类超型表现出强结合力（表S2）。随后，对这些表位的抗原性、免疫原性、非过敏性和无毒性进行了评估。最终，有4个表位（1个来自N蛋白，2个来自G蛋白，1个来自M蛋白）被预测为具有抗原性、免疫原性、非过敏性、无毒性的CTL表位（表3）。随后，通过联合方法确认了免疫表位数据库（IEDB）服务器，其中所选表位被认定为有效的CTL表位。

辅助性T淋巴细胞（HTL）表位预测

根据NetMHCIIpan 4.0版本的预测结果，在抗原蛋白中总共鉴定出14个MHC-II表位（15聚体），它们与至少3种独特的HLA DRB等位基因具有结合亲和力（表S3）。经过对强结合亲和力、抗原性、非致敏性和无毒性的全面评估，最终选择了G蛋白的2个表位，以及NP蛋白和M蛋白各1个表位用于疫苗设计（表4）。这些HLA表位中的大多数都显示出能够诱导干扰素-γ（IFN-γ）和白细胞介素-10（IL-10）等细胞因子的能力。

CTL和HTL表位与HLA等位基因的分子对接

我们使用GalaxyPepDock服务器对CTL和HTL表位与其相关的MHC等位基因进行分子对接，结果显示这些表位与其各自的HLA等位基因之间存在强烈的结合相互作用和热力学稳定性，这体现在负的吉布斯自由能值（ΔG <−7.6千卡·摩尔⁻¹，平均值=−17.5）和较低的解离常数上。图3显示，表位对接在其相应HLA等位基因的结合槽内。

人群覆盖率分析

入围的MHC-I和MHC-II表位在全球范围内的综合人群覆盖率达到96.01%（图3），pc90平均值为1.5，每个HLA抗原的平均命中数为4.08。欧洲和北美的覆盖率最高，分别为98.42%和97.62%，其次是北非，为92.54%（图3）。在被确定为狂犬病影响最严重的非洲地区中，几内亚比绍、佛得角、摩洛哥、突尼斯和津巴布韦的选定表位人群覆盖率超过90%。此外，苏丹、喀麦隆、塞内加尔、南非、乌干达和赞比亚的人群覆盖率超过80%（图4和图S1）。

多表位mRNA候选疫苗序列的构建

我们mRNA疫苗的最终构建体是从N端到C端设计的，如图5所示。

疫苗的编码区由三组具有抗原性、无毒性且无致敏性的CTL、HTL和B细胞候选表位组成，这些表位分别通过AAY、GPGPG和KK连接子依次连接。为增强免疫应答，通过GPGPG连接子添加了佐剂50核糖体蛋白L7/L12。关键的mRNA元件，包括5′m7G帽结构序列、包含AUG密码子的Kozak序列、位于mRNA开放阅读框（ORF）两侧的珠蛋白5′和3′非翻译区（UTRs）以及120-150个碱基的多聚（A）尾，均已整合到疫苗中，以提高翻译效率。此外，为促进翻译后表位的释放，通过GPGPG连接子在N端添加了组织型纤溶酶原激活剂（tPA）分泌信号序列（tPA）。另外，在C端通过EAAK连接子连接了一个MHC I靶向结构域（MITD），以引导CTL表位进入内质网的MHC-I区室。

抗原性、过敏性、毒性和免疫原性的评估

进一步分析发现，该嵌合构建体具有抗原性，VaxiJen v2.0和ANTIGENpro服务器预测的概率得分分别为0.6006和0.869。此外，该多表位疫苗还被检测为非过敏原、无毒、可溶（0.828496）且具有免疫原性（0.29221）。

三级结构建模与优化

通过trRosetta服务器成功模拟出了所设计的多表位疫苗的5种三级结构。预估的TM评分有助于了解两种结构之间的相似性，其值范围为0到1。TM评分高于-0.5表明拓扑结构正确，而评分低于-0.17则意味着存在随机相似性。所选模型的预估TM评分为-0.331。

优化结果以3Drefine分数呈现。保留的模型整体质量优异，全局距离测试高精度（GDT-HA）达到0.9686（[0,1]范围内的分数表明与初始模型的分数相近），RMSD为0.366埃，即优化模型相对于初始模型的Cα原子的均方根偏差分数（Cα-RMSD），molProbity值较低，为1.277，冲突分数为5.2，旋转异构体分数不佳，为0.0。图7A和B展示了保留的优化模型及其序列。

疫苗三维结构验证

通过ProSA-web、RAMPAGE和ERRAT服务器对优化后的三级结构进行了验证。优化后的模型显示，93.2%、5.4%、0.7%和0.7%的残基分别位于最有利区域、有利区域、允许区域、宽松允许区域和异常区域（图7C）。此外，ProSA的Z评分为-5.39（图7D），处于可接受的评分范围内。经ERRAT分析，局部模型质量为94.72（图7E）。所有这些结果表明，该构建体的预测结构是合理且可靠的。

构象B细胞表位的筛选

使用ElliPro预测经过验证的不连续（构象）B细胞表位的3D结构。总共有8个BCL构象表位涉及242个残基，得分范围为0.529至0.871（图8和表S4）。

疫苗蛋白与TLR4的分子对接

疫苗构建体与TLR4之间的分子对接通过ClusPro 2.0服务器完成，该服务器生成了10个对接复合物，并按能量和簇成员进行排名（表S5）。簇7的能量得分最低，为-1600.5千卡/摩尔，有21个簇成员，被选为最佳复合物（图9A）。

PDBsum相互作用分析显示，TLR4的B链与疫苗（A链）之间存在23个氢键和2个盐桥，同时还有302个非键接触（图9B）。

参与氢键形成的氨基酸及其键长列于表S6中。

使用PRODIGY网络服务器进行的高级分析对对接蛋白的吉布斯自由能（ΔG）和解离常数进行了评估，结果表明TLR4疫苗复合物在能量上是可行的，其ΔG值为-18.4千卡/摩尔（负值），解离常数较低，为3.4e-14。

分子动力学模拟

为了评估TLR4-疫苗复合物在分子动力学模拟过程中的整体结构稳定性，我们计算了200纳秒（ns）内主链原子的均方根偏差（RMSD）。如图10A和补充图S2所示，RMSD曲线呈现出以下特征。初始上升阶段（0-50 ns）：RMSD从约0.3 nm急剧上升至约2.5 nm，表明对接过程中施加的初始结构约束迅速松弛。该阶段通常对应于平衡和早期构象调整。中间波动阶段（50-150 ns）：RMSD值在1.5 nm至3.5 nm之间波动，表明复合物具有构象灵活性。这一动态区域可能涉及环运动或结构域重排，这在TLR4等免疫受体中较为常见。稳定阶段（150 ns后）：从约150 ns开始，RMSD值稳定在3.2-3.4 nm左右，表明系统达到了相对稳定的构象状态或局部能量最小值。这种稳定状态反映了模拟的收敛性，并为进一步的结构和能量分析提供了可靠的采样构象。较高的RMSD幅度（>2.0 nm）表明该复合物具有显著的结构可塑性，这可能与受体识别或信号转导的生物学过程相关。尽管存在这些波动，系统仍表现出逐步稳定的趋势，这是成功的长时间尺度分子动力学模拟的标志。

这种RMSD趋势强化了对接构象的质量，并证实TLR4-疫苗相互作用在200纳秒内保持结构稳定，为包括结合自由能计算、氢键稳定性和基本动力学在内的下游分析提供了可靠基础。

为了进一步验证在全局均方根偏差（RMSD）分析中观察到的TLR4-疫苗复合物的结构稳定性，我们对受体与疫苗构建体之间的相互作用界面进行了重点评估。如图10A和B所示，我们专门针对参与结合的界面残基计算了均方根偏差（RMSD）和均方根波动（RMSF）。在200纳秒的模拟过程中，界面的RMSD始终保持在较低水平（约0.12–0.15纳米），表明相互作用稳定，偏差极小。RMSF图谱显示蛋白质整体具有中等灵活性，数值范围在0.04至0.14纳米之间。在残基190（约0.14纳米）、300（约0.135纳米）和340（约0.145纳米）附近观察到较高的波动，这些位置对应于环区或末端区域。大多数残基的波动范围在0.04–0.11纳米内，表明局部灵活性有限，核心区域稳定。

此外，图10C和D分别展示了整个复合物的回转半径（Rg）和溶剂可及表面积（SASA）。回转半径（Rg）起始于约6.0纳米，并在25,000皮秒时稳步下降至约5.2纳米，在整个模拟过程中观察到短暂波动。在100,000皮秒时，其值大幅下降至约4.25纳米，表明系统发生了压缩。从100,000皮秒到150,000皮秒，Rg值在4.25纳米至4.5纳米之间略有波动，并且在模拟结束时降至4.25纳米以下，这表明 compactness有所增加。

溶剂可及表面积（SASA）从约637 nm²稳步下降至约620 nm²，表明其与溶剂的接触减少，这可能是由于三级结构的稳定和结合界面的维持。总之，这些结构描述符为TLR4-疫苗复合物在整体和分子界面上的构象稳定性及紧凑性提供了有力证据，支持其适用于结合自由能和基本动力学等下游分析。主成分分析（PCA）热图补充了自由能景观，突显了TLR4-疫苗复合物如何逐步稳定为窄范围的低能量构象集合。主成分1（PC1）的主导地位（解释最大方差）表明存在协同的结构重排，可能与功能性运动（如结合位点调整）相关。 为捕捉分子动力学（MD）模拟过程中TLR4-疫苗复合物内的主要集体运动，进行了主成分分析（PCA）。图10F显示了模拟轨迹在PC1和PC2上的投影，这两个主成分代表了主要的构象变化。该图展示了轨迹随时间的变化过程，颜色从深红色（模拟开始）渐变至灰色（模拟结束），体现了构象空间的持续且有限的探索。轨迹中多个区域被重新访问，表明在蛋白质-疫苗界面存在稳定和半稳定的构象，且分子间相互作用持续存在。这一观察结果与分子间氢键分析所揭示的模式一致，后者显示在整个模拟过程中氢键数量保持稳定（图10E）。 该图中的模拟快照表明复合物存在稳定且频繁出现的构象状态。相反，较高能量区域对应稳定性较差的瞬时构象，这说明TLR4-疫苗复合物探索的构象空间多样但有限，与稳定的蛋白质-配体系统特征相符，也印证了从RMSD、回转半径（Rg）和氢键分析等其他分析中得出的结论，这些结论共同证实了候选疫苗与TLR4受体之间相互作用的结构稳健性和生物学有利性。密集的红色簇反映了能量上有利且频繁出现的构象盆地，支持了自由能景观（FEL）分析的结果。这些稳定状态表明复合物具有紧凑、稳健的结构行为，与RMSD、Rg、SASA和氢键分析的结果一致。 总体而言，这种动态且有限的运动强调了疫苗构建体与TLR4受体之间具有生物学相关性且持续的相互作用，增强了疫苗构建体的结构完整性。总之，这些结果表明，在200 ns的分子动力学模拟过程中，TLR4-疫苗复合物保持了较强的结构稳定性、稳定的分子间相互作用和紧凑的折叠构象。这些发现验证了候选疫苗在生理条件下有效结合受体的潜力。

疫苗TLR4复合物的NMA

为了进一步分析疫苗与TLR4的复合物，通过iMODS服务器进行了正常模式分析（NMA）。iMODS采用正常模式分析方法，利用正常模式分析（NMA）探究大分子的集体运动。结果如图11所示。

对接复合物的正态模式分析（NMA）（图11A）。复合物的可变形性基于每个残基的特定不准确性，由链铰链表示。图中的峰值与蛋白质的变形区域相关（图11B）。B因子值评分证实了NMA结果（图11C）。与每个正态模式相关的特征值代表运动刚度。其值与变形结构所需的能量直接相关。特征值越低，变形越容易。如图11D所示，复合物的特征值为1.81729e–08，表明该复合物是稳定的。图11E显示了与每个正态模式相关的方差。彩色条显示了个体方差和累积方差。方差与特征值成反比。协方差分析如图11F所示，用于描述复合物分子动态区域中相关、不相关（白色）或反相关（蓝色）的原子运动。弹性网络模型定义了通过弹簧连接的原子对。图中的每个点（图11G）代表相应原子对之间的一个弹簧。这些点根据其刚度着色；颜色越深的灰色表示弹簧越硬。综上所述，这些结果预测了RAB-VAC TLR4复合物的稳定性。

计算机模拟免疫反应

使用C-ImmSim服务器研究了所设计疫苗的有效免疫应答特征。免疫模拟分析预测会产生强烈的体液免疫应答。初次免疫应答显示出高滴度的免疫球蛋白IgM。在后续暴露后，二次和三次免疫应答更为活跃，这可通过IgG1+IgG2和IgM水平的升高（相比IgG+IgM抗体水平较低）得到证明（图12A）。检测到高滴度的（IgG、IgM、IgG1和IgG1+IgG2）抗体，导致抗原浓度降低。此外，还检测到多种模拟的B细胞同种型和记忆B细胞的产生，这表明该疫苗构建体可诱导有效的长效免疫反应，并能更有效地清除抗原（图12B和C）。辅助性T细胞（TH）和细胞毒性T细胞（TC）的数量显著增加（暴露后约50天分别达到11000个/mm³和900个/mm³），并且伴随相应的记忆细胞发育，反应显著（图12D和E）。

Th0型免疫反应的细胞数量较低（少于20000个细胞/mm³），而Th1型免疫反应在第50天左右达到峰值，为120000个细胞/mm³（图12F）。事实上，经过一系列暴露后，Th1细胞达到了持续的显著水平（图12F）。

同时，暴露后自然杀伤细胞（NK细胞）和树突状细胞（DCs）的活性随之增强，同时巨噬细胞的增殖也得到改善（图12G、H和I），这表明其能向T细胞有效呈递抗原。在注射过程中，还可明显观察到干扰素-γ（IFN-γ）、转化生长因子-β（TGF-β）、白细胞介素-23（IL-23）、白细胞介素-10（IL-10）和白细胞介素-12（IL-12）等多种细胞因子的显著诱导（图12J）。根据免疫模拟结果，所设计的嵌合疫苗预计能引发强烈的免疫反应，并确保免疫原性和稳定性。

mRNA疫苗的密码子优化与二级结构预测

通过计算机程序对疫苗构建体的密码子使用情况进行了优化，以确保其在人类宿主中的表达效率。优化后的密码子序列长度为1788个核苷酸，CAI值为0.92（图13A）。改良构建体的平均GC含量和CFD值分别为58.33%和0%。总体而言，这些结果清楚地表明该构建体在目标表达宿主中具有高效表达的巨大潜力。

mRNA二级结构的预测至关重要，因为它是影响翻译起始、延伸以及mRNA生物合成的关键因素。通过RNAfold在线服务获得了与整个mRNA结构相关的自由能。RNA二级结构的最小自由能为ΔG = -634.70千卡/摩尔（图13B）。质心二级结构的最小自由能为-533.34千卡/摩尔（图13C）。最小自由能结构、RNA结构的热力学集合以及质心二级结构的山形图显示出显著的重叠（图13D）。每个位置的 positional entropy 也在图13E中呈现。这些分数反映了宿主内蛋白质翻译的稳定性和效率。mRNA稳定性的增强与表达率的提高直接相关。

局限性

虽然计算机设计的多表位狂犬病疫苗显示出巨大潜力，但它们也存在显著局限性。一个主要问题是表位预测算法的准确性，这些算法往往无法充分反映体内抗原加工和呈递的多层次复杂性。仅靠计算模型可能无法完美预测体内免疫动态，例如先天免疫和适应性免疫反应之间的相互作用、涉及的复杂细胞因子信号通路，以及持久免疫记忆的诱导——所有这些都是疫苗效力的关键因素。

所提出的嵌合疫苗模型在模拟中显示出较强的稳定性和免疫原性，但未经进一步测试，其对狂犬病感染的实际保护效力仍不确定。总体而言，未来的研究应优先开展体外和体内研究，包括动物攻毒实验，以评估其免疫原性和保护效力。

利益冲突

作者声明不存在利益冲突。

补充信息：在线版本包含可在以下网址获取的补充材料。