企业官网建设流程全解析

摘要：本文系统梳理 Amber 分子动力学模拟蛋白-核酸复合物（蛋白质-DNA / 蛋白质-RNA）的完整流程：从 PDB 结构与 pdb4amber 预处理入手，重点说明 ff19SB + OL15/OL21/OL3 + OPC + JC 离子的力场组合；随后演示 tleap 两段式 combine 构建体系的完整命令；接着给出强约束最小化、分段升温、弱约束平衡与生产模拟的可复用输入文件，并总结 pmemd 运行时的常见坑；最后以 cpptraj 实战为验收核心，介绍 nastruct 螺旋参数、糖环 pucker 等核酸特征分析方法，补充离子浓度、Mg²⁺、周期性边界等特殊注意事项及常见错误排查。全文以 Drew-Dickerson 十二聚体为官方教程基线，命令参数可直接照抄复用。

全文约 10 千字 ｜ 预计阅读 25 分钟

本文介绍Amber分子动力学模拟结构准备：蛋白-核酸复合物（蛋白质–DNA / 蛋白质–RNA）的分子动力学模拟操作及注意事项。文中命令与输入文件参数取自 Amber 官方核酸模拟教程（AMBER Hub，Drew-Dickerson 十二聚体案例）与近年力场评测文献，可直接照抄复用。

做转录因子–DNA、Cas9–sgRNA–DNA、核酸酶这类体系的模拟时，很多人直接套用纯蛋白的流程，结果跑出来 DNA 解旋、能量爆炸或者 RMSD 永不收敛。蛋白–核酸复合物与纯蛋白体系有三个本质差异，每一步操作都要围绕它们展开：

• 强负电体系：每两个核苷酸就有一个磷酸基团，一条 12 bp 双链 DNA 净电荷 -22，离子环境直接决定结构稳定性；
• 力场组合：蛋白力场 + 核酸力场必须同时加载，且核酸力场的版本选择比蛋白更讲究（选错会系统性破坏螺旋几何）；
• 构象验证：DNA 有 A/B/Z 三种构象，RNA 以 A-form 为主，糖环 pucker（C2′-endo vs C3′-endo）是否正确是判断模拟可信度的硬指标，纯蛋白流程里没有这一步。

相关教程与核心文献

官方教程

核心文献

一、准备工作

1. 获取初始结构

• 来源：PDB（X 射线或 Cryo-EM 结构），例如1ABC.pdb。
• 检查清单：
  • 是否包含完整蛋白和核酸链；
  • 是否有缺失残基/原子（可用Modeller、ChimeraX补全）；
  • 非标准残基（修饰碱基、甲基化碱基、甲基化氨基酸）需单独参数化（参见非标准残基处理指南）；
  • 删除无关水分子、配体（除非必要）。

✅ 建议保留结合位点附近的结晶水（介导蛋白–核酸识别的水桥），但后续统一用显式水模型重建溶剂环境。

2. pdb4amber 预处理

pdb4amber protein_dna.pdb --dry --reduce --overwrite

• --dry：移除水（可选）；
• --reduce：修复 HIS 命名（HID/HIE），并补氢；
• 自动标准化残基命名。

⚠️实测坑：部分 PDB 条目里核酸链整链是HETATM记录而不是ATOM。--dry只按记录类型删水，但如果后续手动按HETATM清理"杂原子"，会把整条 DNA 一起删掉。预处理后务必grep -c "DA\|DT\|DG\|DC" *.pdb确认核酸链还在。

二、选择合适的力场（核酸模拟的地基）

Amber 官方推荐组合（基于 Amber20/22 推荐实践与 2023 年力场评测）：

三个要点：

