Small molecules in ChEMBL database

Series 1.1.2 - Preprocessing data in Polars dataframe library

Machine learning projects
Polars
Python
Jupyter
ChEMBL database
Cheminformatics
Author

Jennifer HY Lin

Published

January 4, 2023

Modified

July 7, 2025

Import Polars and read parquet file

This second post will mostly be about doing some data preprocessing for the parquet file we’ve saved from the first post, trying to get the data in a reasonable state before moving onto the machine learning model building and training phase. The first step is to import Polars dataframe library and then read the parquet file as a dataframe.

import polars as pl
df_pa = pl.read_parquet("chembl_sm_mols.parquet")
df_pa
shape: (1_831_560, 32)
ChEMBL ID Name Synonyms Type Max Phase Molecular Weight Targets Bioactivities AlogP Polar Surface Area HBA HBD #RO5 Violations #Rotatable Bonds Passes Ro3 QED Weighted CX Acidic pKa CX Basic pKa CX LogP CX LogD Aromatic Rings Structure Type Inorganic Flag Heavy Atoms HBA (Lipinski) HBD (Lipinski) #RO5 Violations (Lipinski) Molecular Weight (Monoisotopic) Molecular Species Molecular Formula Smiles Inchi Key
str str str str i64 f64 i64 i64 f64 f64 i64 i64 i64 i64 str f64 f64 f64 f64 f64 i64 str i64 i64 i64 i64 i64 f64 str str str str
"CHEMBL539070" null null "Small molecule" 0 286.79 1 1 2.28 73.06 6 2 0 5 "N" 0.63 13.84 3.64 2.57 2.57 2 "MOL" -1 17 5 3 0 250.0888 "NEUTRAL" "C11H15ClN4OS" "CCCOc1ccccc1-c1nnc(NN)s1.Cl" "WPEWNRKLKLNLSO-UHFFFAOYSA-N"
"CHEMBL3335528" null null "Small molecule" 0 842.8 2 6 0.18 269.57 18 5 2 17 "N" 0.09 3.2 null 3.31 -0.14 3 "MOL" -1 60 19 5 2 842.2633 "ACID" "C41H46O19" "COC(=O)[C@H](O[C@@H]1O[C@@H](C… "KGUJQZWYZPYYRZ-LWEWUKDVSA-N"
"CHEMBL2419030" null null "Small molecule" 0 359.33 4 4 3.94 85.13 6 1 0 3 "N" 0.66 null null 3.66 3.66 2 "MOL" -1 24 6 1 0 359.0551 "NEUTRAL" "C14H12F3N3O3S" "O=c1nc(NC2CCCC2)sc2c([N+](=O)[… "QGDMYSDFCXOKML-UHFFFAOYSA-N"
"CHEMBL3827271" null null "Small molecule" 0 712.85 1 1 -2.84 319.06 10 11 2 16 "N" 0.07 4.08 10.49 -6.88 -8.95 0 "MOL" -1 50 19 14 3 712.4232 "ZWITTERION" "C31H56N10O9" "CC(C)C[C@@H]1NC(=O)[C@H](CCCNC… "QJQNNLICZLLPMB-VUBDRERZSA-N"
"CHEMBL3465961" null null "Small molecule" 0 319.42 16 22 2.22 50.5 4 1 0 6 "N" 0.87 null 9.38 2.13 -0.44 1 "MOL" -1 23 4 1 0 319.206 "BASE" "C18H26FN3O" "CC(O)CN1CCC(CN(C)Cc2cc(C#N)ccc… "FZEVYCHTADTXPM-UHFFFAOYSA-N"
… … … … … … … … … … … … … … … … … … … … … … … … … … … … … … … …
"CHEMBL2017916" null null "Small molecule" 0 312.35 3 3 2.86 77.0 6 1 0 4 "N" 0.8 8.13 3.49 2.17 2.1 3 "MOL" -1 22 6 1 0 312.0681 "NEUTRAL" "C15H12N4O2S" "COc1ccc(-c2nnc(NC(=O)c3cccnc3)… "XIZUJGDKNPVNQA-UHFFFAOYSA-N"
"CHEMBL374652" null null "Small molecule" 0 403.83 1 1 5.98 36.02 2 2 1 4 "N" 0.42 13.65 null 5.36 5.36 3 "MOL" -1 26 2 2 1 403.0421 "NEUTRAL" "C18H14ClF4NOS" "CC(O)(CSc1ccc(F)cc1)c1cc2cc(Cl… "CRPQTBRTHURKII-UHFFFAOYSA-N"
"CHEMBL1416264" null null "Small molecule" 0 380.41 6 8 3.06 85.07 7 1 0 5 "N" 0.54 13.85 3.86 2.47 2.47 4 "MOL" -1 27 7 1 0 380.0856 "NEUTRAL" "C18H13FN6OS" "O=C(CSc1ccc2nnc(-c3cccnc3)n2n1… "QVYIEKHEJKFNAT-UHFFFAOYSA-N"
"CHEMBL213734" null null "Small molecule" 0 288.26 2 3 2.32 101.7 5 2 0 5 "N" 0.5 7.2 null 2.36 1.95 2 "MOL" -1 21 7 2 0 288.0746 "NEUTRAL" "C14H12N2O5" "O=C(COc1ccccc1)Nc1ccc([N+](=O)… "PZTWAHGBGTWVEB-UHFFFAOYSA-N"
"CHEMBL1531634" null null "Small molecule" 0 320.16 19 21 4.4 29.1 2 1 0 4 "N" 0.67 null null 4.04 4.04 2 "MOL" -1 19 2 1 0 319.0008 "NEUTRAL" "C15H11BrFNO" "O=C(/C=C/Nc1ccc(F)cc1)c1ccc(Br… "DKPWCCDDKFLKEC-MDZDMXLPSA-N"

