Identification and validation of novel microtubule suppressors with an imidazopyridine scaffold through structure-based virtual screening and docking.

Targeting the colchicine binding site of α/β tubulin microtubules can lead to suppression of microtubule dynamics, cell cycle arrest and apoptosis. Therefore, the development of microtubule (MT) inhibitors is considered a promising route to anticancer agents. Our approach to identify novel scaffolds as MT inhibitors depends on a 3D-structure-based pharmacophore approach and docking using three programs MOE, Autodock and BUDE (Bristol University Docking Engine) to screen a library of virtual compounds. From this work we identified the compound 7-(3-hydroxy-4-methoxy-phenyl)-3-(3-trifluoromethyl-phenyl)-6,7-dihydro-3H-imidazo[4,5-b]pyridin-5-ol (6) as a novel inhibitor scaffold. This compound inhibited several types of cancer cell proliferation at low micromolar concentrations with low toxicity. Compound 6 caused cell cycle arrest in the G2/M phase and blocked tubulin polymerization at low micromolar concentration (IC50 = 6.1 ±0.1 μM), inducing apoptosis via activation of caspase 9, increasing the level of the pro-apoptotic protein Bax and decreasing the level of the anti-apoptotic protein Bcl2. In summary, our approach identified a lead compound with potential antimitotic and antiproliferative activity. This journal is © The Royal Society of Chemistry.

Introduction

Microtubules (MTs) consist of α/β tubulin heterodimers[1] and are ubiquitous in all eukaryotic cells being the key component of the cytoskeleton.[2] They play a critical role in many cellular processes, including cell division in which they assemble to make up the mitotic spindle required for segregation of the chromosomes to the spindle poles and consequent daughter cells; cell proliferation, maintenance of cell shape and signal transduction; and MT-motor proteins that transport diverse cellular cargoes.[3] MTs are characterized by their highly dynamic behaviour, as they switch between periods of elongation and shortening.[4] Targeting the process of microtubule dynamics is an excellent strategy for chemotherapy and modulation of MT dynamics is considered to be one of the most successful approaches in the treatment of cancer.[5] Microtubule-targeting agents are classified into microtubule destabilizers and microtubule stabilizers according to the mechanism by which they affect microtubule dynamics.[6] There are four major ligand binding sites identified on the microtubule namely: the vinca and colchicine sites where ligand binding induces microtubule destabilization and the taxane and peloruside/laulimalide sites where binding typically induces microtubule stabilization.[7]

Taxanes, vinca alkaloids, and colchicine all showed potent inhibition of cancer cell lines. However, colchicine showed limitations as an antitumor agent in clinical trials due to its narrow therapeutic window.[8] While vinca alkaloids and taxanes are effective, they are also complex natural products that are difficult to synthesize and generally show poor bioavailability.[9,10] In addition, the emergence of resistance to these drugs has been reported.[11,12] Research has focussed on developing novel colchicine site inhibitors (CSI), since the molecular structure of known colchicine-site inhibitors is less complex than that of taxanes and vinca alkaloids.[13] However, microtubule-destabilizing agents that bind at the colchicine-binding site and reach clinical trials had significant side effects, for example ZD6126 is a phosphate prodrug of N-acetylcolchinol (NAC) releasing the drug after administration and in vivo hydrolysis. NAC binds to the colchicine binding site and inhibits tubulin polymerization, reducing the proliferating immature endothelial cells that line the tumour blood vessels and consequently inducing tumour cell death. ZD6126 reached phase II trials for metastatic renal cell carcinoma and induced necrosis in the tumours causing a large reduction in tumour cell yield after a single dose of ZD6126 but due to its cardiotoxicity, was withdrawn.[14] Similarly, ABT-751 showed antitumor activities against a broad spectrum of cancers including those resistant to conventional chemotherapies. ABT-751 is an anti-tubulin agent with anti-vascular properties that is responsible for the dysfunction of tumour blood vessels. Despite administration of a single dose of ABT-751 (30 mg kg−1, intravenously) disrupting tumour neovascularisation, it was withdrawn from Phase II due to adverse side effects.[15]