1. OL15 是什么。加载leaprc.DNA.OL15时终端会打印一行说明：OL15 force field for DNA (99bsc0-betaOL1-eps-zetaOL1-chiOL4)——即在 ff99 基础上叠加 bsc0 的 α/γ 修正、βOL1、ε/ζOL1、χOL4 四组骨架二面角修正。这些修正专门解决早期力场把 B-DNA 慢慢扭曲成非物理构象的问题。
2. 2023 年之后的最佳实践。Love 等（J. Chem. Theory Comput. 2023）对 Amber 核酸力场做系统 vdW 参数扫描后，推荐 dsDNA 用 OL21 + OPC 水模型。AmberTools 21 起提供leaprc.DNA.OL21；用 AmberTools 20 或更早版本则用 OL15。
3. 蛋白与水模型要配套。ff19SB 的 CMAP 项是配合 OPC 水优化标定的，官方教程（包括核酸教程）现在统一source leaprc.water.opc，OPC 与 JC 离子参数随之自动加载。TIP3P 仍可用（传统组合 ff14SB+TIP3P+JC），但不要让蛋白用 ff19SB+OPC、离子却手动塞 TIP3P 的参数。

⚠️ 不要混用旧力场（如 ff99 + parmbsc0 的原始参数），可能导致能量不一致与螺旋畸变。

三、使用 tleap 构建体系

3.1 纯核酸热身：Drew-Dickerson 十二聚体（官方教程实录）

Amber 官方核酸教程用 d(CGCGAATTCGCG)（Drew-Dickerson dodecamer，DDD）做示范——这是核酸模拟领域被研究最多的短双链，从 A-DNA 构象出发，观察模拟中自发转变为 B-DNA。完整 tleap 会话：

$ tleap > source leaprc.DNA.OL15 > source leaprc.water.opc > dna = loadpdb A-dna.pdb Loading PDB file: ./A-dna.pdb total atoms in file: 758 > check dna Checking 'dna'.... WARNING: The unperturbed charge of the unit: -22.000000 is not zero. check: Warnings: 1 Unit is OK. > addions dna Na+ 0 > solvateoct dna OPCBOX 9.0 > saveamberparm dna dna.prmtop dna.rst7 > quit

逐行解读，这些数字就是你的对照标准：

• 758 个原子：12 bp 双链 DNA 的重原子 + 末端残基标准处理后的规模；
• 净电荷 -22：22 个磷酸基团，check报这个 WARNING 是正常预期，不是错误——看到它说明磷酸骨架完整；如果电荷不是 -2×(链长-1) 的量级，说明缺了核苷酸；
• addions dna Na+ 0：加 Na⁺ 直到净电荷为 0（中和），离子优先放在负电位点附近；
• solvateoct dna OPCBOX 9.0： truncated octahedron 八面体盒子，溶质到盒边至少 9 Å，教程体系约加 5336 个水分子。八面体比立方盒省相当数量的水（省计算量），长核酸体系推荐；
• saveamberparm输出拓扑dna.prmtop+ 坐标dna.rst7。

3.2 蛋白 + B-DNA 复合物构建（两段式 combine）

蛋白–核酸复合物推荐把蛋白链与核酸链拆成独立 PDB、分别加载再 combine（GENESIS 教程用 PDB 3LNQ 演示的同一套思路，Amber 下通用）：

# 加载力场：蛋白 + DNA + 水，三行都要 source leaprc.protein.ff19SB source leaprc.DNA.OL15 source leaprc.water.opc # 若含 RNA：source leaprc.RNA.OL3 # 蛋白与核酸分别加载（残基命名须符合 Amber 规范） # 标准 DNA 残基名：DA5/DA/DA3, DG, DC, DT（5'/3' 端带编号） # 标准 RNA 残基名：A5/A/A3, G, C, U protein = loadpdb protein_clean.pdb dna = loadpdb dna_clean.pdb complex = combine { protein dna } # 检查警告（重点关注 "unknown residue" 或 "missing atoms"） check complex desc complex # 中和净电荷（核酸体系通常加大量 Na+） addIons complex Na+ 0 # 追加生理盐浓度（按盒子体积估算离子对数，如 150 mM NaCl） addIons complex Na+ 12 addIons complex Cl- 12 # 溶剂化：最小水层 10–12 Å solvateoct complex OPCBOX 10.0 # 保存拓扑和坐标 saveamberparm complex system.prmtop system.inpcrd savepdb complex system_solvated.pdb quit

🔍关键检查点：