I’m having a look at the max phase column first to see the distribution or counts of molecules in each max phase.

df_pa.group_by("Max Phase").len()
shape: (5, 2)
Max Phase len
i64 u32
3 773
0 1826345
1 520
4 2870
2 1052

I also want to find out what types of physicochemical properties are there in this dataset.

# Print all column names and data types 
print(df_pa.glimpse())
Rows: 1831560
Columns: 32
$ ChEMBL ID                       <str> 'CHEMBL539070', 'CHEMBL3335528', 'CHEMBL2419030', 'CHEMBL3827271', 'CHEMBL3465961', 'CHEMBL3824158', 'CHEMBL194112', 'CHEMBL2047226', 'CHEMBL1991010', 'CHEMBL195644'
$ Name                            <str> None, None, None, None, None, None, None, None, None, None
$ Synonyms                        <str> None, None, None, None, None, None, None, None, None, None
$ Type                            <str> 'Small molecule', 'Small molecule', 'Small molecule', 'Small molecule', 'Small molecule', 'Small molecule', 'Small molecule', 'Small molecule', 'Small molecule', 'Small molecule'
$ Max Phase                       <i64> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0
$ Molecular Weight                <f64> 286.79, 842.8, 359.33, 712.85, 319.42, 422.48, 366.38, 452.4, 454.05, 375.47
$ Targets                         <i64> 1, 2, 4, 1, 16, 2, 2, 4, 60, 2
$ Bioactivities                   <i64> 1, 6, 4, 1, 22, 4, 3, 8, 60, 3
$ AlogP                           <f64> 2.28, 0.18, 3.94, -2.84, 2.22, 5.09, 4.8, 4.93, 5.18, 4.95
$ Polar Surface Area              <f64> 73.06, 269.57, 85.13, 319.06, 50.5, 109.54, 57.53, 53.08, 40.54, 70.42
$ HBA                             <i64> 6, 18, 6, 10, 4, 6, 3, 5, 3, 4
$ HBD                             <i64> 2, 5, 1, 11, 1, 2, 2, 2, 1, 2
$ #RO5 Violations                 <i64> 0, 2, 0, 2, 0, 1, 0, 0, 1, 0
$ #Rotatable Bonds                <i64> 5, 17, 3, 16, 6, 10, 1, 7, 8, 2
$ Passes Ro3                      <str> 'N', 'N', 'N', 'N', 'N', 'N', 'N', 'N', 'N', 'N'
$ QED Weighted                    <f64> 0.63, 0.09, 0.66, 0.07, 0.87, 0.31, 0.75, 0.53, 0.6, 0.73
$ CX Acidic pKa                   <f64> 13.84, 3.2, None, 4.08, None, 4.59, 8.98, None, 13.88, 9.52
$ CX Basic pKa                    <f64> 3.64, None, None, 10.49, 9.38, 7.99, None, 8.47, 8.48, 3.73
$ CX LogP                         <f64> 2.57, 3.31, 3.66, -6.88, 2.13, 2.49, 4.84, 4.51, 6.34, 3.92
$ CX LogD                         <f64> 2.57, -0.14, 3.66, -8.95, -0.44, 2.42, 4.83, 3.29, 5.22, 3.91
$ Aromatic Rings                  <i64> 2, 3, 2, 0, 1, 2, 1, 3, 2, 2
$ Structure Type                  <str> 'MOL', 'MOL', 'MOL', 'MOL', 'MOL', 'MOL', 'MOL', 'MOL', 'MOL', 'MOL'
$ Inorganic Flag                  <i64> -1, -1, -1, -1, -1, -1, -1, -1, -1, -1
$ Heavy Atoms                     <i64> 17, 60, 24, 50, 23, 31, 26, 29, 31, 28
$ HBA (Lipinski)                  <i64> 5, 19, 6, 19, 4, 7, 3, 5, 3, 4
$ HBD (Lipinski)                  <i64> 3, 5, 1, 14, 1, 2, 2, 2, 1, 2
$ #RO5 Violations (Lipinski)      <i64> 0, 2, 0, 3, 0, 1, 0, 0, 1, 0
$ Molecular Weight (Monoisotopic) <f64> 250.0888, 842.2633, 359.0551, 712.4232, 319.206, 422.1842, 366.1443, 451.1372, 417.2668, 375.1834
$ Molecular Species               <str> 'NEUTRAL', 'ACID', 'NEUTRAL', 'ZWITTERION', 'BASE', 'ACID', 'NEUTRAL', 'NEUTRAL', 'NEUTRAL', 'NEUTRAL'
$ Molecular Formula               <str> 'C11H15ClN4OS', 'C41H46O19', 'C14H12F3N3O3S', 'C31H56N10O9', 'C18H26FN3O', 'C24H26N2O5', 'C20H21F3O3', 'C23H26BrN5', 'C28H36ClNO2', 'C24H25NO3'
$ Smiles                          <str> 'CCCOc1ccccc1-c1nnc(NN)s1.Cl', 'COC(=O)[C@H](O[C@@H]1O[C@@H](C)[C@@H](O)[C@@H](O)[C@@H]1O)[C@@H](O[C@@H]1O[C@H](CO)[C@H](OC(=O)c2ccccc2)[C@H](O[C@H](Cc2ccccc2)C(=O)O)[C@H]1OC(=O)c1ccccc1)C(=O)OC', 'O=c1nc(NC2CCCC2)sc2c([N+](=O)[O-])cc(C(F)(F)F)cc12', 'CC(C)C[C@@H]1NC(=O)[C@H](CCCNC(N)=O)NC(=O)[C@H](CCCCN)NC(=O)[C@H](CC(=O)O)NC(=O)[C@H](CCCCN)NC(=O)CCNC1=O', 'CC(O)CN1CCC(CN(C)Cc2cc(C#N)ccc2F)CC1', 'CCCCCCCNC(C1=C(O)C(=O)c2ccccc2C1=O)c1ccc([N+](=O)[O-])cc1', 'C[C@]12CCC3c4ccc(O)cc4CCC3C1CC(C(=O)C(F)(F)F)=C2O', 'Brc1ccc(CNc2ccnc(NC3CCN(Cc4ccccc4)CC3)n2)cc1', 'CCc1ccc(/C=C/C(=O)C2CN(CC)CCC2(O)/C=C/c2ccc(CC)cc2)cc1.Cl', 'C[C@]12CCC3c4ccc(O)cc4CCC3C1CC(C(=O)c1cccnc1)=C2O'
$ Inchi Key                       <str> 'WPEWNRKLKLNLSO-UHFFFAOYSA-N', 'KGUJQZWYZPYYRZ-LWEWUKDVSA-N', 'QGDMYSDFCXOKML-UHFFFAOYSA-N', 'QJQNNLICZLLPMB-VUBDRERZSA-N', 'FZEVYCHTADTXPM-UHFFFAOYSA-N', 'AXOVDUYYBUYLPC-UHFFFAOYSA-N', 'FIBOSLUEJGPVMK-RYCRIANLSA-N', 'WOAVNWHCIXCOIZ-UHFFFAOYSA-N', 'XJDPAUYFONOZBC-DCPGAFKKSA-N', 'MOBPUUUBXAHZBM-KSAYNYSMSA-N'

