Skip to main content

Proposed CDK changes related to PDBReader and BioPolymer

This is expanding on one of the points that Rajarshi made in his blog (which he followed up here) on the PDB file handling capabilities of the CDK. There are two related topics : reading of PDB format files (the ancient, fixed column-width ATOM files) and the model that these are read into.

PDBReader

The old PDB format is being replaced with mmCIF and/or PDBML formats. Only there are lots of programs that write out this format, so it makes sense to still support it for a while at least.

However, it is a quite nasty format, in some ways. Not so much the fixed column width, but the fact that crystallographers abuse the file format in all sorts of ways. Even simple things like expecting that atom numbers will always increase, may not be true.

So it is not easy making a good reader for PDB files. The current CDK one won't read a file with just ATOM records, for example. Think that's reasonable? Well, tough luck for people that made programs that produce simple files like this.

A more serious problem is the fact that you can't get properly connected ligands from a PDB file. Or easily get at the disordered regions. Or get at the waters. Well, sort of - I suppose that many of these things can be done after reading, with CDK classes.

BioPolymer

In some ways - so long as the atoms are read in - anything can be done to the model post-reading. However, the point of having data model classes like Polymer and BioPolymer is to capture some of the complexity of the macromolecule's organisation.

BioPolymer and the PDBReader do a good job of reading and storing the information in the header files (except that it is not always right!). Apart from calling chains 'Strands', some things are done reasonably well. I don't think that ligand atoms should be stored 'loose' as they are, but probably in referenced atom containers.

The real difficulty with modelling proteins lies with representing the hierarchy. One way is 'PMCA' - protein, model, chain, atom. This misses out secondary structure, but it may be too literal to have objects for every concept; the CDK stores the secondary structure given in the header file as IPDBStructure objects - with insertion codes, which is good to see.

However, there are more secondary structures than helix, turn, and strand. I'm never quite sure what the best way to model the more flexible situation, though.

Integration with biojava

It seems a shame that several open source java projects have very little in the way of integration (that I can see). At least on this topic. Of course, Egon (among others) has done work on both both CDK and Jmol, but it concerns me that incompatible ways of doing things lead to projects drifting apart, that should work together.

For example, here is a class I wrote today (for someone else) that uses both biojava and the CDK. It is very much a hack, but the key point is the method makeLigandsFromGroups that takes both a List of biojava.Group objects and a List of IMolecules, along the way converting biojava Atoms to CDK Atoms.

Clearly, biojava has a better interface to its model as you can get a List of the hetatm groups. The CDK, on the other hand, has better support for determining atom types and setting properties on them.

Comments

Popular posts from this blog

Adamantane, Diamantane, Twistane

After cubane, the thought occurred to look at other regular hydrocarbons. If only there was some sort of classification of chemicals that I could use look up similar structures. Oh wate, there is . Anyway, adamantane is not as regular as cubane, but it is highly symmetrical, looking like three cyclohexanes fused together. The vertices fall into two different types when colored by signature: The carbons with three carbon neighbours (degree-3, in the simple graph) have signature (a) and the degree-2 carbons have signature (b). Atoms of one type are only connected to atoms of another - the graph is bipartite . Adamantane connects together to form diamondoids (or, rather, this class have adamantane as a repeating subunit). One such is diamantane , which is no longer bipartite when colored by signature: It has three classes of vertex in the simple graph (a and b), as the set with degree-3 has been split in two. The tree for signature (c) is not shown. The graph is still bipartite accordin

Király's Method for Generating All Graphs from a Degree Sequence

After posting about the Hakimi-Havel  theorem, I received a nice email suggesting various relevant papers. One of these was by Zoltán Király  called " Recognizing Graphic Degree Sequences and Generating All Realizations ". I have now implemented a sketch of the main idea of the paper, which seems to work reasonably well, so I thought I would describe it. See the paper for details, of course. One focus of Király's method is to generate graphs efficiently , by which I mean that it has polynomial delay. In turn, an algorithm with 'polynomial delay' takes a polynomial amount of time between outputs (and to produce the first output). So - roughly - it doesn't take 1s to produce the first graph, 10s for the second, 2s for the third, 300s for the fourth, and so on. Central to the method is the tree that is traversed during the search for graphs that satisfy the input degree sequence. It's a little tricky to draw, but looks something like this: At the top

1,2-dichlorocyclopropane and a spiran

As I am reading a book called "Symmetry in Chemistry" (H. H. Jaffé and M. Orchin) I thought I would try out a couple of examples that they use. One is 1,2-dichlorocylopropane : which is, apparently, dissymmetric because it has a symmetry element (a C2 axis) but is optically active. Incidentally, wedges can look horrible in small structures - this is why: The box around the hydrogen is shaded in grey, to show the effect of overlap. A possible fix might be to shorten the wedge, but sadly this would require working out the bounds of the text when calculating the wedge, which has to be done at render time. Oh well. Another interesting example is this 'spiran', which I can't find on ChEBI or ChemSpider: Image again courtesy of JChempaint . I guess the problem marker (the red line) on the N suggests that it is not a real compound? In any case, some simple code to determine potential chiral centres (using signatures) finds 2 in the cyclopropane structure, and 4 in the