Despite the great potential of combretastatin and its prodrugs, CA1P and CA4P, these also suffer from drawbacks.[16] As such, there remains an urgent need to design and discover novel MT inhibitors based on the colchicine site of β-tubulin.[17] Therefore, our objective is to identify a novel chemical scaffold to bind the colchicine site that may form the basis of a new lead compound which offers promising antimitotic and antiproliferative activity that is also well-tolerated. Here we have used structure-based pharmacophore virtual screening of a subset of the ZINC15 database.[18] This computational approach represents a quick and efficient method in the identification of novel and diverse scaffolds of CSI. Molecular docking was carried out with three programs:[19] the Bristol University Docking Engine (BUDE),[20,21] AutoDock 4.2[22] and MOE (https://www.chemcomp.com/). Selected hits from this procedure were assessed for antiproliferative activity via their ability to cause mitotic spindle arrest, affect tubulin polymerization and induce apoptosis.

Results and discussion

Identifying colchicine binding site inhibitors

The binding orientation of colchicine was unequivocally established in the structural determination of the tubulin/DAMA-colchicine complex due to the electron density of the sulphur atom (PDB: 1SA0).[23] Hence, a 3D pharmacophore was built based on the analysis of the interaction of colchicine with tubulin in this structure. The derived pharmacophore (Fig. S1†) consists of seven features[24] comprising three H-bond acceptor centres corresponding to the interactions with Cysβ241 and Valα181 (by convention, α and β correspond to subunits A and B PDB labels), one H-bond donor, an aromatic centre, and two hydrophobic centres corresponding to interactions with Leuβ248, Alaβ250, Leuβ255, Asnβ258, Alaβ316, and Valβl318 (zone 2).[25] In order to explore this approach we chose the first one hundred thousand compounds from the current 9.9 million clean, drug-like and purchasable compounds in the ZINC15 database and filtered this set against the 3D pharmacophore using MOE. This process afforded 2476 compounds matching at least 4 points of the pharmacophore. This set of compounds and the native ligand (colchicine) were docked with BUDE, MOE, and AutoDock 4.2 into the colchicine binding site. MOE and Autodock4.2 are used widely in computational studies, BUDE is an in-house docking program using an empirical free energy forcefield to predict ligand affinities and has been used for inhibitor discovery.[26-29] The predicted binding affinities of colchicine using BUDE, MOE, and AutoDock4.2 (−100.27 kJ mol−1, −5.1 kcal mol−1 and −9.52 kcal mol−1) respectively, were used as a cut-off for selecting potential hits. The numbers of compounds passing this filter with a binding score better than colchicine were 188, 107, and 226 respectively. Next, we applied the criterion that compounds must be common to two or three of these docked sets. A total of 99 compounds passed this selection step and were assessed for toxicity risk using Osiris Property Explorer (https://openmolecules.org), resulting in the exclusion of a further 38 compounds (ESI† Properties_99-compounds.xlsx). The remaining 61 compounds were clustered using the Flexophore descriptor implemented in the DataWarrior cheminformatics application (Fig. S2A†) (https://openmolecules.org). A shortlist was generated by sampling from the compound clusters and the final shortlist of 13 compounds was made using actual compound availability and using cost as a proxy for synthetic accessibility. Of the thirteen shortlisted, only compound 3 is reported in the scientific literature[30] according to a search using the Chemical Abstracts service, SciFinder. The similarity of the compound set with the 50 CSI ligands currently present as tubulin complexes[31] in the PDB was assessed using DataWarrior and the results (Fig. S2B†) with the Flexophore descriptor[32] show partial overlap at the 60% level. Likewise, the fundamental geometric and chemical properties of the 50-CSI and our 99-compound sets show extensive overlap, although our compounds from virtual screening are generally somewhat larger (Fig. S3†). The 13 shortlisted compounds all passed the PAINS filters (https://zinc15.docking.org/patterns/home and https://www.cbligand.org/PAINS) and were purchased from MCULE (https://mcule.com) for experimental evaluation (Fig. 1, Tables S1 and S2†).

Fig. 1

Structures of the 13 selected compounds.

Structural characteristics of the shortlisted compounds

In all but one case, the selected compounds are characterised by having core structures comprising two or three fused cyclic rings with either some or extensive aromatic character (Fig. 1). The exception is compound 3 which consists of two substituted phenyl rings coupled by a five-atom long linker. The core structures are fused bicyclic rings: [5,6,0] (4 (pyrrolo[3,2-b]pyridine), 6 (imidazo[4,5-b]pyridine), 10 (pyrazolo[3,4-d]pyrimidine)), [6,6,0] (1 (hexahydroquinoline), 2, 5, 8 and 9 (dihydroquinoline-4-one), 11 (dihydroquinazolin-4-one); [6,7,0] (12 (benzodioxepin)) and fused tricyclic structures with ring sizes 6:6:6 (7 (pyrano[3,2-c]quinolin-5-one)) and 5:6:6 (13 (pyrazolo[4,5-b]quinoline-5-one)) and all are N or O heterocyclic systems including quinolines and quinolones. The core structures are typically elaborated with two phenyl groups, either attached directly or via short linkers to give compounds with moderate flexibility. Four compounds (10, 11, 12, 13) have heterocyclic rings attached to their cores. All but two compounds (10, 12) have at least one methoxy-aryl substituent, four (5, 6, 8, 9) are fluorinated and one (3) has a dichlorophenyl group. Full compound names are included in the ESI† File S1.xlsx.

Biological evaluation

Antiproliferative activity

The anti-proliferation activities of the shortlisted 13 compounds were evaluated against MCF-7, MDA-231, and A549 cancer cells lines using the MTT assay[33] at 10 μM in comparison to paclitaxel (16 nM), employed as a positive control. The results are shown in Fig. 2A and indicate that 6, 8, 9, 13 showed significant anti-proliferative activities against all three cell lines. Compound 6 caused 50, 40, and 60% inhibition of MCF-7, A549, and MDA-231 cell growth, respectively. Similarly, 8 caused 50, 20, and 80% inhibition of MCF-7, A549, and MDA-231 cell growth, respectively. The best activity observed was against MDA-231, a triple-negative cell line that is highly aggressive and resistant to treatment. Compounds 9 and 13 caused 50 and 40% inhibition of all tested cancer cell types respectively, with limited toxicity to normal fibroblast F180 cells (Fig. 2B).

Fig. 2

A) Antiproliferative activity of compounds 1–13 identified from the virtual screening against MCF-7, MDA-231 and A549 cell lines. B) Cytotoxic screening of compounds 6, 8, 9 and 13 against Fibroblast F180 mammalian cell line showed that compound 6 with highest% viability which indicates its low toxicity. Control: Fibroblast F180 mammalian cell line.

The most active four compounds 6, 8, 9 and 13 were titrated against the three cell lines to determine their IC50 values in this cell assay (Fig. 3). Compounds 6, 8, and 9 showed IC50 values in the 9–20 μM range against MCF-7, A549, and MDA-231 while 13 showed IC50 values above 20 μM (Table 1) hence 6, 8 and 9 were chosen for more detailed investigation.

Fig. 3

IC50 of the active compounds 6, 8, 9 and 13 against MCF-7, MDA-231 and A549 cancer cell lines, paclitaxel was used as positive control. Each result is a mean of triplicate experiments, and the mean values and standard error of mean are shown.

IC50 values of antiproliferative activity of compounds against MCF-7, A549, and MDA-231 cancer cell lines

Compound	MCF7	A549	MDA-231
6	9.4 ± 0.1 (μM)	16.1 ± 0.3 (μM)	11.9 ± 0.4 (μM)
8	12.0 ± 0.5 (μM)	13.2 ± 1.1 (μM)	14.1 ± 1.0 (μM)
9	11.0 ± 0.1 (μM)	9.4 ± 0.6 (μM)	19.5 ± 0.7 (μM)
13	23.2 ± 2.4 (μM)	>20 ± 4 (μM)	>23 ± 3 (μM)
Paclitaxel	6.1 ± 1.0 (nM)	15.0 ± 5.0 (nM)	19.3 ± 1.0 (nM)

Inhibition of mitotic spindle formation

An immunofluorescence assay was used to investigate the mechanism of action of 6, 8, and 9 on tubulin organization into mitotic spindles during cell division.[34] All three compounds caused the formation of classical multipolar spindle profiles (Fig. 4), similar to the positive controls paclitaxel and colchicine. These results indicated that compounds 6, 8, and 9 caused modulation of tubulin assembly with irregular morphology, showing typical mitotic arrest.

Fig. 4

The three compounds 6, 8 and 9 disrupt microtubule formation in A549 cells. A549 cells were treated for 24 h with DMSO as control, paclitaxel (14 nM), colchicine (0.1 nM), compound 6 (16 μM), compound 8 (13 μM), compound 9 (9 μM), then fixed, and stained with anti-β-tubulin antibody (green) and with DAPI for DNA (blue) to visualize the microtubules.

Tubulin polymerization assay in vitro

Microtubule polymer solutions scatter light in a concentration-dependent manner.[35,36] This behaviour was used to monitor the effect of ligands on microtubule polymerization (Fig. 5). Compounds 6, 8, and 9 at 15 μM and paclitaxel and colchicine at 3 μM were incubated with unpolymerized tubulin protein at 37 °C. The tubulin polymerization activities were determined by measuring the fluorescence and recording the area under the curve (AUC). Increasing fluorescence indicates increasing polymerization activity while decreasing the fluorescence indicates greater depolymerization activity. Paclitaxel which stabilizes polymerized tubulin caused an increase in the AUC by ∼1.3 fold. On the other hand, the destabilizing compound, colchicine, caused a decrease the AUC by ∼1.3 fold when compared to the negative control, hence inhibiting polymerization of tubulin. Compound 6 showed a modest decrease in fluorescence in comparison to the control while 9 increased the AUC by a similar amount and 8 showed an AUC similar to the control. Although these results are within error, they suggest that 6 inhibits and 9 stabilises tubulin polymerisation (Fig. 5A and B).

Fig. 5

A) Tubulin polymerization of compounds 6, 8 and 9. Tubulin polymerization was monitored by the increase in the fluorescence at 360 nm (excitation) and 420 nm (emission) for 1 h at 37 °C. Paclitaxel and colchicine were used as the positive control while 0.1% DMSO used as negative control. B) Area under the curve for the tested compounds and positive and negative controls.