None


Looking at some of the physicochemical properties in the dataset

I’ve gone through the ChEMBL_31 schema documentation and ChEMBL database website to find out the meanings for some of the column names. The explanations of the selected term are adapted from ChEMBL_31 schema documentation (available as “Release notes” on the website at the time), or if definitions for certain terms are not available, I resort to interpret them myself by going into “Dinstict compounds” section of the ChEMBL database.

The definitions for some of the listed physicochemical properties are:

Max Phase - Maximum phase of development reached for the compound (where 4 = approved). Null is where max phase has not yet been assigned.

Bioactivities - Various biological assays used for the compounds e.g. IC50, GI50, potency tests etc.

AlogP - Calculated partition coefficient

HBA - Number of hydrogen bond acceptors

HBD - Number of hydrogen bond donors

#RO5 Violations - Number of violations of Lipinski’s rule-of-five, using HBA and HBD definitions

Passes Ro3 - Indicating whether the compound passed the rule-of-three (MW < 300, logP < 3 etc)

QED Weighted - Weighted quantitative estimate of drug likeness (Bickerton et al. 2012)

Inorganic flag - Indicating whether the molecule is inorganic (i.e., containing only metal atoms and <2 carbon atoms), where 1 = inorganic compound and 0 = non-inorganic compound (-1 meaning preclinical compound or not a drug)

Heavy Atoms - Number of heavy (non-hydrogen) atoms

CX Acidic pKa - The most acidic pKa calculated using ChemAxon v17.29.0

CX Basic pKa - The most basic pKa calculated using ChemAxon v17.29.0

CX LogP - The calculated octanol/water partition coefficient using ChemAxon v17.29.0

CX LogD - The calculated octanol/water distribution coefficient at pH = 7.4 using ChemAxon v17.29.0

Structure Type - based on compound_structures table, where SEQ indicates an entry in the protein_therapeutics table instead, NONE indicates an entry in neither tables, e.g. structure unknown

Inchi Key - the IUPAC international chemical identifier key

The older version of this post has had a code section for changing the data types of several columns. I’ve found out that this is no longer needed as I’ve taken care of this in the first post when I’ve specifically stated the null_values (e.g. “None” and ““) to be converted to”null” during pl.read_csv(). When the same dataframe df_pa is read here, you’ll notice that the data type for each column should be how it should be, matching the data inside each column.


Dealing with nulls

I’m using null_count() to see the distributions of all null entries in the dataset.

# Alternative code: df_pa.select(pl.all().null_count())
df_pa.null_count()
shape: (1, 32)
ChEMBL ID Name Synonyms Type Max Phase Molecular Weight Targets Bioactivities AlogP Polar Surface Area HBA HBD #RO5 Violations #Rotatable Bonds Passes Ro3 QED Weighted CX Acidic pKa CX Basic pKa CX LogP CX LogD Aromatic Rings Structure Type Inorganic Flag Heavy Atoms HBA (Lipinski) HBD (Lipinski) #RO5 Violations (Lipinski) Molecular Weight (Monoisotopic) Molecular Species Molecular Formula Smiles Inchi Key
u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32 u32
0 1795946 1753154 0 0 0 0 0 26426 26426 26426 26426 26426 26426 26426 26426 783038 691442 26560 26560 26426 0 0 26426 26426 26426 26426 0 43392 0 0 0

The following is basically several different ways to remove null values from the original dataset.

