Hi all,
We’re using the C# wrappers to RDKit for an enumeration engine and we have
run into an issue with stereochemistry in RDKit.
The basic problem is that when we use a reaction which adds a stereocenter
to a molecule, chirality at carbons that were present in the reactants is
kept, but chirality at carbons introduced by the reaction itself is lost.
In the course of trying to figure this out, I created a Python script that
goes through the steps performed by the enumeration engine for a sample
illustrating the issue so that I could try to debug the problem more
easily, only to find that in Python the chirality is retained throughout
the molecule. I then translated that Python script back into a standalone
C# program … and the problem reappears.
The output molecule should have two chiral centers (the carbon attached to
the fluorine and the carbon attached to the chlorine). The Python version
shows both of these as chiral; the C# version has only the fluorine carbon
as chiral.
I’ve attached the Python script and the C# source file. The C# version has
the Python commands included as comments for easy comparison.
Does anyone know what’s going on here? The only significant difference
between the two that I can see is that Python has rxn.Initialize() while C#
has rxn.initReactantMatchers(). Could this be the cause, and if so what C#
method should I be using instead?
Thanks in Advance,
Bob
--
Bob Funchess, Ph.D. Kelaroo,
Inc
Senior Scientist
www.kelaroo.com
bfunch...@kelaroo.com (858)
259-7561 x3
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using System.IO;
using GraphMolWrap;
namespace StereoTest
{
class StereoTest
{
static void Main(string[] args)
{
string molblock_r = @"
Mrv0541 12291411382D
6 5 0 0 1 0 999 V2000
-2.8875 0.7145 0.0000 Br 0 0 0 0 0 0 0 0 0 6 0 0
-2.0625 0.7145 0.0000 C 0 0 0 0 0 0 0 0 0 5 0 0
-1.6500 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 4 0 0
-2.0625 -0.7145 0.0000 C 0 0 1 0 0 0 0 0 0 1 0 0
-1.6500 -1.4289 0.0000 C 0 0 0 0 0 0 0 0 0 3 0 0
-2.8875 -0.7145 0.0000 F 0 0 0 0 0 0 0 0 0 2 0 0
1 2 1 0 0 0 0
2 3 1 0 0 0 0
3 4 1 0 0 0 0
4 5 1 0 0 0 0
4 6 1 1 0 0 0
M END
";
// molblock_r = """ <same heredoc? """
string molblock_p = @"
Mrv0541 12291411382D
9 8 0 0 1 0 999 V2000
3.1821 0.6188 0.0000 C 0 0 0 0 0 0 0 0 0 3 0 0
3.8966 0.2063 0.0000 C 0 0 2 0 0 0 0 0 0 1 0 0
3.8966 -0.6187 0.0000 F 0 0 0 0 0 0 0 0 0 2 0 0
4.6111 0.6188 0.0000 C 0 0 0 0 0 0 0 0 0 4 0 0
5.3256 0.2062 0.0000 C 0 0 0 0 0 0 0 0 0 5 0 0
6.0400 0.6187 0.0000 C 0 0 0 0 0 0 0 0 0 7 0 0
6.7545 0.2062 0.0000 C 0 0 1 0 0 0 0 0 0 8 0 0
7.4690 0.6187 0.0000 C 0 0 0 0 0 0 0 0 0 9 0 0
6.7545 -0.6188 0.0000 Cl 0 0 0 0 0 0 0 0 0 10 0 0
1 2 1 0 0 0 0
2 3 1 1 0 0 0
2 4 1 0 0 0 0
4 5 1 0 0 0 0
5 6 1 0 0 0 0
7 6 1 0 0 0 0
7 8 1 0 0 0 0
7 9 1 1 0 0 0
M END
";
// molblock_p = """ <same heredoc> """
ChemicalReaction rxn = new ChemicalReaction();
// rxn = Chem.rdChemReactions.ChemicalReaction()
ROMol template_r = RWMol.MolFromMolBlock(molblock_r, true, false);
// template_r = Chem.MolFromMolBlock(molblock_r,True,False)
ROMol reactant_H = template_r.addHs(true);
// reactant_H = Chem.AddHs(template_r,True)
rxn.addReactantTemplate(reactant_H);
// rxn.AddReactantTemplate(reactant_H)
ROMol template_p = RWMol.MolFromMolBlock(molblock_p, true, false);
// template_p = Chem.MolFromMolBlock(molblock_p,True,False)
ROMol product_H = template_p.addHs(true);
// product_H = Chem.AddHs(template_p,True)
rxn.addProductTemplate(product_H);
// product_H = Chem.AddHs(template_p,True)
rxn.initReactantMatchers();
// rxn.Initialize()
string molblock_r0 = @"
Mrv0541 12291411372D
6 5 0 0 1 0 999 V2000
-2.8875 0.7145 0.0000 Br 0 0 0 0 0 0 0 0 0 6 0 0
-2.0625 0.7145 0.0000 C 0 0 0 0 0 0 0 0 0 5 0 0
-1.6500 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 4 0 0
-2.0625 -0.7145 0.0000 C 0 0 1 0 0 0 0 0 0 1 0 0
-1.6500 -1.4289 0.0000 C 0 0 0 0 0 0 0 0 0 3 0 0
-2.8875 -0.7145 0.0000 F 0 0 0 0 0 0 0 0 0 2 0 0
1 2 1 0 0 0 0
2 3 1 0 0 0 0
3 4 1 0 0 0 0
4 5 1 0 0 0 0
4 6 1 1 0 0 0
M END
";
// molblock_r0 = """ <same heredoc> """
RWMol r0 = RWMol.MolFromMolBlock(molblock_r0);
// r0 = Chem.MolFromMolBlock(molblock_r0)
r0.sanitizeMol();
// Chem.SanitizeMol(r0)
ROMol r0_H = r0.addHs(false);
// r0_H = Chem.AddHs(r0,False)
ROMol_Vect reactants = new ROMol_Vect();
reactants.Add(r0_H);
ROMol_Vect_Vect ps = rxn.runReactants(reactants);
// ps = rxn.RunReactants((r0_H,))
ROMol prod = ps[0][0].removeHs(false);
// prod = Chem.RemoveHs(ps[0][0],False)
prod.compute2DCoords();
// Chem.rdDepictor.Compute2DCoords(prod)
string molblock = prod.MolToMolBlock(true);
// print Chem.MolToMolBlock(prod,True)
string smiles = prod.MolToSmiles(true);
// print Chem.MolToSmiles(prod,True)
Console.WriteLine(molblock);
Console.WriteLine(smiles);
Console.Read();
}
}
}
from rdkit import Chem
from rdkit.Chem import AllChem
molblock_r = """
Mrv0541 12291411382D
6 5 0 0 1 0 999 V2000
-2.8875 0.7145 0.0000 Br 0 0 0 0 0 0 0 0 0 6 0 0
-2.0625 0.7145 0.0000 C 0 0 0 0 0 0 0 0 0 5 0 0
-1.6500 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 4 0 0
-2.0625 -0.7145 0.0000 C 0 0 1 0 0 0 0 0 0 1 0 0
-1.6500 -1.4289 0.0000 C 0 0 0 0 0 0 0 0 0 3 0 0
-2.8875 -0.7145 0.0000 F 0 0 0 0 0 0 0 0 0 2 0 0
1 2 1 0 0 0 0
2 3 1 0 0 0 0
3 4 1 0 0 0 0
4 5 1 0 0 0 0
4 6 1 1 0 0 0
M END"""
molblock_p = """"
Mrv0541 12291411382D
9 8 0 0 1 0 999 V2000
3.1821 0.6188 0.0000 C 0 0 0 0 0 0 0 0 0 3 0 0
3.8966 0.2063 0.0000 C 0 0 2 0 0 0 0 0 0 1 0 0
3.8966 -0.6187 0.0000 F 0 0 0 0 0 0 0 0 0 2 0 0
4.6111 0.6188 0.0000 C 0 0 0 0 0 0 0 0 0 4 0 0
5.3256 0.2062 0.0000 C 0 0 0 0 0 0 0 0 0 5 0 0
6.0400 0.6187 0.0000 C 0 0 0 0 0 0 0 0 0 7 0 0
6.7545 0.2062 0.0000 C 0 0 1 0 0 0 0 0 0 8 0 0
7.4690 0.6187 0.0000 C 0 0 0 0 0 0 0 0 0 9 0 0
6.7545 -0.6188 0.0000 Cl 0 0 0 0 0 0 0 0 0 10 0 0
1 2 1 0 0 0 0
2 3 1 1 0 0 0
2 4 1 0 0 0 0
4 5 1 0 0 0 0
5 6 1 0 0 0 0
7 6 1 0 0 0 0
7 8 1 0 0 0 0
7 9 1 1 0 0 0
M END
"""
rxn = Chem.rdChemReactions.ChemicalReaction()
template_r = Chem.MolFromMolBlock(molblock_r,True,False)
reactant_H = Chem.AddHs(template_r,True)
rxn.AddReactantTemplate(reactant_H)
template_p = Chem.MolFromMolBlock(molblock_p,True,False)
product_H = Chem.AddHs(template_p,True)
rxn.AddProductTemplate(product_H)
rxn.Initialize()
molblock_r0 = """"
Mrv0541 12291411372D
6 5 0 0 1 0 999 V2000
-2.8875 0.7145 0.0000 Br 0 0 0 0 0 0 0 0 0 6 0 0
-2.0625 0.7145 0.0000 C 0 0 0 0 0 0 0 0 0 5 0 0
-1.6500 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 4 0 0
-2.0625 -0.7145 0.0000 C 0 0 1 0 0 0 0 0 0 1 0 0
-1.6500 -1.4289 0.0000 C 0 0 0 0 0 0 0 0 0 3 0 0
-2.8875 -0.7145 0.0000 F 0 0 0 0 0 0 0 0 0 2 0 0
1 2 1 0 0 0 0
2 3 1 0 0 0 0
3 4 1 0 0 0 0
4 5 1 0 0 0 0
4 6 1 1 0 0 0
M END
"""
r0 = Chem.MolFromMolBlock(molblock_r0)
Chem.SanitizeMol(r0)
r0_H = Chem.AddHs(r0,False)
ps = rxn.RunReactants((r0_H,))
prod = Chem.RemoveHs(ps[0][0],False)
Chem.rdDepictor.Compute2DCoords(prod)
print Chem.MolToMolBlock(prod,True)
print Chem.MolToSmiles(prod,True)
------------------------------------------------------------------------------
Dive into the World of Parallel Programming! The Go Parallel Website,
sponsored by Intel and developed in partnership with Slashdot Media, is your
hub for all things parallel software development, from weekly thought
leadership blogs to news, videos, case studies, tutorials and more. Take a
look and join the conversation now. http://goparallel.sourceforge.net
_______________________________________________
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
