如何将化学分子SMILES字符串转化为Pytorch图数据结构——ESOL分子水溶性数据集解析

硬核原创,转载请注明出处:
https://leytton.blog.csdn.net/article/details/130406553

一、前言

在用Pytorch图神经网络对化学分子进行数据分析的时候,经常使用现有的数据集。看到自动处理完毕的数据结构,里面的特征值让我们一脸懵逼,不知道代表的是什么含义。本文将带大家分析这些数据结构的来龙去脉。

二、数据原始特征

在使用图神经网络(GNN)对化学分子进行水溶性预测的实验中,加载了MoleculeNetESOL数据集。我们打开原始的csv文件,结构是这样的(非专业翻译,有误恳请留言纠正):

表头含义示例
Compound ID化合物ID2-pyrrolidone
ESOL predicted log solubility in mols per litreESOL预测对数溶解度(mol/L)0.243
Minimum Degree最小度1
Molecular Weight分子量85.10600000000001
Number of H-Bond Donors氢键供体数1
Number of Rings环数1
Number of Rotatable Bonds可旋转键数0
Polar Surface Area极性表面积29.1
measured log solubility in mols per litre测量对数溶解度(mol/L)1.07
smiles分子SMILES字符串O=C1CCCN1

从smiles中可以看到,这个分子有OCN、H(一般省略)四种元素,除去H有六个原子。

三、分析预处理数据

加载数据

如下面代码,加载ESOL数据集后将其打印出来:

from torch_geometric.datasets import MoleculeNet

data = MoleculeNet(root="data", name="ESOL")

print("Dataset Size:", len(data))
print("Dataset classes:", data.num_classes)
print("Dataset features:", data.num_features)
Dataset Size: 1128
Dataset classes: 734
Dataset features: 9

从结果可以看到,有1128个分子样本,734种类型,每个分子有9个特征。

分析数据

我们选择第11个分子(smiles比较短)进行分析:

print("Sample:", data[10])
print("Sample y:", data[10].y)
Sample: Data(x=[6, 9], edge_index=[2, 12], edge_attr=[12, 3], smiles='O=C1CCCN1', y=[1, 1])
Sample y: tensor([[1.0700]])

可以看到x、edge_index、edge_attr 是二维数组,y可以看成一个值(水溶性)。

关于水溶性参考《理化性质|(logSw和logP)小分子化合物水溶性和脂溶性指标

画出分子图

根据SMILES字符串,将其分子图画出来:

from rdkit import Chem
from rdkit.Chem import Draw

molecule = Chem.MolFromSmiles(data[10]["smiles"])
# Draw.MolToFile(molecule, "mol.png")
Draw.MolToImage(molecule)

在这里插入图片描述

edge_index数据分析

edge_index数组打印出来:

print(data[10].edge_index.T)
tensor([[0, 1],
        [1, 0],
        [1, 2],
        [1, 5],
        [2, 1],
        [2, 3],
        [3, 2],
        [3, 4],
        [4, 3],
        [4, 5],
        [5, 1],
        [5, 4]])

这是O、C、N六个原子的连接关系。

x数据分析

x数组打印出来:

print(data[10].x.shape)
print(data[10].x)
torch.Size([6, 9])
tensor([[8, 0, 1, 5, 0, 0, 3, 0, 0],
        [6, 0, 3, 5, 0, 0, 3, 0, 1],
        [6, 0, 4, 5, 2, 0, 4, 0, 1],
        [6, 0, 4, 5, 2, 0, 4, 0, 1],
        [6, 0, 4, 5, 2, 0, 4, 0, 1],
        [7, 0, 3, 5, 1, 0, 3, 0, 1]])

这就不太看得懂了,看起来像是描述6个原子9个特征的二维数组。

四、真相:从SMILES字符串得到

作者查阅资料无果,那么久只能去分析MoleculeNet中的代码了,到底对原始数据进行了怎样的处理,x中的数据是怎样来的。

点进去看到一个process函数应该是处理数据的,我对其进行了注释:

# Format: name: [display_name, url_name, csv_name, smiles_idx, y_idx]
    names = {
        'esol': ['ESOL', 'delaney-processed.csv', 'delaney-processed', -1, -2],
        'freesolv': ['FreeSolv', 'SAMPL.csv', 'SAMPL', 1, 2],
        'lipo': ['Lipophilicity', 'Lipophilicity.csv', 'Lipophilicity', 2, 1],
        'pcba': ['PCBA', 'pcba.csv.gz', 'pcba', -1,
                 slice(0, 128)],
        'muv': ['MUV', 'muv.csv.gz', 'muv', -1,
                slice(0, 17)],
        'hiv': ['HIV', 'HIV.csv', 'HIV', 0, -1],
        'bace': ['BACE', 'bace.csv', 'bace', 0, 2],
        'bbbp': ['BBPB', 'BBBP.csv', 'BBBP', -1, -2],
        'tox21': ['Tox21', 'tox21.csv.gz', 'tox21', -1,
                  slice(0, 12)],
        'toxcast':
        ['ToxCast', 'toxcast_data.csv.gz', 'toxcast_data', 0,
         slice(1, 618)],
        'sider': ['SIDER', 'sider.csv.gz', 'sider', 0,
                  slice(1, 28)],
        'clintox': ['ClinTox', 'clintox.csv.gz', 'clintox', 0,
                    slice(1, 3)],
    }
    