Tubulin polymerization inhibition mechanism in cells

An ELISA assay was used to measure tubulin polymerization in MCF7 cells in the presence of compounds 6, 8, and 9. The results are shown in Fig. 6 and indicate that 6 behaves like colchicine as a suppressor of microtubule polymerization, while 8 and 9 enhance polymerization. These results are in reasonable accordance with the in vitro tubulin polymerization assay results and are more reliable, showing differences well outside the standard error. IC50 values of compounds 6, 8, and 9 required to modulate tubulin polymerization were 6.1 ± 0.1, 13.1 ± 0.3, and 12.8 ± 0.2 μM respectively (Fig. 6, Table 2). Compound 6 was selected for further investigation due to its low toxicity on normal cells while affording good potency against cancer cells (particularly MCF7 cells). It also showed the best activity in the tubulin polymerization assays. Although the IC50 of compound 6 is a little higher than colchicine (6.1 ± 0.1 μM and IC50 1.4 ± 0.02 μM respectively, Table 2) the low cytotoxicity of this compound suggests that it is a promising lead-molecule and scaffold for drug development.

Fig. 6

ELISA binding assay. A) MCF7 cell line treated with compounds 6, 8 and 9 at their IC50 for 24 h, paclitaxel and colchicine were used as positive controls. The tubulin polymer was extracted, and the quantities of monomeric and polymeric tubulin were measured using ELISA. Compound 6 showed the same effect as colchicine. B) IC50 of the active compounds 6, 8, 9 and colchicine in tubulin polymerization assay. Each result is a mean of triplicate experiments, and the mean values and standard error of mean are shown.

IC50 value of compounds affecting tubulin polymerization in MCF-7 cells

Compounds	IC50 (μM)a
6	6.1 ± 0.1
8	13.1 ± 0.3
9	12.8 ± 0.2
Colchicine	1.4 ± 0.02

IC50 values are the mean of three replicate experiments ± SD.

Compound 6 inhibited cell cycle progression at G2/M and induced apoptosis

Cell cycle analysis was performed to determine at which phase compound 6 exerted its antimitotic effect. MCF7 cells were treated with compound 6 at its IC50 concentration (6 μM) and 0.1% DMSO as control and incubated for 24 h, followed by measuring cell cycle distribution by flow cytometry. The results showed that the antimitotic activity of compound 6 was through apoptosis. Compound 6 caused a 10-fold increase in cell populations at the G2/M phase compared to control (Fig. 7A–C, Table S3†). The apoptotic activity of compound 6 was further evaluated by propidium iodide (PI) and annexin-V-FITC labeling assay on MCF-7 cells and using flow cytometry analysis.[37] Compound 6 caused an 84-fold increase in the late stage of cellular apoptosis compared to the negative control (Fig. 7D–F). Furthermore, the compound also caused a decrease in early apoptosis by 11-fold compared to the negative control (Fig. 7D–F).

