#!/usr/bin/env python3 # -*- coding: utf-8 -*- """MD 结果页生成：动画 + RMSD/RMSF + 方法与边界用法：python build_md_report.py <蛋白名> <作用> <域说明> """ import base64 import io import json import os import sys sys.path.insert(0, '/Users/liuxinyuan/Desktop/Hermes输出-工作类/脚本') import brand_kit as BK # noqa: E402 W = '/Users/liuxinyuan/Desktop/Hermes输出-工作类' d, name, uniprot, pdbid, role, domdesc = sys.argv[1:7] s = json.load(open(os.path.join(d, 'summary.json'), encoding='utf-8')) def b64(p, mime): return 'data:%s;base64,%s' % (mime, base64.b64encode(open(p, 'rb').read()).decode()) gif = b64(os.path.join(d, 'md.gif'), 'image/gif') vid = b64(os.path.join(d, 'md.mp4'), 'video/mp4') curve = b64(os.path.join(d, 'md_curves.png'), 'image/png') mp4 = os.path.join(d, 'md.mp4') import csv as _csv _rr = list(_csv.DictReader(open(os.path.join(d, 'rmsd.csv')))) _vals = [float(x['rmsd_A']) for x in _rr] import statistics as _st _a = _vals[int(len(_vals) * 0.8):int(len(_vals) * 0.9)] _b = _vals[int(len(_vals) * 0.9):] _plateau = abs(_st.mean(_b) - _st.mean(_a)) < 0.3 _gap = abs(_st.mean(_b) - _st.mean(_a)) _drift = _vals[-1] - _vals[int(len(_vals) * 0.2)] _drg = s['rg_end_A'] - s['rg_start_A'] _interp = [] if _drg < -0.4: _interp.append('回转半径由 %.1f Å 降到 %.1f Å（%.1f%%），结构在模拟中明显压实。这是隐式溶剂（GB）模型的已知偏差——缺少显式水与离子，疏水塌缩被放大；不能据此判断酶在真实缓冲体系中会收缩。' % (s['rg_start_A'], s['rg_end_A'], abs(_drg) / s['rg_start_A'] * 100)) elif _drg > 0.4: _interp.append('回转半径由 %.1f Å 增至 %.1f Å，结构略舒展。' % (s['rg_start_A'], s['rg_end_A'])) else: _interp.append('回转半径基本稳定（%.1f → %.1f Å），整体紧凑度未见明显变化。' % (s['rg_start_A'], s['rg_end_A'])) if _plateau: _interp.append('最后 1 ns 的 RMSD 均值 %.2f Å，与前 0.5 ns 仅差 %.2f Å，已进入平台段；但 5 ns 仍不足以证明充分平衡，只能说明这一段没有持续单向漂移。' % (_st.mean(_b), _gap)) else: _interp.append('最后 1 ns 的 RMSD 均值较前一区间仍变化 %.2f Å（末帧 %.2f Å），未达平台；此时任何“稳定性”结论都不成立，只能用于观察局部柔性。' % (_gap, _vals[-1])) _interp.append('柔性最大的区域集中在 N/C 两端（%s）；末端涨落在多数酶里属正常，不应解读为活性位点不稳。' % '、'.join('%s%d' % (x['name'], x['resSeq']) for x in s['top_flex_residues'][:4])) _interp.append('本体系未包含 DNA 底物、Mg²⁺ 与结晶水：UNG 的催化构象需要底物结合才被组织起来，无底物模拟只能看 apo 态的呼吸运动。') _interp_html = ''.join('

%s

' % x for x in _interp) flex = ''.join('%s%d%.2f' % (x['name'], x['resSeq'], x['rmsf']) for x in s['top_flex_residues']) page = """ %(name)s 分子动力学模拟

%(name)s 分子动力学模拟（隐式溶剂）

UniProt %(uniprot)s ｜起始结构 PDB %(pdbid)s ｜ OpenMM 8.6.1 · AMBER14 + GBSA-OBC2 · 300 K ｜截至 2026-10-17

%(ns).2f ns模拟时长

%(atoms)s主链原子数（含氢）

%(rmsd).2f Å主链 RMSD 均值

%(rmsf).2f ÅCA RMSF 均值

轨迹动画每帧 %(frame_ps).1f ps，共 %(frames)d 帧

若视频不自动播放，点一下播放键；源文件 数据/md/%(dirname)s/md.mp4

RMSD / RMSF 曲线

左：轨迹回放（结构按残基编号着色）｜右：上=主链 RMSD，下=各残基 CA 涨落

这条蛋白为什么值得看动力学 %(role)s

%(domdesc)s

柔性最大的残基（RMSF 前 8）：

%(flex)s

残基	RMSF (Å)

结果解读自动按本次轨迹生成，不预设结论

%(interp)s

方法与边界必须一起讲

这是隐式溶剂 MD，不是全原子显式水模拟。

溶剂用 GBSA-OBC2 连续介质近似，省去水分子与离子，速度快 10–50 倍，但不能用来判断结合水、离子桥、界面水网络；
时长 %(ns).2f ns，只能看局部构象涨落与早期形变，不能据此判断热稳定性排序或货架期；
起始结构为晶体结构（含缺失段由 pdbfixer 补齐），预测段与补齐段的运动学解释力弱；
单条轨迹、无重复，结论属“观察”，不属“统计显著”。

销售怎么用只说观察，不说保证

讲构象不是静止的：酶在 300 K 下持续涨落，活性位点附近的柔性决定它“容不容易被温度/离子条件带偏”；
指到具体残基：RMSF 高的环区往往是突变设计、耐热改造的着力点；
边界话术：客户问“能不能证明更耐热”，回答用实测（半衰期、Tm、加速老化），不用 MD 动画代替。

%(footer)s

""" % { 'name': name, 'uniprot': uniprot, 'pdbid': pdbid, 'role': role, 'domdesc': domdesc, 'ns': s['duration_ns'], 'atoms': s['atoms'], 'rmsd': s['rmsd_mean_A'], 'rmsf': s['rmsf_mean_A'], 'frames': s['frames'], 'frame_ps': s['duration_ns'] * 1000 / max(1, s['frames'] - 1), 'gif': gif, 'vid': vid, 'dirname': os.path.basename(d), 'curve': curve, 'flex': flex, 'footer': '__FOOTER__', 'interp': _interp_html, } page = BK.inject(page, theme='dark', tag='%s 分子动力学' % name, footer_extra='内部资料 · 请勿外发｜隐式溶剂 MD，结论仅限构象涨落') page = page.replace('__FOOTER__', '') out = os.path.join(W, '销售分析', '分子动力学_%s_%s.html' % (name, pdbid)) io.open(out, 'w', encoding='utf-8').write(page) print('✅', out, '%.1f MB' % (len(page.encode()) / 1048576))