A while ago I had a project which needed PMI descriptors (specifically NPR1
and NPR2) which were not available in the main branch of RDKit at the time.
At the time I used the fork by 'hahnda6' which provided the
calcPMIDescriptors() function, and this worked well. Now that PMI
descriptors are available in the main RDKit distrubution I thought I'd
rewrite my code to use the official version.

Building the new RDKit was no problem, but things went downhill shortly
after that. There's every chance that I've missed the relevant
documentation (I hope someone can point me in the right direction if so)
and done something stupid!

The issues are -
1) I can't find any documentation of the C++ API - the only reference to
PMI in the online RDKit documentation appears to be to the PMI.h file
2) Having written a program using the PMI[123] and/or NPR[12] functions, I
couldn't get it to compile until I added the  -DRDK_BUILD_DESCRIPTORS3D
directive -
g++ -o sdf_pmi_blob sdf_pmi.cpp -I/packages/rdkit/include/rdkit
-L/packages/rdkit/lib -lDescriptors -lGraphMol -lFileParsers
This seems a bit odd...
3) Is it necessary to make separate calls to the individual PMI() and/or
NPR() functions? Surely this results in duplication of some of the heavier
calculations? I can't find any equivalent of calcPMIDescriptors() which
returned a 'Moments' struct containing all the PMI and NPR values in one go.
4) The big one! The returned results look very odd. They appear to relate
more to the dimensions of the molecule than the moments of inertia. For a
rod-like molecule (dimethylacetylene) I'd expect two large and one small
PMI (e.g. PMI1: 6.61651   PMI2: 150.434   PMI3: 150.434  NPR1: 0.0439828
NPR2: 0.999998) but actually get PMI1: 0.061647  PMI2: 0.061652  PMI3:
25.3699  NPR1: 0.002430  NPR2: 0.002430.
For disk-like (benzene) the result should be one large and two medium (e.g.
PMI1: 89.1448  PMI2: 89.1495  PMI3: 178.294  NPR1: 0.499987  NPR2:
0.500013) but get PMI1: 2.37457e-10  PMI2: 11.0844  PMI3: 11.0851  NPR1:
2.14213e-11  NPR2: 0.999933.
Finally for a roughly spherical molecule (neopentane) the NPR values look
reasonable (no great surprise) but the absolute PMI values may be too
small: old program - PMI1: 114.795  PMI2: 114.797  PMI3: 114.799
NPR1: 0.999966  NPR2: 0.999988, new program - PMI1: 6.59466  PMI2: 6.59488
PMI3: 6.59531  NPR1: 0.999902  NPR2: 0.999935

As I say, it's entirely likely that I'm doing something stupid here so any
pointers will be gratefully received. FWIW, the core of my program is -
mol = MolBlockToMol(ctab, true, false);
double pmi1 = RDKit::Descriptors::PMI1(*mol);
double pmi2 = RDKit::Descriptors::PMI2(*mol);
double pmi3 = RDKit::Descriptors::PMI3(*mol);
double npr1 = RDKit::Descriptors::NPR1(*mol);
double npr2 = RDKit::Descriptors::NPR2(*mol);

Thanks for any help!
Check out the vibrant tech community on one of the world's most 
engaging tech sites, SlashDot.org! http://sdm.link/slashdot
Rdkit-discuss mailing list

Reply via email to