Fig. 7

Cell cycle distribution. MCF-7 cells were treated with A) 0.1% DMSO (B), compound 6 (6 μM) for 24 h. Next, the cells were harvested and stained with propidium iodide, and flow cytometry cell cycle analysis was used to evaluate the cell cycle progression. C) Percentages of cells in the different phases indicate compound 6 arrests the cell cycle at (G2/M). Cell apoptosis D) MCF7 treated with compound 6 for 24 h showed late apoptosis. E) MCF7 treated with 0.1% DMSO as control. F) Percentage of MCF7 cell apoptosis showed that compound 6 increased both stages of apoptosis with a large increase in late apoptosis.

Compound 6 acts by a dual apoptosis mechanism

The apoptotic effect of compound 6 was determined by measuring caspase 9 levels. Caspase 9 initiates apoptosis by activating a cascade of intracellular events.[38,39] MCF7 cancer cells were treated with compound 6 at 6 μM for 24 h. The level of caspase 9 was then measured by Enzyme-Linked Immuno-Sorbent assay (ELISA) analysis. Compound 6 showed an increase in the level of caspase-9 in treated cells by 7.6-fold compared to control-vehicle-treated cells (Table 3). The Bcl2 family of proteins plays a crucial role in regulating mitochondrial apoptosis by either increasing the level of Bcl2 anti-apoptotic proteins or downregulating the level of Bax, a pro-apoptotic protein. Hence, cancer cells may develop resistance to apoptosis by changing the level of the Bcl2 and Bax protein expression. The effect of compound 6 on the balance of Bcl2/Bax proteins in MCF7 cancer cells was investigated. The results showed that compound 6 increased the level of the pro-apoptotic protein Bax by 7.2-fold compared to the control, and at the same time decreased the level of the anti-apoptotic protein Bcl2 by 2.7-fold (Table 3).

Determination of caspase 9, Bcl2, and Bax levels in the MCF7 cells treated with compound 6

Compounds	Caspase9 (ng ml−1)	Bcl2 (ng ml−1)	Bax (pg ml−1)
6	13.7 ± 0.7	2.1 ± 0.1	207.6 ± 6.3
Control	1.8 ± 0.3	5.7 ± 0.1	28.8 ± 1.6

Compound verification

1H NMR spectra of compounds 6, 8 and 9 (Fig. S6–S11†) were measured to confirm the chemical structures. Likewise, molecular weights were confirmed by mass spectroscopy for compound 6 MS m/z = 404.36 [M+ + 1] and compound 8 MS m/z 497.37 [M+ + 1] (Fig. S2 and S4† respectively).

Computational modelling

Molecular docking

To further elucidate the interactions between compounds 6, 8 and 9 and tubulin, the binding mode of the three compounds within the colchicine binding site was investigated. The results of docking 6, 8, and 9 are summarized in Table S4.† The docking study revealed that the compounds fitted well within the hydrophobic pocket of β-tubulin, making hydrophobic interactions with the hydrophobic residues lining zone 2 of the colchicine binding pocket, namely Leuβ242, Leuβ248, Alaβ250, Leuβ255, Valβl315, Alaβ316, and Ileβ378. The trifluoro phenyl ring of compound 6 occupied the same position as the trimethoxy phenyl ring of colchicine, facilitating hydrophobic interactions with key residues (Fig. 8A). In the case of compound 8, the 4-fluoro phenyl ring occupied the position of the trimethoxy phenyl ring of colchicine and the dimethoxy quinoline ring adopted an analogous position as the tropone ring of colchicine (Fig. 8B). However, compound 9 showed a different binding mode with the dimethoxy quinoline ring positioned where the trimethoxy phenyl ring of colchicine and the 4-fluoro phenyl group of tropone (Fig. 8C) reside. The three compounds showed strong H-bonds, as the OH of compound 6 served as H donor and formed an H-bond with Thrα179 and OCH3 formed an H-bond with NH-Alaβ250 (Fig. 9A, and Table S4†). Both compounds 8 and 9 formed H-bonds with Serα178, and Cysβ241. (Fig. 9B and C, and Table S4†). The trifluoro phenyl group of compound 6 faced SH-Cysβ241, allowing the lone pair of electrons of the sulphur to interact with the pi cloud of the aromatic ring (SH---π).[40] Analyses of the docking results indicate that the hydrophobic interactions of compounds 6, 8 and 9 predominate while H-bonds help in the proper orientation of the compounds within the binding pocket.

Fig. 8

Comparison of the proposed binding modes of compounds 6, 8 and 9 into α/β interface of tubulin (PDB: 1SA0) with DAMA–colchicine binding mode (dark green sticks). A) Compound 6 (purple sticks), B) compound 8 (purple sticks) C) compound 9 (purple sticks).

Fig. 9

The interactions of compounds 6, 8, and 9 with tubulin dimer (PDB: 1SA0). A) Compound 6, B) compound 8, C) compound 9, hydrogen bonds are shown as dotted lines.

Molecular dynamics

