📚 课程目录第4章 RDKit分子操作

4、RDKit分子操作：读取、写入、显示分子

做分子信息学，绕不开 RDKit。我个人习惯把它当作分子的「瑞士军刀」——读取、改造、输出，一条龙。这一章我们从最基础的操作讲起：怎么把 SMILES 变成 RDKit 能懂的对象，怎么把它写回成各种格式，以及怎么在屏幕上看到分子长什么样。

你想想看，SMILES 说到底就是一串字符串，但分子是有三维结构的。RDKit 这个中间层就起到了「翻译」作用。我刚开始用的时候，经常忘记检查分子是否读成功，结果后续代码一顿报错。嗯，这里要注意——读取一定要判空。

4.1 分子读取：从 SMILES 到 RDKit 对象

RDKit 读取 SMILES 的核心函数是 Chem.MolFromSmiles。它接受一个字符串，返回一个 RDKit Mol 对象。如果 SMILES 不合法，返回 None。我在项目中遇到过一个坑：用户提交的 SMILES 里混了空格和换行，直接传进去就炸了。所以读之前先 strip 一下。

from rdkit import Chem

smiles = "CCO"
mol = Chem.MolFromSmiles(smiles)
if mol is not None:
    print("分子读取成功")
else:
    print("无效的 SMILES")

你也可以一次读多个，比如从文件里加载。很多数据集里的 SMILES 带标题行，用 pandas 读出来后逐行转就行。

记得对每个 MolFromSmiles 的结果做 None 检查。我曾经因为没检查，后续调用 mol.GetNumAtoms() 直接报 AttributeError，排查了半小时。

4.2 分子写入：输出为不同格式

读进来之后，往往要写回去。RDKit 提供了丰富的写入函数，最常用的三种：

函数	输出格式	典型用途
`Chem.MolToSmiles`	标准 SMILES	规范化、去重、存储
`Chem.MolToMolBlock`	MOL (SDfile 格式)	保留坐标信息，对接其他工具
`Chem.MolToSDF`	SDF (多个分子)	批量写入，含属性

举个例子，我想把乙醇的 SMILES 写回标准化的形式：

smiles_in = "CCO"
mol = Chem.MolFromSmiles(smiles_in)
canonical_smiles = Chem.MolToSmiles(mol)
print(canonical_smiles)  # 输出: CCO

如果分子有坐标，写成 MOL 文件就可以用其他软件打开查看：

mol_block = Chem.MolToMolBlock(mol)
with open("ethanol.mol", "w") as f:
    f.write(mol_block)

注意事项： MolToSmiles 默认会做规范化（如调整芳香性表示）。如果你需要保留原始的 SMILES 写法，用 MolToSmiles(mol, canonical=False)。

4.3 分子显示与可视化

干了这么久，我还是觉得「亲眼看到分子」最踏实。RDKit 内置了基于 rdkit.Chem.Draw 的可视化工具。在 Jupyter 里直接用 Draw.MolToImage 就能出图。

from rdkit.Chem import Draw

mol = Chem.MolFromSmiles("CCO")
img = Draw.MolToImage(mol)
img.show()  # 在 Jupyter 里会自动显示

如果你没有图形界面，可以用 Draw.MolToFile 保存成 PNG：

Draw.MolToFile(mol, "ethanol.png")

当然，显示的样式可以调整。比如我习惯把原子编号显示出来，方便 debug：

from rdkit.Chem.Draw import rdMolDraw2D
d = rdMolDraw2D.MolDraw2DSVG(300, 300)
opts = d.drawOptions()
opts.addAtomIndices = True
d.DrawMolecule(mol)
d.FinishDrawing()
svg = d.GetDrawingText()

这个 SVG 字符串可以直接嵌入网页或拿去打印。我在做报告时经常这么干，观众能直接看到原子编号。

在非 Jupyter 环境下，img.show() 可能不会弹出窗口。建议用 img.save("mol.png") 离线查看。

4.4 操作链路：读取 → 修改 → 显示

实际工作中，读取、修改、输出很少分开做。我举个例子：从数据文件里读一批 SMILES，给每个分子加上氢原子，然后显示并统计原子数。

smiles_list = ["CCO", "c1ccccc1", "CC(=O)O"]
for s in smiles_list:
    mol = Chem.MolFromSmiles(s)
    if mol is None:
        continue
    # 加氢
    mol_with_H = Chem.AddHs(mol)
    # 显示
    Draw.MolToFile(mol_with_H, f"{s.replace('/', '_')}.png")
    print(f"{s} -> {mol_with_H.GetNumAtoms()} 个原子（含H）")

这个链路非常典型。我个人习惯把读取、验证、转换包装成一个函数，后面调起来方便。

下面这张图，概括了整个 RDKit 分子操作的核心链条：

这个流程里，我特别想强调一点：不要以为读出来的 Mol 对象就是「干净的」。它可能没有氢原子、没有立体化学信息，甚至坐标缺失。在实际项目中，我一般会在读取后立即调用 Chem.SanitizeMol 做一次完整性校验。代码就像这样：

mol = Chem.MolFromSmiles(smiles)
if mol is not None:
    Chem.SanitizeMol(mol)
    print("分子已规范化")

养成“读取即校验”的习惯，能省掉后面 80% 的调试时间。我自己的 pipeline 里，读取函数会返回一个 (mol, error_message) 元组，把异常抛在明面上。

好了，这一章就到这里。学会读取、写入、显示，你就能在 RDKit 的世界里自由穿梭了。后面的章节会在此基础上，教你如何操作原子、键、子结构，一步步搭出真正的分子处理管线。