• 若 PDB 中 DNA 残基名为ADE,THY等旧命名，需重命名为DA,DT，否则 tleap 报 unknown residue；
• combine之后 tleap 会把两条链重排成连续编号（蛋白 1..N，核酸 N+1..M），后面写 restraintmask 时以重排后的编号为准，用desc complex核对；
• 5′/3′ 末端残基类型（DA5/DA3）没给对时 tleap 会在末端缺磷酸或多余氧，check阶段就能暴露。

四、能量最小化与逐步平衡

官方教程对核酸体系用的是"强约束最小化 → 约束下缓慢升温 → 弱约束 NPT 平衡"三步，比纯蛋白流程更保守（核酸初始结构常有实验几何与力场不匹配的张力）。三个输入文件可直接复用：

mini.in —— 25 kcal/mol 强约束下最小化 DNA：

Minimzing the system with 25 kcal/mol restraints on DNA, 500 steps of steepest descent and 500 of conjugated gradient &cntrl imin=1, ntx=1, irest=0, ntpr=50, ntf=1, ntb=1, cut=9.0, ntr=1, maxcyc=1000, ncyc=500, ntmin=1, restraintmask=':1-24', restraint_wt=25.0, &end &ewald ew_type = 0, skinnb = 1.0, &end

heat.in —— 约束下 100 K → 300 K 分段升温（NVT）：

Heating the system with 25 kcal/mol restraints on DNA, V=const &cntrl imin=0, ntx=1, ntpr=500, ntwr=500, ntwx=500, ntwe=500, nscm=5000, ntf=2, ntc=2, ntb=1, ntp=0, nstlim=50000, t=0.0, dt=0.002, cut=9.0, tempi=100.0, ntt=1, ntr=1, nmropt=1, restraintmask=':1-24', restraint_wt=25.0, &end &wt type='TEMP0', istep1=0, istep2=5000, value1=100.0, value2=300.0, &end &wt type='TEMP0', istep1=5001, istep2=50000, value1=300.0, value2=300.0, &end &wt type='END', &end

eq.in —— 0.5 kcal/mol 弱约束下 NPT 平衡 500 ps：

Equilibrating the system with 0.5 kcal/mol restraints on DNA, during 500 ps, at constant T=300 K, P=1 ATM &cntrl imin=0, ntx=5, ntpr=500, ntwr=500, ntwx=500, ntwe=500, nscm=5000, ntf=2, ntc=2, ntb=2, ntp=1, tautp=0.2, taup=0.2, nstlim=25000, t=0.0, dt=0.002, cut=9.0, ntt=1, ntr=1, irest=1, restraintmask=':1-24', restraint_wt=0.5, &end &ewald ew_type = 0, skinnb = 1.0, &end

三个文件的要点：

• 约束从 25 → 0.5 kcal/mol 阶梯释放：最小化和升温阶段把 DNA（:1-24为教程体系 24 个核苷酸残基）位置约束住，只让水和离子松弛，避免溶剂化冲击破坏晶体几何；平衡阶段降到 0.5 让核酸开始弛豫但不上蹿下跳。生产前再把ntr=0放开跑一段无约束平衡；
• 升温用 nmropt=1 + TEMP0 分段：前 5000 步从 100 K 线性升到 300 K，之后恒温——比一步到位加热温和得多；
• ntc=2/ntf=2（SHAKE）：从升温阶段就启用，配合 dt=0.002 ps。

运行顺序（注意每步的-ref都要给）：

pmemd -O -i mini.in -o mini.out -p dna.prmtop -c dna.rst7 -r mini.rst -ref dna.rst7 pmemd -O -i heat.in -o heat.out -p dna.prmtop -c mini.rst -r heat.rst -x heat.nc -e heat.mden -ref mini.rst pmemd -O -i eq.in -o eq.out -p dna.prmtop -c heat.rst -r eq.rst -x eq.nc -e eq.mden -ref heat.rst

⚠️实测坑：pmemd 16+ 运行时一律加-O（overwrite），否则输出文件已存在时直接报Unit NN Error on OPEN终止；设置了ntwe的输入文件必须配-e mden文件名，否则同样报 OPEN 错误。教程网页上的命令行没写这两点，照抄会翻车。