Simulations of the tubulin–ligand complexes (colchicine, 6, 8, 9, 14 and 15) were performed for 3 repeats of 500 ns with different initial starting velocities to assess their stability and persistence over this short period. The Root Mean Square Deviation (RMSD) of the protein and the ligand with respect to their initial (time 0 ns) positions are shown in Fig. S4† and average values in Table S5.† In all cases the RMSD of the protein stays between 2.0 and 3.0 Å for 99.15% of the 9 μs of total simulation time (maximum 3.4 Å). Colchicine and 6 also stay close to their original positions in all simulations with average RMD values of 1.5 Å (maximum 3.2 Å) and 1.4 Å (maximum 2.8 Å) respectively, indicating that the conformation and poses of these complexes are close to their original docked positions. The behaviour of 8 is more labile with the ligand RMSD values typically ranging between 3.0 and 4.0 Å. In these simulations, the 4-fluorophenyl group reorients closer to the deep-binding site (by Y224 and V238) the naphthyl group tilts by some 45° in two simulations while remaining close to the original pose in run 3. The greatest variability in pose and conformation occurs in the 3,4-difluorophenyl group at the dimer interface. Compound 9 stays close to its original docked position in runs 1 and 3 while in run 2 it rapidly slips the 4-methoxyphenyl ring into the deep-binding pocket causing the whole ligand to shift by some 3 Å into the B subunit, accounting for the large but stable rise in RMSD to 4–5 Å in this case. Compounds 14 and 15 are versions of 6 redesigned to occupy the deep-binding pocket (see below) and their plots of RMSD versus time indicate little movement from the original pose as seen for colchicine and 6. Ligand-residue contact data are provided for the 500 ns structures from the 12 simulations of tubulin with 6, 8, 9 and colchicine in Fig. S5.†

Rational design and synthetic accessibility

Many of the protein ligand interactions between 6 and tubulin are mediated by the two aryl groups attached to the imidazopyridine core. A simple synthetic route to varying these substituents has been described[41] and is shown here.

The mild conditions of this route and wide variety of substituted aniline and benzaldehyde starting materials available commercially make this a very attractive goal for medicinal chemistry to explore SAR and rational design for improving affinity and the drug-likeness of this potential lead. As an exemplar, we describe the design and modelling of 14 with variations to exploit the deep-binding pocket (X1 = 3-hydroxy; X2 = 4-pyrimidine) and improve interactions at the dimer interface (Y = 3,5-bistrifluoromethylbenzene).

We also investigated 15, equivalent to 14 with X1 = 3-hydroxy removed, since this is cheaper to synthesise and in conjunction with 14 would allow exploration of any contribution of the binding affinity of this group that hydrogen bonds with the carbonyl oxygen of V238 in our modelling.

Conclusion

In this paper we have used computational methods to identify inhibitors and modulators of tubulin polymerization. Firstly, 7-feature pharmacophore matching was used to filter out some 98% of compounds from a library of 100 000 available compounds. Consensus docking of these 2746 matched compounds was performed using three different methods (MOE, BUDE and Autodock) and passed 99 compounds. Cheminformatics was used to remove compounds with toxicological risk and to cluster compounds by similarity. A total of 13 compounds were chosen from the clustering for experimental investigation. The antiproliferative activities of the 13 compounds were evaluated against three cancer cell lines (MCF-7, MDA-231, and A549) and four compounds (6, 8, 9 and 13) showed significant results. Compounds 6, 8 and 9 had IC50 values ≤20 μM and these three compounds disrupted spindles in the mitotic cells giving a phenotype similar to colchicine. Compounds 6, 8 and 9 modulated tubulin polymerisation in vitro and in three cancer cell lines with minimal toxicity. Molecular dynamics simulations into the microsecond regime support the predicted binding modes. Compound 6 showed the lowest IC50 (6.1 ± 0.1 μM) inhibiting tubulin polymerization in MCF-7 cells using ELISA. FACS cell cycle analysis showed that 6 arrests the cell cycle at the G2/M phase and induces late apoptosis via upregulating caspase-9 and Bax while downregulating Bcl2. The activity of 6 against cancer cells, its low cytotoxicity to normal fibroblasts and ease of synthesis of variants provides a start-point for medicinal chemistry development. The results of docking and simulation of 6 were used to suggest elaborations to exploit the deep-binding site in the colchicine pocket and the persistence of binding of these derivatives (14 and 15) was demonstrated by further molecular dynamics simulations. The proposed synthetic route should facilitate production of a large number of derivatives based on the imidazopyridine scaffold of 6 to explore SAR and improve the drug-likeness of this series.

Experimental section

Chemistry

1H NMR spectra were recorded on a Bruker spectrometer at 500 MHz. Chemical shifts are expressed in parts per million (ppm) relative to tetramethylsilane and coupling constants (J values) are represented in Hertz (Hz) and the signals are designated as follows: s, singlet; d, doublet; t, triplet; m, multiple. Mass spectroscopic data were obtained through electrospray ionization (ESI) mass spectrometry.

7-(3-Hydroxy-4-methoxy-phenyl)-3-(3-trifluoromethyl-phenyl)-6,7-dihydro-3H-imidazo[4,5-b] pyridin-5-ol (compound 6). 1H NMR (DMSO-d6) δ: 2.72 (m, 1H, H̲–CH), 3.12 (m, 1H, H–C–H̲), 3.75 (s, 3H, OCH3), 4.23 (m, 1H, CH), 6.25 (m, 2H, ArH), 6.80 (d, 1H, J = 8, ArH), 7.54 (d, 1H, J = 6.5, ArH), 7.82 (m, 4H, ArH, imidazole CH̲), 8.60 (s, 1H, OH), 10.2 (s, 1H, OH). 13C NMR (DMSO-d6): δ 35.59 ppm (CH of the pyridine ring), δ 38.44 ppm (CH2 of pyridine ring), δ 56.79 ppm (OCH3), 124.14(CF3), δ 135.44 (CH imidazole), 113.22–147.46 (aromatic carbons). MS analysis for C20H16F3N3O3 Calcd mass 403.11, found (m/z, ESI+) (M+ + 1): 404.36.

