Array-like structure for MolBlock IO

[cwognum]

Hi everyone,

Context: For Polaris, I would like to minimize the number of different file types the storage backend supports. I am therefore exploring whether Zarr could serve as a general purpose data format for saving drug discovery datasets. So far so good, but I'm now running into some challenges w.r.t. saving SDF files, specifically the Atom and Bond Block. For this application, I am not interested in supporting query molecules so that should simplify things.

---

The Atom and Mol Block lend themselves quite naturally to be saved as an array. For example:

Say we have the following SDF file:

     RDKit          3D

 14 15  0  0  0  0  0  0  0  0999 V2000
   -2.1392   -3.1722   -1.1488 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.1302   -2.4867   -0.3430 N   0  0  0  0  0  0  0  0  0  0  0  0
   -0.4692   -2.9563    0.7299 C   0  0  0  0  0  0  0  0  0  0  0  0
    0.3815   -2.0279    1.2013 N   0  0  0  0  0  0  0  0  0  0  0  0
    0.2834   -0.9242    0.4287 C   0  0  0  0  0  0  0  0  0  0  0  0
   -0.6660   -1.2364   -0.5269 C   0  0  0  0  0  0  0  0  0  0  0  0
   -0.9777   -0.2738   -1.4725 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.8413   -0.5546   -2.3437 O   0  0  0  0  0  0  0  0  0  0  0  0
   -0.3474    0.9114   -1.4182 N   0  0  0  0  0  0  0  0  0  0  0  0
   -0.6901    1.9046   -2.4192 C   0  0  0  0  0  0  0  0  0  0  0  0
    0.5817    1.1832   -0.4629 C   0  0  0  0  0  0  0  0  0  0  0  0
    1.1511    2.3087   -0.4465 O   0  0  0  0  0  0  0  0  0  0  0  0
    0.9066    0.2690    0.4698 N   0  0  0  0  0  0  0  0  0  0  0  0
    1.9089    0.5829    1.4879 C   0  0  0  0  0  0  0  0  0  0  0  0
  1  2  1  0
  2  3  1  0
  3  4  2  0
  4  5  1  0
  5  6  2  0
  6  7  1  0
  7  8  2  0
  7  9  1  0
  9 10  1  0
  9 11  1  0
 11 12  2  0
 11 13  1  0
 13 14  1  0
  6  2  1  0
 13  5  1  0
M  END
$$$$

We can save this as two arrays:

# For the atom block, we convert the atom symbol to its atomic number
atom_block = [
    [-2.1392, -3.1722, -1.1488, 6, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    [-1.1302, -2.4867, -0.3430, 7, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
    ...
    [1.9089, 0.5829, 1.4879, 6, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
]

bond_block = [
    [1, 2, 1, 0],
    [2, 3, 1, 0],
    ...
    [13, 5, 1, 0]
]

As long as we save the Chiral Flag from the counts line, I believe we would not lose any information this way?

Does RDKit provide an easy way to get all the properties that make up the atom and bond block? I only found rdkit.Chem.rdmolfiles.MolFromMolBlock and rdkit.Chem.rdmolfiles.MolToMolBlock in the Python binding, but these operate on strings. I could:

• Parse the SDF string myself, but that seems hacky.
• I could also work my through the RDKit source code to duplicate the necessary code, but that too does not feel like good practice.

Curious to hear if anyone in the community has any thoughts on this.

[cwognum]

Alternatively, we could first convert the mol to an entirely different molecular representation and save this alongside the 3D coordinates separately? But I don't want to lose the information stored in the atom and mol block flags. Maybe rdkit.Chem.rdchem.Mol.ToBinary, although it seems this could change between RDKit versions (see #4919) ? Maybe Inchi?

[greglandrum]

Hi @cwognum, if I were doing this I would have a single column in which I store the V3000 mol block for each molecule.
If you want to enable more efficient processing with the RDKit, you could also have a column with the output of Mol.ToBinary() (that can change from version to version, but we maintain backwards compatibility).

Storing arrays of atoms and bonds is limiting, adds complication, requires you to write your own code to serialize/deserialize the molecules into the array format, and delivers very little additional value - there aren't very many use cases where you'd want to query those features and if you need that, you can always add code that creates those arrays while still storing the molecules themselves in more portable ways

[cwognum]

Hi @greglandrum, thank you for the quick and thoughtful answer!

This all makes sense to me. Just to make sure I understood correctly: There's no loss of information when using Mol.ToBinary()? We can create the exact same Mol object from it?

One reason I was curious to explore arrays is because I was wondering whether it could lead to more efficient storage usage due to the compressors that Zarr supports (e.g. since a lot of the flags in the atom block are 0, any good compression algorithm would compress that drastically). But come to think of it, I assume Mol.ToBinary() has already been quite effectively optimized. I will run some tests to compare storage requirements! Thanks again!

[greglandrum]

Hi @greglandrum, thank you for the quick and thoughtful answer!

This all makes sense to me. Just to make sure I understood correctly: There's no loss of information when using Mol.ToBinary()? We can create the exact same Mol object from it?

One reason I was curious to explore arrays is because I was wondering whether it could lead to more efficient storage usage due to the compressors that Zarr supports (e.g. since a lot of the flags in the atom block are 0, any good compression algorithm would compress that drastically). But come to think of it, I assume Mol.ToBinary() has already been quite effectively optimized. I will run some tests to compare storage requirements! Thanks again!

This all makes sense to me. Just to make sure I understood correctly: There's no loss of information when using Mol.ToBinary()? We can create the exact same Mol object from it?

That is the intent. It's not a direct binary dump of the molecule, but going to/from that binary representation should be mostly lossless. There are some details depending on whether or not you choose to preserve atom/bond/mol properties (and whether or not you have properties you need to preserve), but in general you should be fine.

One reason I was curious to explore arrays is because I was wondering whether it could lead to more efficient storage usage due to the compressors that Zarr supports (e.g. since a lot of the flags in the atom block are 0, any good compression algorithm would compress that drastically). But come to think of it, I assume Mol.ToBinary() has already been quite effectively optimized. I will run some tests to compare storage requirements! Thanks again!

Though you may win with storage, it's unlikely to be dramatic and you are definitely going to lose with performance. Serializing/Deserializing molecules using the binary format is going to be MUCH faster.

[cwognum]

Thank you once more! This has been very helpful!

[cwognum]

There are some details depending on whether or not you choose to preserve atom/bond/mol properties (and whether or not you have properties you need to preserve), but in general you should be fine.

@greglandrum One more question came up. The properties you're referring to here ☝️ , I found this bit of code to specify this, but cannot find any documentation on what this means exactly. I already have a way to save any properties stored in the SDF, e.g.:

>  <smiles>  (1) 
COC(=O)[C@@]1(Cc2ccc(OC)cc2)[C@H]2c3cc(C(=O)N(C)C)n(Cc4ccc(OC(F)(F)F)cc4)c3C[C@H]2CN1C(=O)c1ccccc1

If I want to preserve as much info as possible, is this necessary?

bytes_data = [mol.ToBinary(Chem.PropertyPickleOptions.AllProps) for mol in mols]

[greglandrum]

I would suggest pulling out the molecular properties and storing them as a separate column (or cell or whatever) in the database. It’s easy to imagine want to query on those properties.

If the system you’re using supports any kind of key-value data type (or even just JSON), this should be quite easy

[cwognum]

Yes, we would support that. Given that I'm converting directly from a SDF file, would the mol properties be just the ones defined in the SDF file (after the bond block)?