蛋白–核酸复合物在 eq 之后建议再加一步完全无约束的短平衡（ntr=0，50–100 ps），观察蛋白–核酸界面接触是否维持，再进生产。

五、生产模拟

官方教程的 md.in（NPT，200 ps 量级示例）：

production &cntrl imin=0, ntx=5, ntpr=1000, ntwr=1000, ntwx=1000, ntwe=1000, nscm=1000, ntf=2, ntc=2, ntb=2, ntp=1, tautp=5.0, taup=5.0, nstlim=100000, t=0.0, dt=0.002, cut=9.0, ntt=1, irest=1, iwrap=1, ioutfm=1, &end &ewald ew_type = 0, skinnb = 1.0, &end

pmemd -i md.in -o md.out -p dna.prmtop -c eq.rst -r md.rst7 -x md.nc -inf mdinfo

生产实践建议：

• 教程的nstlim=100000 × dt=0.002 = 200 ps只是演示量级；蛋白–核酸复合物的科研模拟建议 100 ns 起步（蛋白–核酸界面水合与离子再分布的弛豫时间比蛋白–蛋白界面长）；
• 大体系上 GPU：pmemd.cuda直接替换pmemd，输入文件不用改；
• tautp/taup=5.0比平衡阶段（0.2）放松，避免过度耦合扭曲动力学；
• ioutfm=1输出 NetCDF 轨迹，长模拟必须（ASCII mdcrd 会大到没法处理）。

六、核酸特征分析（cpptraj 实战）

这是蛋白–核酸流程区别于纯蛋白流程的核心验收环节。官方教程用 cpptraj 完成三类分析：

6.1 轨迹预处理

parm dna.prmtop trajin md.nc autoimage strip :WAT trajout md-nowater.nc netcdf

autoimage处理周期性边界，strip :WAT去水减小轨迹体积。去水后要用 parmed 生成对应的新拓扑（strip :WAT+parmout dna-nowater.prmtop）才能在 VMD 里看。

6.2 螺旋参数：nastruct

cpptraj 的nastruct命令一次输出核酸构象的三类参数文件：

判断标准（B-DNA 常见值域）：Rise ≈ 3.3–3.4 Å，Twist ≈ 34–36°，Slide 为小幅负值。若 Twist 漂到 30° 以下或 Rise 超过 3.6 Å，通常意味着力场/离子环境出了问题。教程里用一段 awk 提取中心碱基步的 Twist 随时间演化：

sed -e '/^$/d' BPstep.data_6-7 > BPstep.data_6-7.new awk '{print $1 "\t" $9}' BPstep.data_6-7.new > twist_6-7.dat

DDD 案例的示范结果：A-DNA 起点在模拟中自发向 B-DNA 转变，Twist 在约第 180 帧附近抬升，与 RMSD 收敛时间点对应——构象转变事件要同时用 RMSD + 螺旋参数交叉确认，只看 RMSD 会漏判。

6.3 糖环 pucker

cpptraj 的pucker命令计算五元糖环的伪旋转相位角（Altona & Sundaralingam 方法）。注意 mask 写法：一条:1-24@C1' ...的范围 mask 会被 cpptraj 静默截成只算第一个残基（每列 mask 只取首个匹配原子），必须逐残基写：

pucker p1 :1@C1' :1@C2' :1@C3' :1@C4' :1@O4' out pucker.dat pucker p2 :2@C1' :2@C2' :2@C3' :2@C4' :2@O4' out pucker.dat ...（逐残基写全）

C2′-endo 对应伪旋转相位角约 110–180°（B-DNA 特征），C3′-endo 约 0–40°（A-DNA / RNA 特征）。

实测示例（Drew-Dickerson 十二聚体，A-DNA 起点，OL15+OPC，按本节流程跑 200 ps 生产）：全链 24 个残基 × 100 帧中仅 31% 采样落在 C2′-endo 区，残基 2–6 和 15–19（两条链的 5′ 端半段）超过 80% 的帧仍停留在 C3′-endo 区——A→B 转变在 200 ps 内只完成了一半。这与官方教程"~800 帧太少、延长模拟分布才平滑"的提示一致：糖环 pucker 的构象转变是慢变量，判断 B-form 是否达成需要 ns 级采样，200 ps 只够看到转变开始。RNA 体系判据反过来：主峰应落在 C3′-endo。