N-(3,4-Difluorophenyl)-2-[3-(4-fluorobenzoyl)-6,7-dimethoxy-4-oxoquinolin-1-yl]-acetamide (compound 8). 1H NMR (DMSO-d6) δ: 3.85 (s, 6H, 2 OCH3), 5.23 (s, 2H, CH2) 7.02 (s, 1H, A̲r̲H), 7.32 (m, 3H, ArH), 7.45 (q, 1H, ArH), 7.62 (s, 1H, ArH), 7.83 (m, 3H, Ar–H), 8.33 (s, 1H, quinolone 2-CH), 10.80 (s, 1H, NH of acetamide). 13C NMR (DMSO-d6): δ 56.79 ppm: 2 singlet signals for 2 OCH3, δ 35.19 ppm for C̲H2–CO, δ 104.96–165.72 ppm for aromatic carbons, δ 169.84, 174.96, 188.76 ppm for 3 CO. MS analysis for C26H19F3N2O5 calcd mass 496.12, found (m/z, ESI+) (M+ + 1): 497.37.

N-(4-Fluorophenyl)-2-[6,7-dimethoxy-3-(4-methoxybenzoyl)-4-oxoquinolin-1-yl]-acetamide (compound 9). 1H NMR (DMSO-d6) δ: 3.85 (s, 9H, 3 OCH3), 5.24 (s, 2H, CH2) 7.02 (m, 3H, A̲r̲H), 7.2 (t, 2H, ArH), 7.60 (m, 3H, ArH), 7.74 (d, 2H, ArH), 8.3 (s, 1H, quinolone 2-CH), 10.62 (s, 1H, NH of acetamide). 13C NMR (DMSO-d6): δ 53.19–56.79 ppm (3 singlet signals for 3 OCH3), δ 56.79 ppm (CH2–CO), δ 104.96–163.49 ppm (aromatic carbons), δ 169.84, 174.56, 188.76 ppm for 3 CO. MS analysis for C27H23FN2O6 calcd mass 490.48, found (m/z, ESI+) (M+ + 1): 491.48.

Computational methods

Protein and database preparation: the X-ray crystal structure of α/β tubulin in complex with colchicine-DAMA (PDB code: 1SA0) was downloaded from protein data bank https://www.rcsb.org/ and used for virtual screening. Chains C and D were removed, the protein structure was prepared by inserting the missing loop regions using MODELLER[42]via the UCSF Chimera graphical interface.[43] Hydrogen atoms were added, and water molecules were removed using MOE. A set comprising the first 100 000 compounds was selected from the clean, druglike subset (9.9 M compounds) of the ZINC15 database to use for virtual screening. The compounds were saved as mdb format files by MOE. We constructed a 3D structure-based pharmacophore from the tubulin-colchicine complex (1SA0) using the protein–ligand interaction fingerprints (PLIF) application implemented in MOE. The pharmacophore model comprises seven features: three H-bond acceptors and one H-bond donor (F1, F2, F7 and F5) respectively; two hydrophobic centres (F4 and F6), and an aromatic centre (F3). These represent i) the two acceptors F1 and F7 corresponding to interaction with Cysβ241, ii) the third acceptor (F2) corresponding to the interaction of Val α181, iii) aromatic centre (F3), iv) two hydrophobic centres (F4 and F6) corresponding to hydrophobic interaction with Leuβ248, Alaβ250, Leuβ255, Asnβ258, Alaβ316, and Valβl318. The pharmacophore model was employed as a search query using MOE to identify commercial compounds targeting the colchicine binding site, matching at least 4 of the 7 pharmacophore features. This process afforded 2476 compounds for virtual screening by docking described below.

Molecular docking was performed using three programs: BUDE, AutoDock and MOE. Validation of the screening model and applied protocol was carried out by re-docking the native ligand (colchicine) in the binding site of the tubulin protein (1SA0). The RMSD value was then calculated with respect to the co-crystallized ligands. An RMSD value ≤1.0 Å between the X-ray structure and the best-scored conformations of the native ligand, the docking process was considered successful.[44] The RMSD between the re-docked and co-crystal ligand is less than 1 Å for the three programs, indicating a consensus of the methods and consistency of pose prediction. The predicted binding energies for colchicine re-docked by BUDE, AutoDock and MOE were −100.27 kJ mol−1, −9.52 kcal mol−1 and −5.1 kcal mol−1 respectively.

Virtual screening with BUDE

Molecular docking was performed using the Bristol University Docking Engine (BUDE) for the compounds filtered by the 3D pharmacophore into the tubulin–colchicine binding site. The BUDE search area was defined as grid centred on the native ligand (X = 42.848, Y = 52.376, Z = −8.531). BUDE is a rapid rigid-body docking program, hence ligand flexibility is achieved by docking multiple conformers of each ligand. All compounds with a predicted binding energy better than colchicine were selected (177 compounds with binding energy ≤−100 kJ mol−1).

Virtual screening with AutoDock 4.2

The native ligand was removed from the crystal structure (PDB: 1SA0) and AutoDock.4.2 used to convert both protein structure and the native ligand separately into PDBQT format. Polar hydrogen atoms and Kollman charges were assigned to the protein. Gasteiger partial charges were assigned to the ligand and non-polar hydrogen atoms merged with their heavy atoms Rotatable bonds in the ligand were defined using an AutoDock utility, AutoTors. The grid box was placed at the centroid of native ligand (X = 42.845, Y = 52.376 and Z = −8.531), the box size was 100 × 100 × 100 Å with a 0.2 A grid spacing and the grid map was calculated using Autogrid tool and saved as a gpf file. Docking was performed using the Lamarckian genetic algorithm, each docking experiment was performed 100 runs, the configuration file was saved as dpf format. Raccoon was used to prepare all the 2476 ligands to perform a docking with Autodock 4.2, the Raccoon software splits the multi-structure files of the ligands to separate PDBQT input files and generate configuration files and scripts for both and Autodock. The results were sorted according to the lowest predicted binding energy. 226 compounds had predicted binding energies ≤−9.52 kcal mol−1 (the colchicine binding energy calculated with AutoDock).