def process(self):
        with open(self.raw_paths[0], 'r') as f: #读取原始数据文件
            dataset = f.read().split('\n')[1:-1] #按行分割,并去掉第一行
            dataset = [x for x in dataset if len(x) > 0]  # 去掉空行

        data_list = []
        for line in dataset:  #遍历每行
            line = re.sub(r'\".*\"', '', line)  # 去掉".*"字符串
            line = line.split(',') #逗号分隔

            smiles = line[self.names[self.name][3]] #获取到smiles字符串
            ys = line[self.names[self.name][4]] #获取到y值
            ys = ys if isinstance(ys, list) else [ys] #将y值统一成数组形式

            ys = [float(y) if len(y) > 0 else float('NaN') for y in ys] #将y转成float类型
            y = torch.tensor(ys, dtype=torch.float).view(1, -1) #将y转成torch.float类型
			
			# 重点:获取x、edge_index、edge_attr数据,需要查看from_smiles函数
            data = from_smiles(smiles)  
            data.y = y  #y处理完毕

            if self.pre_filter is not None and not self.pre_filter(data):
                continue

            if self.pre_transform is not None:
                data = self.pre_transform(data)

            data_list.append(data)

        torch.save(self.collate(data_list), self.processed_paths[0])

从上面分析可以知道,原来x、edge_index、edge_attr数据都是通过将smile字符串传递到from_smiles函数获取到的!

from_smiles函数如下:

def from_smiles(smiles: str, with_hydrogen: bool = False,
                kekulize: bool = False) -> 'torch_geometric.data.Data':
    # 太多了省略。。。
    return Data(x=x, edge_index=edge_index, edge_attr=edge_attr, smiles=smiles)

这下可以参考这个函数的代码进一步分析了。

我们直接指定smiles进行分析:

smiles='O=C1CCCN1'
from rdkit import Chem

mol = Chem.MolFromSmiles(smiles)
for atom in mol.GetAtoms():
    print(f'原子序号:{atom.GetAtomicNum()}, 手性信息:{atom.GetChiralTag()}, 度:{atom.GetTotalDegree()}, 电荷:{atom.GetFormalCharge()}, 连接氢原子数:{atom.GetTotalNumHs()}, 自由基:{atom.GetNumRadicalElectrons()}, 杂化类型:{atom.GetHybridization()}, 芳香性:{atom.GetIsAromatic()}, 是否在环上:{atom.IsInRing()}')
原子序号:8, 手性信息:CHI_UNSPECIFIED, 度:1, 电荷:0, 连接氢原子数:0, 自由基:0, 杂化类型:SP2, 芳香性:False, 是否在环上:False
原子序号:6, 手性信息:CHI_UNSPECIFIED, 度:3, 电荷:0, 连接氢原子数:0, 自由基:0, 杂化类型:SP2, 芳香性:False, 是否在环上:True
原子序号:6, 手性信息:CHI_UNSPECIFIED, 度:4, 电荷:0, 连接氢原子数:2, 自由基:0, 杂化类型:SP3, 芳香性:False, 是否在环上:True
原子序号:6, 手性信息:CHI_UNSPECIFIED, 度:4, 电荷:0, 连接氢原子数:2, 自由基:0, 杂化类型:SP3, 芳香性:False, 是否在环上:True
原子序号:6, 手性信息:CHI_UNSPECIFIED, 度:4, 电荷:0, 连接氢原子数:2, 自由基:0, 杂化类型:SP3, 芳香性:False, 是否在环上:True
原子序号:7, 手性信息:CHI_UNSPECIFIED, 度:3, 电荷:0, 连接氢原子数:1, 自由基:0, 杂化类型:SP2, 芳香性:False, 是否在环上:True

如上所示,这9个特征就是x变量中每个原子的含义,对其进行一些编码变换就构造成了x变量。具体的原子更多的属性,可以参考 RDKit 文档

接下来我们分析edge_attredge_index变量含义:

for bond in mol.GetBonds(): #便利所有的键
    i = bond.GetBeginAtomIdx()
    j = bond.GetEndAtomIdx()
    print(f'连接:{i,j},{j,i}')
    print(f'键的类型:{bond.GetBondType()}, Stereo:{bond.GetStereo()}, 是否共轭:{bond.GetIsConjugated()}')
连接:(0, 1),(1, 0)
键的类型:DOUBLE, Stereo:STEREONONE, 是否共轭:True
连接:(1, 2),(2, 1)
键的类型:SINGLE, Stereo:STEREONONE, 是否共轭:False
连接:(2, 3),(3, 2)
键的类型:SINGLE, Stereo:STEREONONE, 是否共轭:False
连接:(3, 4),(4, 3)
键的类型:SINGLE, Stereo:STEREONONE, 是否共轭:False
连接:(4, 5),(5, 4)
键的类型:SINGLE, Stereo:STEREONONE, 是否共轭:False
连接:(5, 1),(1, 5)
键的类型:SINGLE, Stereo:STEREONONE, 是否共轭:True

这就是分子SMILES字符串转化成图数据结构的过程,可以看到只用到了原始数据里的SMILES字符串水溶性结果

在Pytorch官网找了半天没找到数据集的说明资料,等我分析完后,才发现,这里已经有大佬发表了相关文章。不过,如果不知道图结构数据是从SMILES字符串分析得到,很难通过关键字找到这些资料。
How to turn a SMILES string into a molecular graph for Pytorch Geometric

  • 8
    点赞
  • 33
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值