## ---Drop all rows with null entries---
#df_pa.drop_nulls()

## Number of rows reduced to 2,645 - seems way too much...

## ---Restricting drop nulls to only subsets of nulls in strings ---
# import polars.selectors as cs
# df_pa.drop_nulls(subset=cs.string())

## Number of rows reduced to 17,735 - hmm... try other ways?

## ---Drop a row if all values are null---
# df_pa.filter(~pl.all_horizontal(pl.all().is_null()))

## No change in number of rows - meaning no one row contains all null values

## ---Restricting drop nulls to a specific column---
# df_pa.drop_nulls(subset="CX LogP")

## Number of rows reduced to 1,873,678

Initially, I’ve tried to remove all nulls first, then I realise that removing nulls above may not be the best way to prepare the data as all max phase 4 compounds are actually also removed completely! So what I’m going to do instead is to keep them by using fill_null() to replace all “null” as “0” (this’ll only apply to all the integers or floats in the data).

df_pa = df_pa.fill_null(0)
df_pa
shape: (1_831_560, 32)
ChEMBL ID Name Synonyms Type Max Phase Molecular Weight Targets Bioactivities AlogP Polar Surface Area HBA HBD #RO5 Violations #Rotatable Bonds Passes Ro3 QED Weighted CX Acidic pKa CX Basic pKa CX LogP CX LogD Aromatic Rings Structure Type Inorganic Flag Heavy Atoms HBA (Lipinski) HBD (Lipinski) #RO5 Violations (Lipinski) Molecular Weight (Monoisotopic) Molecular Species Molecular Formula Smiles Inchi Key
str str str str i64 f64 i64 i64 f64 f64 i64 i64 i64 i64 str f64 f64 f64 f64 f64 i64 str i64 i64 i64 i64 i64 f64 str str str str
"CHEMBL539070" null null "Small molecule" 0 286.79 1 1 2.28 73.06 6 2 0 5 "N" 0.63 13.84 3.64 2.57 2.57 2 "MOL" -1 17 5 3 0 250.0888 "NEUTRAL" "C11H15ClN4OS" "CCCOc1ccccc1-c1nnc(NN)s1.Cl" "WPEWNRKLKLNLSO-UHFFFAOYSA-N"
"CHEMBL3335528" null null "Small molecule" 0 842.8 2 6 0.18 269.57 18 5 2 17 "N" 0.09 3.2 0.0 3.31 -0.14 3 "MOL" -1 60 19 5 2 842.2633 "ACID" "C41H46O19" "COC(=O)[C@H](O[C@@H]1O[C@@H](C… "KGUJQZWYZPYYRZ-LWEWUKDVSA-N"
"CHEMBL2419030" null null "Small molecule" 0 359.33 4 4 3.94 85.13 6 1 0 3 "N" 0.66 0.0 0.0 3.66 3.66 2 "MOL" -1 24 6 1 0 359.0551 "NEUTRAL" "C14H12F3N3O3S" "O=c1nc(NC2CCCC2)sc2c([N+](=O)[… "QGDMYSDFCXOKML-UHFFFAOYSA-N"
"CHEMBL3827271" null null "Small molecule" 0 712.85 1 1 -2.84 319.06 10 11 2 16 "N" 0.07 4.08 10.49 -6.88 -8.95 0 "MOL" -1 50 19 14 3 712.4232 "ZWITTERION" "C31H56N10O9" "CC(C)C[C@@H]1NC(=O)[C@H](CCCNC… "QJQNNLICZLLPMB-VUBDRERZSA-N"
"CHEMBL3465961" null null "Small molecule" 0 319.42 16 22 2.22 50.5 4 1 0 6 "N" 0.87 0.0 9.38 2.13 -0.44 1 "MOL" -1 23 4 1 0 319.206 "BASE" "C18H26FN3O" "CC(O)CN1CCC(CN(C)Cc2cc(C#N)ccc… "FZEVYCHTADTXPM-UHFFFAOYSA-N"
… … … … … … … … … … … … … … … … … … … … … … … … … … … … … … … …
"CHEMBL2017916" null null "Small molecule" 0 312.35 3 3 2.86 77.0 6 1 0 4 "N" 0.8 8.13 3.49 2.17 2.1 3 "MOL" -1 22 6 1 0 312.0681 "NEUTRAL" "C15H12N4O2S" "COc1ccc(-c2nnc(NC(=O)c3cccnc3)… "XIZUJGDKNPVNQA-UHFFFAOYSA-N"
"CHEMBL374652" null null "Small molecule" 0 403.83 1 1 5.98 36.02 2 2 1 4 "N" 0.42 13.65 0.0 5.36 5.36 3 "MOL" -1 26 2 2 1 403.0421 "NEUTRAL" "C18H14ClF4NOS" "CC(O)(CSc1ccc(F)cc1)c1cc2cc(Cl… "CRPQTBRTHURKII-UHFFFAOYSA-N"
"CHEMBL1416264" null null "Small molecule" 0 380.41 6 8 3.06 85.07 7 1 0 5 "N" 0.54 13.85 3.86 2.47 2.47 4 "MOL" -1 27 7 1 0 380.0856 "NEUTRAL" "C18H13FN6OS" "O=C(CSc1ccc2nnc(-c3cccnc3)n2n1… "QVYIEKHEJKFNAT-UHFFFAOYSA-N"
"CHEMBL213734" null null "Small molecule" 0 288.26 2 3 2.32 101.7 5 2 0 5 "N" 0.5 7.2 0.0 2.36 1.95 2 "MOL" -1 21 7 2 0 288.0746 "NEUTRAL" "C14H12N2O5" "O=C(COc1ccccc1)Nc1ccc([N+](=O)… "PZTWAHGBGTWVEB-UHFFFAOYSA-N"
"CHEMBL1531634" null null "Small molecule" 0 320.16 19 21 4.4 29.1 2 1 0 4 "N" 0.67 0.0 0.0 4.04 4.04 2 "MOL" -1 19 2 1 0 319.0008 "NEUTRAL" "C15H11BrFNO" "O=C(/C=C/Nc1ccc(F)cc1)c1ccc(Br… "DKPWCCDDKFLKEC-MDZDMXLPSA-N" 
# Summary statistics for df_dn dataset
df_pa.describe()
shape: (9, 33)
statistic ChEMBL ID Name Synonyms Type Max Phase Molecular Weight Targets Bioactivities AlogP Polar Surface Area HBA HBD #RO5 Violations #Rotatable Bonds Passes Ro3 QED Weighted CX Acidic pKa CX Basic pKa CX LogP CX LogD Aromatic Rings Structure Type Inorganic Flag Heavy Atoms HBA (Lipinski) HBD (Lipinski) #RO5 Violations (Lipinski) Molecular Weight (Monoisotopic) Molecular Species Molecular Formula Smiles Inchi Key
str str str str str f64 f64 f64 f64 f64 f64 f64 f64 f64 f64 str f64 f64 f64 f64 f64 f64 str f64 f64 f64 f64 f64 f64 str str str str
"count" "1831560" "35614" "78406" "1831560" 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 "1805134" 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 "1831560" 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 1.83156e6 "1788168" "1831560" "1831560" "1831560"
"null_count" "0" "1795946" "1753154" "0" 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 "26426" 0.0 0.0 0.0 0.0 0.0 0.0 "0" 0.0 0.0 0.0 0.0 0.0 0.0 "43392" "0" "0" "0"
"mean" null null null null 0.008967 418.62922 6.499244 9.722575 3.456912 80.668031 5.18864 1.549416 0.398376 5.616697 null 0.54128 5.269918 3.389479 3.156906 2.541117 2.441586 null -0.995465 28.001124 6.187292 1.69613 0.443761 414.111544 null null null null
"std" null null null null 0.177118 186.287522 14.273599 49.829009 1.897507 42.587594 2.445579 1.480168 0.710891 3.680179 null 0.225859 5.342831 3.640044 2.115381 2.433513 1.247167 null 0.067517 9.239756 2.86568 1.700808 0.780102 184.202916 null null null null
"min" "CHEMBL1" "(+) NEOMENTHOL" "'CLOPRA-''YELLOW'''|AHR-3070-C… "Small molecule" 0.0 4.0 1.0 1.0 -14.26 0.0 0.0 0.0 0.0 0.0 "N" 0.0 -19.61 0.0 -20.95 -29.96 0.0 "MOL" -1.0 0.0 0.0 0.0 0.0 4.0026 "ACID" "Ag+" "B.CC(=O)OC1CN2CCC1CC2" "AAAAEENPAALFRN-UHFFFAOYSA-N"
"25%" null null null null 0.0 326.42 1.0 2.0 2.34 54.35 4.0 1.0 0.0 3.0 null 0.38 0.0 0.0 1.96 1.29 2.0 null -1.0 22.0 4.0 1.0 0.0 323.2209 null null null null
"50%" null null null null 0.0 393.49 3.0 4.0 3.45 75.27 5.0 1.0 0.0 5.0 null 0.56 3.89 2.13 3.21 2.72 2.0 null -1.0 27.0 6.0 1.0 0.0 389.2355 null null null null
"75%" null null null null 0.0 471.73 6.0 9.0 4.59 99.38 6.0 2.0 1.0 7.0 null 0.73 10.61 6.58 4.43 4.02 3.0 null -1.0 33.0 8.0 2.0 1.0 467.1879 null null null null
"max" "CHEMBL99999" "xerophilusin G" "zygosporamide" "Small molecule" 4.0 12546.32 1334.0 17911.0 22.57 568.39 32.0 25.0 4.0 67.0 "Y" 0.95 14.0 61.9 24.88 22.99 30.0 "MOL" 1.0 79.0 35.0 31.0 4.0 8214.3786 "ZWITTERION" "Zn+2" "n1onc2c1NC1Nc3nonc3NC1N2" "ZZZZVQYIUQRCGN-UHFFFAOYSA-N"