Virtual screening with MOE

The same protein (1SA0) and set of ligands from the 3D pharmacophore filtration and converted to mdb format. MOE was used to add hydrogen atoms to the ligands and energy minimised until the gradient of energy with respect to coordinates fell below 0.05 kcal mol−1 Å−1 under the was MMFF94X force field. The binding site was defined as the colchicine site. Ligands were docked using the Triangle Matcher method with the London ΔG scoring function. Refinement was performed using the rescoring affinity ΔG method. The lowest energy pose was chosen for each docked compound yielding 188 compounds with a binding energy better than colchicine (−5.1 kcal mol−1).

Selection of compounds for testing

Next, the requirement for a compound to be present in at least two docking search results was applied, giving 99 compounds (ESI† File S1.xlsx). These the 99 compounds were re-docked with MOE to allow a consistent set for visualization with Pymol.[39] Further selection by inspection was performed as described in Results. This process gave a shortlist of 13 compounds (Tables S1 and S2†). The shortlisted hits were screened for pan assay interference compounds (PAINS) using the online PAINS filters at https://zinc15.docking.org/patterns/home/ and https://www.cbligand.org/PAINS/. The 13 compounds passed both filters and were purchased from MCULE (USA) for experimental testing.

Molecular dynamics (MD) simulations

The recent structure of a tubulin tetramer with a stathmin-4 domain was used as the basis for the MD simulations (6F7C; 2.0 Å) firstly, subunits C and D were removed and stathmin truncated at N91. Crystallographic waters, cofactors GTP and GDP and metal ions associated with subunits A and B were retained apart from waters 653A and 691B at the colchicine site. Ligand coordinates of the six docked complexes of tubulin with colchicine, 6, 8, 9, 14 and 15 were transferred to this model and simulations performed for 3 repeats of 500 ns each using GROMACS 2019.4 and 2021.2[46] as follows. Pdb2gmx was used to add hydrogen atoms to the protein consistent with pH 7 and generate a topology file under the Amber99-SB-ildnb force field.[47] Acpype[48] was used to generate topology files of the compounds 6, 8 and 9 under the GAFF force field.[49] The ligand and protein complexes were centred in a triclinic box with a minimum margin of 1.5 nm and filled with TIP3P water. The system was neutralized by adding sodium and chloride ions to give an ionic strength of 0.15 M. The energy minimization (5000 steps) was conducted using steepest descents. All simulations were performed as NPT ensembles at 310 K under periodic boundary conditions. The Particle Mesh Ewald (PME) method was used for calculating long range electrostatics, and Van der Waals (VdW) interactions. The cut-off distance for the short-range VdW and Coulombic interactions was set to 1.2 nm.[50] Pressure was controlled by the Parrinello–Rahman barostat and temperature by the V-rescale thermostat. The simulations were integrated with a leap-frog algorithm over a 2 fs time step, constraining h-bond vibrations with the P-LINCS method. Molecular dynamics simulations were carried out for 500 ns on BlueCrystal, the University of Bristol's high-performance computing machine and the GW4 tier-2 machine Isambard. Simulation analyses were carried out using GROMACS tools, Xmgrace and gnuplot were used for plotting data, molecular graphics manipulations and visualizations were performed using Chimera v1.14,[43] VMD v1.9.4 (ref. 51) and OpenPymol v1.8.[45]

Biological methods

Cell culture and maintenance

Breast cancer cell lines including breast cancer cell lines (MCF-7), triple-negative breast cancer (MDA-231), adenocarcinomic human alveolar basal epithelial cells (A549), and normal fibroblast cells (F180) were cultured in Roswell Park Memorial Institute media (RPMI, Sigma-Aldrich, St. Louis, MO, USA) supplemented with 10% fetal bovine serum, and 1% penicillin/streptomycin. All cell lines were purchased from the European Collection of Cell Cultures (ECACC, UK). Cell lines were incubated at 37 °C in a humidified incubator containing 5% CO2.

Antiproliferative and cytotoxic activity

The antiproliferative activity of compounds 1–13 were evaluated by testing their effects on the aforementioned cell lines 3-(4,5-dimethyl thiazolyl-2)-2,5-diphenyltetrazolium bromide (MTT) assay as described before.[17] The cancer cell lines were seeded as 1 × 104 cells per well in 96-well flat bottom plates for 24 h. The cells were then treated with 10 μM of compounds 1–13 and incubated for 48 h at 37 °C in a humidified incubator containing 5% CO2. Paclitaxel and 0.1% DMSO were employed as positive and negative controls, respectively. The culture media was then removed, washed, and then incubated with 200 μl culture media containing 0.5 mg ml−1 MTT for 2 h. The blue formazan crystals, converted from the yellow MTT by viable cells, were dissolved by adding 200 μl DMSO and measured spectrophotometrically at 570 nm in a microplate reader (Thermo-Scientific, Vantaa, Finland). Cell viability was calculated, using the following formula: % of living cells = (OD experimental)/(OD control) × 100, while % of cell death was calculated by subtracting the living cells from the total number of cells. To determine the IC50 of the active compounds, the cells were treated with different concentrations 2.5, 5, 10, 20 μM of compounds (6, 8 and 9).

