Predicting Molecular Properties With Machine Learning

Property prediction is regression on molecules, and the modelling is rarely the hard part. A four-term linear model published in 2004 is still a respectable baseline for aqueous solubility, which tells you where the difficulty actually lives: in the split, the label noise and the domain.

Descriptors: the smallest useful input

A descriptor is a number computed from the structure. The classical set is small and interpretable: molecular weight, computed logP (the octanol-water partition coefficient, a lipophilicity measure), topological polar surface area, counts of hydrogen-bond donors and acceptors, rotatable bonds, ring count, fraction of sp3 carbons. RDKit computes a couple of hundred of them.

They are not independent, and that matters for anything linear. Molecular weight, heavy-atom count and molar refractivity are close to collinear across most libraries, so a regression given all three distributes one effect across three unstable coefficients and the resulting signs are uninterpretable. Either select a small non-redundant set, as the model below does, or use a method that tolerates correlated inputs and stop reading the coefficients as chemistry.

Descriptors and fingerprints answer different questions and are usually concatenated. A fingerprint says which substructures are present; a descriptor says what the whole molecule is like. A model given only fingerprints has to reconstruct molecular weight from bit counts, badly.

Worked: solubility of benzene from four numbers

John Delaney fitted a linear model to 2,874 measured aqueous solubilities and published it as ESOL in the Journal of Chemical Information and Computer Sciences in 2004. Four descriptors survived: computed logP, molecular weight, the number of rotatable bonds, and the aromatic proportion — aromatic heavy atoms divided by total heavy atoms.

log S  =  0.16
        - 0.63  * clogP
        - 0.0062 * MW
        + 0.066  * rotatable bonds
        - 0.74   * aromatic proportion

benzene:  clogP 2.13,  MW 78.11,  rotatable bonds 0,  AP = 6/6 = 1.0

  0.16
- 0.63 * 2.13    = -1.342
- 0.0062 * 78.11 = -0.484
+ 0.066 * 0      =  0.000
- 0.74 * 1.0     = -0.740
                 ---------
  log S          = -2.41       (S in mol/L)

measured: benzene dissolves at about 1.8 g/L
          1.8 / 78.11 = 0.023 mol/L  ->  log S = -1.64

error = 0.77 log units, i.e. the prediction is about 6x too insoluble

That gap is representative rather than embarrassing: Delaney reports average absolute errors around three-quarters of a log unit for this equation on held-out sets. Four numbers, no training run, and you are within a factor of ten of the truth on a property that spans ten orders of magnitude. Any modern model has to be judged against that, not against zero.

Random splits lie

Chemical datasets are full of analogue series: twenty compounds from one paper differing by a substituent. Split them randomly and the test set contains near-neighbours of training compounds, so the model is rewarded for interpolating within a series — which is not what you will ask it to do.

The standard alternative is the scaffold split, which groups molecules by their Bemis-Murcko framework — the ring systems plus the linkers between them, with side chains stripped, as defined by Guy Bemis and Mark Murcko in 1996 — and assigns whole scaffold groups to train or test. No scaffold appears on both sides. The reported error rises, sometimes by a factor of two, and the higher number is the one that predicts how the model behaves on a new chemical series.

MoleculeNet, the benchmark suite from Zhenqin Wu and colleagues published in Chemical Science in 2018, packages ESOL, FreeSolv, Lipophilicity, BBBP, Tox21 and others with recommended splits. Use its splits when you want to compare against published numbers, and state which split you used when you report yours; a random-split number and a scaffold-split number for the same model are not comparable and are routinely quoted as if they were.

The error floor is experimental

Aqueous solubility measurements disagree between laboratories. David Palmer and John Mitchell examined this directly in “Is experimental data quality the limiting factor in predicting the aqueous solubility of druglike molecules?” (Molecular Pharmaceutics, 2014), and concluded that the reproducibility of the experiment itself sets a floor of the same order as the best models’ errors.

Two things follow. A model reporting a test RMSE well below the inter-laboratory spread of its own label is not better than the experiment; it is fitting the idiosyncrasies of one data source, which usually means the split leaked. And chasing the last 0.1 log units of accuracy is chasing noise — the return is in widening the domain the model covers, not in tightening the residual on what it already covers.

Report RMSE in the units of the label alongside the label’s own spread. A coefficient of determination on its own hides both, because it is a ratio against the variance of your particular test set: assemble a test set spanning fifteen log units and a mediocre model scores well.

The same argument decides what is worth predicting at all. Endpoints measured cleanly and often — solubility, lipophilicity, permeability in a standard assay — have enough consistent data for a model to be useful. Endpoints assembled from heterogeneous literature, where the protocol varies between rows, have a label whose noise is larger than the effect being modelled, and no architecture recovers a signal that the aggregation destroyed. Check the provenance of the label column before choosing the model class.

Applicability domain

A regression model always returns a number. It has no mechanism for saying “this molecule is unlike anything I was trained on” unless you build one.

• Distance to training set. Compute the maximum Tanimoto similarity from the query to any training molecule. Below roughly 0.3 to 0.4 on ECFP4, treat the prediction as an extrapolation and say so in the output.
• Ensemble spread. Train five models on different splits or with different seeds and report the standard deviation of their predictions. Wide spread flags a molecule the training data does not constrain.
• Descriptor range. A molecule whose molecular weight or logP falls outside the training range is outside the domain of a linear model by construction, whatever the model outputs.

This matters most where the prediction feeds a decision. A solubility model that quietly extrapolates is a nuisance; a toxicity or ADME-endpoint model that quietly extrapolates is a claim about safety that nobody made. Predictions on those endpoints are hypotheses for experimental and expert review, and regulatory use of computational toxicology carries its own validation requirements that a model score does not satisfy.