GROMACS FRET 距离与取向因子演示

下载本组输入、真实参考结果和绘图脚本，解压后进入资料包目录。inputs 是可以修改的输入，reference 是鸢算科技集群已经计算得到的参考结果。每次运行在 work 下建立新目录，路径写入 LAST-RUN.txt。计算请在已经分配的资源中执行；整组可用 sbatch run.sbatch 提交。

本例能回答什么#

FRET 不只取决于距离，也与供体、受体的跃迁偶极取向有关。这里使用真实泛素轨迹中的两对 Cα 作几何代理，说明 dyecoupl 如何选端点、计算距离、κ² 与模型效率。没有化学染料、光谱标定或真实跃迁偶极，因此结果不能作为该蛋白的实验 FRET 效率预测。

bash run.sh protein
cd "$(cat LAST-RUN.txt)/protein-analysis"
cat fret-proxy.ndx

供体代理是 23、24 残基的 Cα，受体代理是 33、34 残基的 Cα；每组恰好两个原子，顺序确定向量方向。各组位置在本例按相同质量 Cα 的中点计算。

1. 建立索引#

手动执行时先复制 inputs/protein 到独立目录并加载普通 GROMACS 模块。可以用原子选择生成两组：

gmx select -s protein.tpr -select 'resnr 23 24 and name CA' 'resnr 33 34 and name CA' -on fret-proxy.ndx
cat fret-proxy.ndx

核对每组两个原子和编号；资料包 fret.py 也会按真实 GRO 中的残基名称和顺序生成该文件。对化学染料，应根据染料模型和光谱知识定义端点，不能把任意两个原子连线视为真实跃迁偶极。

2. 运行原生工具#

printf '0\n1\n' | gmx dyecoupl -f protein.xtc -n fret-proxy.ndx -ot proxy-r-kappa.xvg -oe proxy-efficiency.xvg -rhist proxy-r-hist.xvg -khist proxy-kappa-hist.xvg -R0 2 -norm

本例人为设定 Förster 半径 R0=2 nm，只是公式演示。输出包含距离/取向时间序列、模型效率，以及距离与取向直方图；-norm 归一化直方图，不会把几何代理转变成真实染料。

3. 看真实轨迹里的结果#

这组代理向量的取向因子较小，不能只用距离套入各向同性平均公式后解释图中的效率。

4. 独立复算的公式#

令单位向量 d、a 分别沿两组偶极方向，rhat 沿两个中心连线，则：

κ² = [d·a − 3(d·rhat)(a·rhat)]²
E  = 1 / [1 + (R/R0)^6 × (2/3)/κ²]

这里 R0 对应各向同性取向因子 2/3 的约定，实际 κ² 在上式中另行代入。零取向因子按模型极限得到零效率；计算脚本应避免除零。不要重复乘入取向修正，或无条件把 κ² 固定成 2/3。

资料包逐帧从坐标独立复算，与工具的最大差别：距离约 3.9e-7 nm、κ² 约 3.0e-7、效率约 1.6e-6，与原生文本输出位数相符。结果表在 independent-fret.csv。

head -30 proxy-r-kappa.xvg
head -30 proxy-efficiency.xvg
head independent-fret.csv

真实 FRET 研究还需要染料连接臂、染料构象采样、实验对应的 R0 和光谱信息；100 ps 蛋白几何示例不替代这些步骤。参数见 GROMACS dyecoupl 官方说明。