Looking at max phase and QED weighted

To explore some of the physicochemical and molecular properties in the dataframe, “Max Phase” is one of the first few that I want to have a look. Each ChEMBL compound will have a max phase number from 0 to 4 (for this particular dataset only, you’ll notice that this can start from -1 in other ChEMBL datasets), where 4 means the compound is approved (e.g. a prescription medicine).

While looking at max phases, I’m thinking the compounds with max phase 0 could be a testing set for prediction of their max phase outcomes, while the max phase 4 compounds could be the training set for building a machine learning (ML) model. There are of course some obvious and potential flaws of splitting a dataset like this, please treat this as an exercise for demonstration only.

Below I’m checking that I do have max phases of molecules spanning across the entire range of 0 to 4.

## Alternative code
# df_pa.group_by("Max Phase", maintain_order = True).agg(pl.len())
df_pa.group_by("Max Phase").len()
shape: (5, 2)
Max Phase len
i64 u32
3 773
0 1826345
1 520
4 2870
2 1052

One of the other parameters I’m interested in is “QED Weighted” (Bickerton et al. 2012). It’s a measure of druglikeness for small molecules based on the concept of desirability, which is based on a total of 8 different molecular properties. These molecular properties include molecular weight, ALogP, polar surface area, number of hydrogen bond acceptors, number of hydrogen bond donors, number of rotatable bonds, number of aromatic rings and structural alerts. It’s normally recorded as a number ranging from 0 to 1, where 0 is the least druglike and 1 being the most druglike.