6.4 其他常规分析

七、特殊注意事项

1. 离子浓度设置

• 蛋白–核酸复合物带强负电（磷酸骨架），需足够阳离子屏蔽。先addIons ... 0中和，再按盒子体积追加 100–150 mM NaCl；
• 高盐体系的离子数要自己算：n ≈ 0.15 × 0.6022 × V(ų) 个离子对，别拍脑袋；
• 可用parmed查看残基电荷核对净电荷：

parmed system.prmtop parmed> printDetails :1-100

2. Mg²⁺ 处理

• 若结构中有 Mg²⁺（常见于核酸酶、核酶、聚合酶）：
• Amber 的 JC 离子模型支持 Mg²⁺，但对配位键的描述精度有限——12-6 LJ 球形模型无法重现八面体配位几何；
• 对关键结构/催化 Mg²⁺，建议用Li/Merz 12-6-4 非键模型（在 12-6 基础上加诱导偶极项，配位几何明显改善）或MCPB.py键合模型（成键参数化，最准但最贵）；
• 弱结合的游离 Mg²⁺（只是电荷屏蔽用）用 JC 模型即可。

3. 周期性边界与盒子大小

• 核酸较长且各向异性（长条形），旋转扩散快，需确保任意方向水层 ≥10 Å，避免双链首尾穿过周期边界与自身相互作用；
• solvateOct（八面体盒子）对棒状溶质特别省水——教程 DDD 体系 9 Å 水层只用约 5336 个水，立方盒要多出一截。

4. 验证核酸构象（验收清单）

平衡与生产后逐项检查 DNA/RNA 是否保持目标构象（B/A/Z 或 A-form）：

• nastruct螺旋参数落在目标构象值域（见 6.2）；
• 糖环 pucker 主峰位置正确（见 6.3）；
• RMSD 收敛后检查是否发生 A↔B 构象漂移（DDD 是有意让它转，你的体系如果不是研究对象就该警惕）。

八、常见错误排查

九、结语

蛋白–核酸复合物模拟的完整心法一句话：力场组合（蛋白 ff19SB + 核酸 OL15/OL21/OL3 + OPC 水 + 匹配离子）→ tleap 两段式 combine → 强约束起步的阶梯式平衡 → 用 nastruct 和 pucker 验收构象。前三点决定模拟能不能稳定跑下去，最后一点决定结果能不能信。

如提供具体体系（如 p53–DNA 复合物、Cas9–sgRNA–DNA 三元复合物），可进一步定制流程（如处理 sgRNA 修饰、多链组装、Mg²⁺ 配位等）。

关键字

蛋白-核酸复合物、分子动力学模拟、AMBER、tleap、OL15力场、OPC水模型、nastruct、糖环pucker

参考文献与延伸阅读

1. AMBER Hub 官方核酸模拟教程（Drew-Dickerson dodecamer 全流程，本文 tleap/mdin/nastruct 命令来源）: https://amberhub.chpc.utah.edu/analisis-of-nucleic-acid-simulation/
2. OL15 力场系列综述之一（J. Chem. Theory Comput., 2016）: DOI 10.1021/acs.jctc.6b00186
3. 核酸力场变体系统评估（Galindo-Murillo et al., Nucleic Acids Res. 2017, 45(7): 4217）: DOI 10.1093/nar/gkx136
4. Love et al. dsDNA 力场 vdW 参数扫描，提出 OL21+OPC 推荐（J. Chem. Theory Comput., 2023）: DOI 10.1021/acs.jctc.3c01164
5. OPC 水模型（Izadi, Anandakrishnan, Onufriev, J. Phys. Chem. Lett. 2014, 5, 3863）: DOI 10.1021/jz501780a
6. OL15/OL21 官方主页（Olomouc 力场组）: http://ffol.upol.cz