Immunofluorescence assay

A549 cells (5 × 104/well) were plated on coverslips in 6-well plates and treated with compounds 6, 8 and 9 at concentrations of 16, 13, and 9 μM, respectively for 24 h. Paclitaxel (14 nM) and Colchicine (0.1 μM) used as positive controls while 0.1% DMSO used as a negative control. The cells were then rinsed twice with PBS, fixed with 3.7% paraformaldehyde, and permeabilized with 0.1% Triton X-100. The cells were blocked with 1% BSA in PBS for 1 h prior to incubation with anti-β-tubulin mouse monoclonal antibody (#86298, Cell Signaling, San Francisco, CA, USA) overnight at 4 °C. The cells were washed with PBS for 1 h in the dark, and then incubated with Alexa Fluor® 488 secondary antibodies (Abcam). The cellular microtubules were observed with a fluorescence microscope (Olympus BX43, Japan).

In vitro tubulin polymerization assay

Tubulin polymerization was analyzed in vitro using a Tubulin Polymerization Assay kit (Cytoskeleton, Denver, CO). Briefly, 2 mg ml−1 porcine tubulin was dissolved in buffer 1 (80 mM PIPES, 2 mM MgCl2, 0.5 mM EGTA pH 6.9, 10 μM fluorescent reporter, 1 mM GTP, 15% glycerol). This solution was transferred to a pre-warmed 96-well plate and treated with 15 μM of the test compounds 6, 8 and 9. Colchicine and paclitaxel at 3 μM were employed as positive controls, while 0.1% DMSO treated-cells were employed as a negative control. Tubulin polymerization was monitored at 37 °C for 60 min using fluorescence microscopy. The reading speed was programmed at 1 cycle per min with excitation and emission wavelengths of 360 and 450 nm, respectively, using the Varioskan Flash spectral scanning multimode reader (Thermo Fisher Scientific). The 100% polymerization value was defined as the area under the curve (AUC) of the untreated control.

Extraction of soluble and polymerized tubulin fractions

Extraction was done using the protocol previously reported.[52] After treatment with drugs for 24 h, medium containing cells in suspension was recovered and pooled with adherent cells scraped in PBS pre-warmed at 37 °C. After centrifugation 5 min at 400 × g, wash with PBS, cells were extracted for 5 min with pre-warmed at 37 °C microtubule-stabilizing buffer (0.1 M PIPES pH 6.9, 14.5% glycerol, 0.5% Triton X-100, 20 mM EGTA and 5 mM MgCl2) containing Complete (Sigma FAST Protease Inhibitor Cocktail Tablet) and 10 ng ml−1 paclitaxel (Sigma). After centrifugation at 20 000 × g for 10 min at 25 °C, supernatants containing soluble fractions were transferred to a new tube, while polymerized fractions in pellets were recovered by incubation in RIPA buffer for 45 min on ice followed by centrifugation for 10 min at 20 000 × g. The total protein content in the samples was determined using DC protein assay kit (Bio-Rad, Hercules, CA, USA). Equivalent aliquots from polymeric fractions were measured using Human Beta-tubulin ELISA Kit (Abcam) following the manufacturer's instructions.

Determination of the IC50 of compound 6, 8 and 9

The IC50 of compounds 6, 8 and 9 required to inhibit tubulin polymerization was measured using ELISA kit (cat. # BK011P, Cytoskeleton, Denver, CO) with tubulin protein (cat. #T240-DX, Cytoskeleton, Denver, CO) according to the manufacturer's instructions.[53]

Cell cycle analysis

MCF7 cancer cells were treated with compounds 6 at 6 μM for 24 h. DMSO was employed as vehicle control. The cells were harvested, centrifuged, and the cell pellets were fixed with 70% ethanol on ice for 15 min. The fixed pellets were incubated with propidium iodide (Sigma-Aldrich, St. Louis, MO, USA), staining solution (50 mg mL−1 PI, 0.1 mg mL−1 RNaseA, 0.05% Triton X-100) for 1 h at room temperature. Cell cycle was assessed by Gallios flow cytometer (Beckman Coulter, Brea, CA, USA).

Apoptosis assay

The apoptosis assay was carried out using the FITC Annexin-V/PI commercial kit (Becton Dickenson, Franklin Lakes, NJ, USA), MCF-7 cells were treated with compound 6 at 6 μM for 24 h. 0.1% DMSO was used as negative control. Treated and control cells were stained using V/PI apoptosis kit. Samples were analysed by fluorescence-activated cell sorting (FACS) using flow cytometer run over one hour. Data were analysed using Kaluzav 1.2 (Beckman Coulter).[54]

Caspase-9 assay

MCF-7 cells were incubated without and with 6 μM of compound 6 for 24 h. The cells were washed in phosphate buffered saline and cell lysates were collected and level of caspase 9 was determined using ELISA kit (cat. # EIA-4860, Invitrogen, Carlsbad, CA, USA), and according to manufacturer instruction https://www.thermofisher.com/elisa/product/Caspase-9-Human-ELISA-Kit/BMS2025.

Determination of the effect compound 6 on BAX and Bcl-2 protein levels

MCF-7 cells, which were grown in RPMI1640 containing 10% fetal bovine serum. The cells were treated with compound 6 at 6 μM for 24 h. The cells were then lysed using cell extraction buffer. The collected lysate was diluted in standard diluent buffer and the levels of Bax and Bcl2 were measured as previously reported.[55]

Statistical analysis

Data were plotted using GraphPad Prism (5.04, La Jolla, CA, USA). Results are shown as the mean ± SEM of three independent replicates.

Conflicts of interest

There is no conflict of interest to declare.