Prepare data prior to running ML model

The goal of this ML model is to answer this question - which physicochemical features may be useful to predict whether a compound may enter into max phase 4?

A rough plan at this stage is:

  • to filter out max phase 4 and 0 compounds from the df_pa dataset

  • to use “Max Phase” as the target y variable for a logistic regression (LR) model because ultimately stakeholders are more likely to be interested in knowing which candidate compounds have the most likely chance to reach the final approved phase during a drug discovery project (please bear in mind that this ML model building is more like a demo only and may not reflect real-life scenarios…)

  • to use “Polar Surface Area”, “HBA”, “HBD”, “#RO5 Violations”, “QED Weighted”, “CX LogP”, “CX LogD” and “Heavy atoms” as the training features for building the ML model (these molecular features are selected randomly and ideally I should include more… the first three are newly added in this updated version)

  • to use max phase 0 compounds as the testing set since they are the ones not assigned with any max phase category yet

  • to use max phase 4 compounds as the training set since they’ve reached the approved phase

  • to build a confusion matrix in the end to see if the features selected may generate a somewhat reasonable model for predicting the outcomes of these small molecules

Downsides of this plan are (there are many…):

  • the dataset used here may include some inorganic compounds (e.g. molecules containing metals)

  • the dataset is likely not very-chemically-diverse

  • there are only a small number of max phase 4 compounds for model training (especially comparing with the “much larger” amount of max phase 0 compounds present)

  • LR may not be the best choice of ML model to start off with (I just wanted to find out how it works at the time and because I’m updating the post now so I’m trying to keep it the same…)

  • other exploratory data analysis should be done really… e.g. how chemically diverse is the dataset?

  • this is only an exercise so will not be taking into account many other important aspects of drug discovery e.g. what disorders or illnesses are we targeting? (is it an experimentally validated and actually a suitable target?) and what have been done so far in the literatures? etc.

Firstly, I’m filtering out all the max phase 0 compounds and also all the max phase 4 compounds then save them in separate dataframes.

# Selecting all max phase 0 molecules with their training features
df_0 = df_pa.filter((pl.col("Max Phase") == 0)).select([
  "ChEMBL ID", 
  "Max Phase",
  "Polar Surface Area",
  "HBA",
  "HBD",
  "#RO5 Violations", 
  "QED Weighted", 
  "CX LogP", 
  "CX LogD", 
  "Heavy Atoms"
  ])
df_0  #1,826,345 mols
shape: (1_826_345, 10)
ChEMBL ID Max Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms
str i64 f64 i64 i64 i64 f64 f64 f64 i64
"CHEMBL539070" 0 73.06 6 2 0 0.63 2.57 2.57 17
"CHEMBL3335528" 0 269.57 18 5 2 0.09 3.31 -0.14 60
"CHEMBL2419030" 0 85.13 6 1 0 0.66 3.66 3.66 24
"CHEMBL3827271" 0 319.06 10 11 2 0.07 -6.88 -8.95 50
"CHEMBL3465961" 0 50.5 4 1 0 0.87 2.13 -0.44 23
… … … … … … … … … …
"CHEMBL2017916" 0 77.0 6 1 0 0.8 2.17 2.1 22
"CHEMBL374652" 0 36.02 2 2 1 0.42 5.36 5.36 26
"CHEMBL1416264" 0 85.07 7 1 0 0.54 2.47 2.47 27
"CHEMBL213734" 0 101.7 5 2 0 0.5 2.36 1.95 21
"CHEMBL1531634" 0 29.1 2 1 0 0.67 4.04 4.04 19 
# Selecting all max phase 4 molecules with their training features
df_4 = df_pa.filter((pl.col("Max Phase") == 4)).select([
      "ChEMBL ID", 
      "Max Phase",
      "Polar Surface Area",
      "HBA",
      "HBD",
      "#RO5 Violations", 
      "QED Weighted", 
      "CX LogP", 
      "CX LogD", 
      "Heavy Atoms"
      ])
