mol2graph and graph2mol #rdkit #igraph

I posted about mol to graph object before.
In the blog post, I wrote example that convert RDKit mol object to igraph object. It was one way. There was no method igraph to rdkit mol object.
So I wrote very simple converter from graph to molecule.

First, import modules.

import numpy as np
import pandas as pd
import igraph
from rdkit import Chem
from rdkit.Chem.rdchem import RWMol
from rdkit.Chem import Draw
from rdkit.Chem import rdmolops
from rdkit.Chem.Draw import IPythonConsole
IPythonConsole.ipython_useSVG = True

Then define two way function, mol2graph and graph2mol. It is very simple.I did not sanitize process because I could not handle some compounds. RWMol method is very useful to do this work.

def mol2graph(mol):
    atoms_info = [ (atom.GetIdx(), atom.GetAtomicNum(), atom.GetSymbol()) for atom in mol.GetAtoms()]
    bonds_info = [(bond.GetBeginAtomIdx(), bond.GetEndAtomIdx(), bond.GetBondType(), bond.GetBondTypeAsDouble()) for bond in mol.GetBonds()]
    graph = igraph.Graph()
    for atom_info in atoms_info:
        graph.add_vertex(atom_info[0], AtomicNum=atom_info[1], AtomicSymbole=atom_info[2])
    for bond_info in bonds_info:
        graph.add_edge(bond_info[0], bond_info[1], BondType=bond_info[2], BondTypeAsDouble=bond_info[3])
    return graph

def graph2mol(graph): 
    emol = RWMol()
    for v in graph.vs():
        emol.AddAtom(Chem.Atom(v["AtomicNum"]))
    for e in graph.es():
        emol.AddBond(e.source, e.target, e['BondType'])
    mol = emol.GetMol()
    return mol

Finally, I checked my function on jupyter notebook. And It worked well.

All code is uploaded my repo and can check from following URL. https://nbviewer.jupyter.org/github/iwatobipen/chemo_info/blob/master/rdkit_notebook/mol2g_g2mol.ipynb

Advertisements

Convert rdkit molecule object to igraph graph object.

Molecules are often handled as graph in chemoinformatics.
There are some libraries for graph analysis in python. Today, I wrote a sample script that convert from molecule to graph.
I used python-igraph and rdkit.
RDkit has method to get adjacency matrix from molecule so, I used the method.
Code is following.

from rdkit import Chem
from rdkit.Chem import rdmolops
import igraph
import numpy as np
# mol2graph function can convert molecule to graph. And add some node and edge attribute.
def mol2graph( mol ):
    admatrix = rdmolops.GetAdjacencyMatrix( mol )
    bondidxs = [ ( b.GetBeginAtomIdx(),b.GetEndAtomIdx() ) for b in mol.GetBonds() ]
    adlist = np.ndarray.tolist( admatrix )
    graph = igraph.Graph()
    g = graph.Adjacency( adlist ).as_undirected()
    for idx in g.vs.indices:
        g.vs[ idx ][ "AtomicNum" ] = mol.GetAtomWithIdx( idx ).GetAtomicNum()
        g.vs[ idx ][ "AtomicSymbole" ] = mol.GetAtomWithIdx( idx ).GetSymbol()
    for bd in bondidxs:
        btype = mol.GetBondBetweenAtoms( bd[0], bd[1] ).GetBondTypeAsDouble()
        g.es[ g.get_eid(bd[0], bd[1]) ][ "BondType" ] = btype
        print( bd, mol.GetBondBetweenAtoms( bd[0], bd[1] ).GetBondTypeAsDouble() )
    return g

Now test it.

# Oseltamivir ( tamiful )
mol2 = Chem.MolFromSmiles( "CCOC(=O)C1=C[C@@H](OC(CC)CC)[C@H](NC(C)=O)[C@@H](N)C1" )
g2=mol2graph(mol2)
(0, 1) 1.0
(1, 2) 1.0
(2, 3) 1.0
(3, 4) 2.0
(3, 5) 1.0
(5, 6) 2.0
(6, 7) 1.0
....

Seems work fine. ;-)
Next, get bond infromation.

for e in g2.es:
    print(e)
igraph.Edge(<igraph.Graph object at 0x10ae999a8>, 0, {'BondType': 1.0})
igraph.Edge(<igraph.Graph object at 0x10ae999a8>, 1, {'BondType': 1.0})
igraph.Edge(<igraph.Graph object at 0x10ae999a8>, 2, {'BondType': 1.0})
......

OK next, I tried to plot tamiful as Graph.
Igraph has plot function.
I added some options for plotting.

from igraph import GraphBase
layout = g2.layout_graphopt()
color_dict = {"C": "gray", "N": "blue", "O" :  "white" }
my_plot = igraph.Plot()
my_plot.add(g2,layout=layout,bbox=(400,400),
             margin=90,
             vertex_color = [color_dict[atom] for atom in g2.vs["AtomicSymbole"]],
             vertex_size = [ v.degree()*10 for v in g2.vs ])
my_plot.show()

tamiful

Ummmm @_@ not so good view.
But fun!