df_4  # 2,870 mols
shape: (2_870, 10)
ChEMBL ID Max Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms
str i64 f64 i64 i64 i64 f64 f64 f64 i64
"CHEMBL1200675" 4 12.47 2 0 1 0.31 6.27 4.89 29
"CHEMBL1200436" 4 46.53 3 1 0 0.7 2.95 2.95 22
"CHEMBL1096882" 4 186.07 10 5 0 0.31 -1.97 -5.12 24
"CHEMBL256997" 4 76.22 4 1 0 0.8 3.92 0.68 21
"CHEMBL2023898" 4 174.64 8 4 2 0.14 4.18 4.16 54
… … … … … … … … … …
"CHEMBL1619785" 4 128.03 8 2 0 0.49 2.09 1.86 34
"CHEMBL4297142" 4 0.0 0 0 0 0.0 0.0 0.0 0
"CHEMBL2364639" 4 74.02 6 1 0 0.68 3.65 2.3 30
"CHEMBL1232131" 4 94.83 4 3 0 0.44 1.2 -1.18 12
"CHEMBL3989949" 4 95.92 6 1 0 0.9 1.66 1.66 18 
df_4.describe()
shape: (9, 11)
statistic ChEMBL ID Max Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms
str str f64 f64 f64 f64 f64 f64 f64 f64 f64
"count" "2870" 2870.0 2870.0 2870.0 2870.0 2870.0 2870.0 2870.0 2870.0 2870.0
"null_count" "0" 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
"mean" null 4.0 74.951927 4.574216 1.780488 0.340418 0.54947 1.966551 0.841592 23.863066
"std" null 0.0 53.875964 3.06589 1.762048 0.703626 0.235063 2.61846 3.015392 11.064927
"min" "CHEMBL1000" 4.0 0.0 0.0 0.0 0.0 0.0 -10.11 -17.92 0.0
"25%" null 4.0 37.61 2.0 1.0 0.0 0.39 0.34 -0.63 17.0
"50%" null 4.0 65.78 4.0 1.0 0.0 0.58 2.29 1.12 23.0
"75%" null 4.0 100.9 6.0 2.0 0.0 0.74 3.73 2.75 30.0
"max" "CHEMBL99946" 4.0 359.42 20.0 15.0 4.0 0.94 11.62 11.62 69.0

I have to go back to the first post at this point to change the data restriction criteria for the parquet file as I’ve found out in this step above that if I restrict data via inorganic flags, then it’ll remove a lot of compounds with calculated physicochemical properties (because most of the ones in the data are preclinical compounds with inorganic flag of -1 or max phase 0 compounds; inorganic flag of 1 is inorganic ones and 0 should be the non-inorganic ones - I’ve updated its definition above too after re-checking the ChEMBL schema). For now, I’m restricting data via removing compounds with zero targets rather than via the inorganic flag, and this leaves us a set of max phase 4 compounds with their physicochemical properties calculated in ChEMBL (which is needed for training the model in the next post).


Re-sampling via under-sampling - dealing with imbalanced dataset

Because of the presence of a large number of max phase 0 compounds (1,826,345 molecules) and the relatively smaller number of max phase 4 compounds (2,870 molecules) in the original dataset, I’m going to randomly sample 2800 small molecules from the max phase 0 group to balance out the ratio of max phase 0 versus max phase 4 compounds. This will allow us to have a similar amount of data in each group to avoid a very imbalanced dataset.

Note

For an alternative and likely a better way to deal with imbalanced datasets, please see my later post that attempts to use a “generalized threshold shifting” method.

df_0 = df_0.sample(n = 2800, shuffle = True, seed = 0)
df_0
shape: (2_800, 10)
ChEMBL ID Max Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms
str i64 f64 i64 i64 i64 f64 f64 f64 i64
"CHEMBL1420527" 0 66.81 4 1 0 0.47 3.94 3.94 32
"CHEMBL1878707" 0 62.55 3 1 0 0.93 3.38 3.38 25
"CHEMBL3971884" 0 73.86 5 1 2 0.12 9.34 9.34 40
"CHEMBL3438846" 0 84.22 4 2 0 0.76 2.01 -0.19 26
"CHEMBL1819378" 0 40.46 4 0 0 0.62 4.0 4.0 26
… … … … … … … … … …
"CHEMBL2145354" 0 84.67 5 1 0 0.68 3.51 3.51 30
"CHEMBL3926224" 0 87.15 5 1 2 0.33 2.51 2.51 38
"CHEMBL3448203" 0 86.37 4 1 0 0.88 1.7 1.7 26
"CHEMBL1478078" 0 37.38 2 0 0 0.74 2.9 2.9 20
"CHEMBL3105607" 0 64.99 5 1 0 0.68 4.15 4.15 28

Since the plan is to use LR method for the ML model, the y variable I’m interested in is going to be a binary categorical variable - meaning it needs to be 0 (not approved) or 1 (approved). To do this, I’m going to add a new column with a new name of “Max_Phase” and replace “4” as “1” by dividing the whole column by 4 to reach this new label in the max phase 4 dataset.

Then I’m changing the data type of “Max_Phase” from float to integer, so that the two different dataframes can be concatenated (which will only work if both are of the same data types).

df_4 = df_4.with_columns((pl.col("Max Phase") / 4).alias("Max_Phase")).cast({"Max_Phase": pl.Int64})
df_4
shape: (2_870, 11)
ChEMBL ID Max Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms Max_Phase
str i64 f64 i64 i64 i64 f64 f64 f64 i64 i64
"CHEMBL1200675" 4 12.47 2 0 1 0.31 6.27 4.89 29 1
"CHEMBL1200436" 4 46.53 3 1 0 0.7 2.95 2.95 22 1
"CHEMBL1096882" 4 186.07 10 5 0 0.31 -1.97 -5.12 24 1
"CHEMBL256997" 4 76.22 4 1 0 0.8 3.92 0.68 21 1
"CHEMBL2023898" 4 174.64 8 4 2 0.14 4.18 4.16 54 1
… … … … … … … … … … …
"CHEMBL1619785" 4 128.03 8 2 0 0.49 2.09 1.86 34 1
"CHEMBL4297142" 4 0.0 0 0 0 0.0 0.0 0.0 0 1
"CHEMBL2364639" 4 74.02 6 1 0 0.68 3.65 2.3 30 1
"CHEMBL1232131" 4 94.83 4 3 0 0.44 1.2 -1.18 12 1
"CHEMBL3989949" 4 95.92 6 1 0 0.9 1.66 1.66 18 1

I’m also going to create a new column with the same name of “Max_Phase” in the max phase 0 dataframe, so that the two dataframes can be combined.

df_0 = df_0.with_columns((pl.col("Max Phase")).alias("Max_Phase"))
df_0
shape: (2_800, 11)
ChEMBL ID Max Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms Max_Phase
str i64 f64 i64 i64 i64 f64 f64 f64 i64 i64
"CHEMBL1420527" 0 66.81 4 1 0 0.47 3.94 3.94 32 0
"CHEMBL1878707" 0 62.55 3 1 0 0.93 3.38 3.38 25 0
"CHEMBL3971884" 0 73.86 5 1 2 0.12 9.34 9.34 40 0
"CHEMBL3438846" 0 84.22 4 2 0 0.76 2.01 -0.19 26 0
"CHEMBL1819378" 0 40.46 4 0 0 0.62 4.0 4.0 26 0
… … … … … … … … … … …
"CHEMBL2145354" 0 84.67 5 1 0 0.68 3.51 3.51 30 0
"CHEMBL3926224" 0 87.15 5 1 2 0.33 2.51 2.51 38 0
"CHEMBL3448203" 0 86.37 4 1 0 0.88 1.7 1.7 26 0
"CHEMBL1478078" 0 37.38 2 0 0 0.74 2.9 2.9 20 0
"CHEMBL3105607" 0 64.99 5 1 0 0.68 4.15 4.15 28 0

Then I’m combining df_0 (dataframe with max phase 0 compounds) and df_4 (dataframe with max phase 4 compounds).

df_concat = pl.concat([df_0, df_4])
df_concat
shape: (5_670, 11)
ChEMBL ID Max Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms Max_Phase
str i64 f64 i64 i64 i64 f64 f64 f64 i64 i64
"CHEMBL1420527" 0 66.81 4 1 0 0.47 3.94 3.94 32 0
"CHEMBL1878707" 0 62.55 3 1 0 0.93 3.38 3.38 25 0
"CHEMBL3971884" 0 73.86 5 1 2 0.12 9.34 9.34 40 0
"CHEMBL3438846" 0 84.22 4 2 0 0.76 2.01 -0.19 26 0
"CHEMBL1819378" 0 40.46 4 0 0 0.62 4.0 4.0 26 0
… … … … … … … … … … …
"CHEMBL1619785" 4 128.03 8 2 0 0.49 2.09 1.86 34 1
"CHEMBL4297142" 4 0.0 0 0 0 0.0 0.0 0.0 0 1
"CHEMBL2364639" 4 74.02 6 1 0 0.68 3.65 2.3 30 1
"CHEMBL1232131" 4 94.83 4 3 0 0.44 1.2 -1.18 12 1
"CHEMBL3989949" 4 95.92 6 1 0 0.9 1.66 1.66 18 1

This df_concat dataset is checked to see if it has all compounds in max phases 0 and 4 only.

Note: max phase 4 (approved) compounds are re-labelled as Max_Phase = 1.

df_concat.group_by("Max_Phase").len()
shape: (2, 2)
Max_Phase len
i64 u32
0 2800
1 2870

So here we have the final version of the dataset, renamed to df_ml to avoid confusion from the previous dataframes, before entering into the ML phase.

df_ml = df_concat.select([
  "Max_Phase", 
  "Polar Surface Area",
  "HBA",
  "HBD",
  "#RO5 Violations", 
  "QED Weighted", 
  "CX LogP", 
  "CX LogD", 
  "Heavy Atoms"
  ])
df_ml
shape: (5_670, 9)
Max_Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms
i64 f64 i64 i64 i64 f64 f64 f64 i64
0 66.81 4 1 0 0.47 3.94 3.94 32
0 62.55 3 1 0 0.93 3.38 3.38 25
0 73.86 5 1 2 0.12 9.34 9.34 40
0 84.22 4 2 0 0.76 2.01 -0.19 26
0 40.46 4 0 0 0.62 4.0 4.0 26
… … … … … … … … …
1 128.03 8 2 0 0.49 2.09 1.86 34
1 0.0 0 0 0 0.0 0.0 0.0 0
1 74.02 6 1 0 0.68 3.65 2.3 30
1 94.83 4 3 0 0.44 1.2 -1.18 12
1 95.92 6 1 0 0.9 1.66 1.66 18

I’m just checking for any nulls in the final dataset (should be none as I’ve filled all nulls in one of the steps above).

df_ml.null_count()
shape: (1, 9)
Max_Phase Polar Surface Area HBA HBD #RO5 Violations QED Weighted CX LogP CX LogD Heavy Atoms
u32 u32 u32 u32 u32 u32 u32 u32 u32
0 0 0 0 0 0 0 0 0

Also, I’m checking the data types in df_ml dataset to make sure they’re all integers or floats for scikit-learn algorithms to work.

# also shown in the dataframe below the column names
df_ml.dtypes
[Int64, Float64, Int64, Int64, Int64, Float64, Float64, Float64, Int64]

Then finally I’m saving this tidied dataset as a separate .csv file for use in the next post.

df_ml.write_csv("df_ml.csv", separator = ",")
Back to top

References

Bickerton, G. Richard, Gaia V. Paolini, Jérémy Besnard, Sorel Muresan, and Andrew L. Hopkins. 2012. “Quantifying the Chemical Beauty of Drugs.” Nature Chemistry 4 (2): 90–98. https://doi.org/10.1038/nchem.1243.