Combined Molecular Docking, 3D-QSAR, and Pharmacophore Model: Design of Novel Tubulin Polymerization Inhibitors by Binding to Colchicine-binding Site

Article abstract of DOI:10.1111/cbdd.12545

Interference with dynamic equilibrium of microtubule-tubulin has proven to be a useful tactics in the clinic. Based on investigation into the structure-activity relationship (SAR) studies of tubulin polymerization inhibitors obtained from several worldwide groups, we attempted to design 691 compounds covering several main heterocyclic scaffolds as novel colchicine-site inhibitors (CSIs). Evaluated by a series of combination of commonly used computer methods such as molecular docking, 3D-QSAR, and pharmacophore model, we can obtain the ultimate 16 target compounds derived from five important basic scaffolds in the field of medicinal chemistry. Among these compounds, compound A-132 with in silico moderate activity was synthesized, and subsequently validated for preliminary inhibition of tubulin polymerization by immunofluorescence assay. In additional, the work of synthesis and validation of biological activity for other 15 various structure compounds will be completed in our laboratory. This study not only developed a hierarchical strategy to screen novel tubulin inhibitors effectively, but also widened the spectrum of chemical structures of canonical CSIs.

Full text of DOI:10.1111/cbdd.12545

Received Date : 11-Oct-2014
Revised Date : 18-Jan-2015
Accepted Date : 18-Feb-2015
Article type : Research Article
Combined Molecular Docking, 3D-QSAR, and
Pharmacophore Model: Design of Novel Tubulin
Polymerization Inhibitors by Binding to Colchicine-
binding Site
a
‡
*
b
‡
b
b
b*
Dong-Dong Li , Ya-Juan Qin , Xin Zhang , Yong Yin , Hai-Liang Zhu
a
College of Chemical Engineering, Nanjing Forestry University, Nanjing, 210073, P. R. China.
State Key Laboratory of Pharmaceutical Biotechnology, School of Life Sciences, Nanjing
b
University, Nanjing, 210093, P. R. China.
Corresponding Author
*
Dong-Dong Li, Tel.: +86 025 85427962; E-mail: lidon@njfu.edu.cn. Hai-Liang Zhu, Tel.:
86-25-8359 2572; Fax: +86 025 83592672; E-mail address: zhuhl@nju.edu.cn. Lin-Guo Zhao,
+
E-mail: lgzhao@njfu.edu.cn.
This article has been accepted for publication and undergone full peer review but has not been
through the copyediting, typesetting, pagination and proofreading process which may lead to
differences between this version and the Version of Record. Please cite this article as an
'Accepted Article', doi: 10.1111/cbdd.12545
This article is protected by copyright. All rights reserved.

Author Contributions
‡
These authors contributed equally.
ABSTRACT
Interference with dynamic equilibrium of microtubule-tubulin has proven to be a useful tactics
in the clinic. Based on investigation into the structure–activity relationship (SAR) studies of
tubulin polymerization Inhibitors obtained from several worldwide groups, we attempted to
design 691 compounds covering several main heterocyclic scaffolds as novel colchicine site
inhibitors (CSIs). Evaluated by a series of combination of commonly used computer methods
such as molecular docking, 3D-QSAR, and pharmacophore model, we can obtain the ultimate 16
target compounds derived from five important basic scaffolds in the field of medicinal chemistry.
Among these compounds, compound A-132 with in silico moderate activity was synthesized,
and subsequently validated for preliminary inhibition of tubulin polymerization by
immunofluorescence assay. In additional, the work of synthesis and validation of biological
activity for other 15 various structure compounds will be completed in our lab. This study not
only developed a hierarchical strategy to screen novel tubulin inhibitors effectively, but also
widened the spectrum of chemical structures of canonical CSIs.
KEYWORDS: Microtubules; Antitumor agents; Colchicine-binding site; 3D-QSAR;
Pharmacophore.
This article is protected by copyright. All rights reserved.

INTRODUCTION
Targeting the process of microtubule dynamics contributes much for cancer therapies, due to it plays an crucial
role in the maintenance of cell shape, signal transduction, and chromosome segregation during mitosis (1). So far,
microtubule-targeting agents could be simply classified into microtubule stabilizer and microtubule distabilizer
according to the mechanism by interfering with microtubule dynamics. Further divided, there are four major binding
sites provided by the microtubule: the taxane site and the laulimalide/ peloruside A site for microtubule stabilizer,
and the vinca site and the colchicine site for microtubule destabilizer (2).
The motivations for developing colchicine site inhibitors (CSIs) can be concluded as follows. Firstly, by
comparison of the two successful antitumor drugs taxanes and vinca alkaloids, the first CSI colchicine (Figure 1),
while inhibiting cancer cell proliferation potently, has limited its applications in clinical trials as antitumor agent due
to its narrow therapeutic index (3). Much possibility for searching novel scaffolds into the CSI chemical space will
be provided. Currently the CSI-like compound most close to the drug approved by FDA is CA-4P) (Zybrestat,
Figure 1), which is being developed under clinical trials as a treatment for solid tumors. Second , a large number of
small molecules binding to colchicine site have been reported, which possess significant structural diversity.
Meanwhile consideration in the development of CSIs as vascular disrupting agents (3) (VDA) that obviously differs
from cytotoxic agents such as taxanes and vinca alkaloids, This will stimulate the exploitation of CSIs in the further
step. Thirdly, molecular complexity of colchicine site inhibitors are less than that of taxol and vinca site binders,
which could be in more favor of exploring the CSI chemical space efficiently (5). In addition, unsolved problems
such as multiple drug resistance (MDR) and secondary toxicity induced by common cytotoxic agents in the clinical
use might be circumvented by means of colchicine site agents (6, 7). Thus, there is an urgent need to design and
synthesize novel tubulin inhibitors based on the colchicine site at the β-tubulin.
To date, there are dozens of CSIs (Figure 1) being evaluated under clinical investigation, and even more in
preclinical studies, mainly including combretastatin analogs, indole-, quinolone- and thiophene-based compounds,
chalcone compounds, and sulfonanilide compounds (3). Such structural diversity of the CSIs seems to be closely
involved with the inherent flexibility of the binding site, which can be demonstrated by several tubulin protein
crystal structures (PDB ID: 1SA0 (8), 1SA1 (9), 3HKC (10), 3HKE (10), 3HKD (10), 3N2K (11), and 3N2G (11)).
This article is protected by copyright. All rights reserved.

Based on the increased crystal structures related to the colchicine site, an increasing number of computer modeling
studies of CSIs have been reported (12-14). However, that all the available crystal structures of tubulin protein
complex have low resolution (3.5~4.0 Å, Table S1) leads to less success in structure-based drug discovery. The way
by employing combination of 3D-QSAR, docking studies and molecular dynamics simulations actually improved
accuracy of computer model (15). Additionally, there are very few reports about discovering novel tubulin inhibitors
using curated design and a series of structure-based virtual screening (VS) methods. Therefore, in this manuscript,
we would attempt to design novel CSIs based on the widely-used tubulin crystal complex 1SA0 as the drug target
template.
In order to discover novel tubulin inhibitors, the detailed investigation of the compounds that were collected
from recently published papers for structure–activity relationship (SAR), together with binding site analyses, was
displayed in Chart 1. On the basis of the above results, the curated process of designed compounds was performed
on the seven basic scaffolds (Figure S1). As to these compounds with structural diversity, we developed a
hierarchical strategy for structure-based virtual screening (VS) by integrating different computational methods
(
Figure 2). Briefly, the whole procedure could be depicted in the successive five steps: (1) On the basis of the SAR
analysis, designing 691 novel small molecules covering seven basic scaffolds; (2) Predicting the binding poses of
ligands using the CDOCKER docking program (16); (3) Optimizing and calculating the binding energy of the
docked ligand in protein binding pocket by molecular mechanics (MM) force field in combination with generalized
born surface area (GB/SA) implicit solvent (17); (4) With the help of those 115 reported active tubulin inhibitors
(
Figure S2) extracted from recent papers, building 3D-QSAR models located in the colchicine site using the
Discovery Studio suite 3.5 (18); (5) The establishment of pharmacophore modeling based on ten docked ligands
with high activity in the binding site, also built by the Discovery Studio suite 3.5. In each step, the predefined filters
were applied to eliminate the undesired hits, and it was expected that such hierarchical protocol can maximally
reduce the false positives produced by the manual experience-based judgment, which finally led to that the top 16
novel molecules were obtained. In the next section, we would elaborate this series of processes above in detailed.
This article is protected by copyright. All rights reserved.

RESULTS AND DISCUSSION
Investigation into the SAR of Compounds with Various Scaffolds
As said above, many recently published good works had been selected for optimizing our strategy, particularly
continuous studies from the same group. We retrieved the term “tubulin” on the J. Med. Chem. website, and curated
sixteen important medicinal chemistry papers (19-34) carefully for analysis of the structure-activity relationship
(
SAR) and the binding pose of these ligands (Chart 1). From the year 2004 to 2012, Silvestri group (19-24) reported
a series of colchicine site inhibitors with arylthioindole skeleton, and recently their attentions were shift to the
replacement of the indole ring with thiazole (25, 26). Kremmidiotis group (27) depicted a class of tubulin
polymerization inhibitors derived from the benzofuran scaffold, in which compound 95 (BNC105, Figure 1) had
been investigated under the clinical trial. Chang group (28, 29) reported that compounds containing 5-amino-2-
aroylquinolines skeleton possessed the inhibitory activities towards tubulin polymerization. In recent years Li group
(
33,34) reported that a class of 4-Aryl-2-benzoyl-imidazoles (ABI-III) were developed as as Tubulin Polymerization
Inhibitors. In addition, several kinds of CSIs with various aromatic ring system as basic skeleton (30-32) were also
reported.
Due to the inherent flexibility of the binding site, colchicine-binding pocket could accommodate various structure
ligands. In aim to probe the role of structural diversity of these compounds as CSIs, the molecular docking studies of
these optimized compounds from the above literatures were performed based on the common tubulin protein crystal
complex (PDB code: 1SA0), and subsequently the reasonable conformations of these compounds were manually
selected. Referring to the obtained poses located at the binding site, we artificially divided the modification space
into four small-size regions (A: the modifications of the linker; B and C: the modifications of substituents of the
parent nucleus; D: the modifications of the trimethoxybenzoyl “headgroup”). While this kind of simple classification
appeared to be arbitrary, it actually was in much favor of exploring the SARs of these compounds.
It is worthwhile noting that these compounds seemed to possess much of the same structural features when
focusing the binding pose at the colchicine-binding site (Chart 1). Firstly, almost all the optimized compounds
contained trimethoxybenzoyl moiety, which could be inserted into the same position of the active pocket. It is the
prerequisite for inhibition of tubulin polymerization. Focusing on D region, the methoxy group at the para position
This article is protected by copyright. All rights reserved.

can largely influence tubulin inhibitory activity (20), and a halo or a nitro group at the ortho position can moderately
increase potency (30). Second , how to carry on the modifications of B and C region heavily depended on the
structure of the parent nucleus (indol, thiophene, thiazole, or otherelse). The common electron-donating group such
as hydroxy, methoxy, amino or halo substituent can be employed at this site for enhancing potency (19-21, 23, 27-
3
0). The difference of modification on B and C area between these lead compounds was merely the various position
of the parent core. Another interest point was that the B region could accommodate a large-size group, no matter
hydrophobic or hydrophilic substituent (27). Besides, the alternation of the parent nucleus could circumvent tumor
resistance (31). The compounds with mono-aromatic ring system (thiazole or thiophene) might possess much of
variation on alternative binding pose (26, 27), thus this needed to be discussed carefully in the design process of
novel CSIs. Thirdly, the SAR studies (19-24) on the linker of a series of arylthioindole compounds revealed that A
region tended towards the rigid moiety, such as ketone or amide, while the methylene can be tolerated (23).
However, as to compounds with mono-aromatic ring as the parent nucleus, the analysis of their binding pose
displayed that the mono-aromatic ring could act as a role of the linker group (27, 33, 34).
Design CSIs based on the SAR Study and Limited Laboratory Conditions
In this section, we delineated our approach for the design of novel CSIs (Figure 3). It started with the employment
of seven basic scaffolds (A-series: indol; B-series: thionaphthene; C-series: benzofuran; D-series: imidazole; E-
series: parazole; F-series: benzimidazole; G-series: 1,3,4-Oxadiazole). A-, B, and C-series were the common basic
scaffolds for developing tubulin inhibitors, while the rest could be regarded as new basic entities. Besides, for E- and
G-series, there were a large number of works (35-40) reported in our group, which can provide the practical basis for
synthesizing these two kinds of compounds in the future.
As to A-series, we retained the bromine or nitro group at the ortho position of trimethoxybenzoyl headgroup (D
region), the bromine or methoxy group at the C-5 position of the indol ring (C region), and methoxy/ethoxycarbonyl
group at position 2 of the indole nucleus (B region) on the basis of the SAR studies. Considering that the rigidity of
the linker bridge (A region) and the aromatic moiety at B region would enhance oral bioavailability of compounds,
we introduced the four various length of rigid linkages (Figure 3, A-series) for bridging the indole nucleus and
trimethoxybenzoyl headgroup, and attached the 1,3,4-Oxadiazole aromatic system that was easy to be synthesized
This article is protected by copyright. All rights reserved.

to the C-2 position of the parent core. After carefully checking whether this kind of molecules had been reported
previously, the number of the molecules belong to A-series reached up to 143. Compared B and C series to A series,
the only difference focused on additional amino group that was closely related to the parent core (21). The number
of B and C series was 152 and 153, respectively. For the mono-aromatic ring (D-, E-, G-series), the designed
process seemed to be slightly complicated, due to docking study indicated that the binding conformation of these
classes of compounds in the active pocket would vary with the modification at the different position. As to D-series,
this sort of molecules can be further classified into three blocks (D1: trimethoxybenzoyl group at the C4-position
and amino group at the C-2 position of the imidazole nucleus; D2: trimethoxybenzoyl group at the C4-position and
amino group at the C-5 position; D3: trimethoxybenzoyl group at the C2-position), according to the grafted position
of trimethoxybenzoyl and amino group on the imidazole ring. In addition, owing to the small size of the imidazole
core, it was surmised that this scaffold could spin in the deep active pocket (33, 34). Thus, it was a bit difficult to
differentiate the modification at B or C region clearly. This stimulated us to consider the design of compounds
1
32~143, mainly due to that it was not sure at which position of these indole-pyrazole analogs the
trimethoxybenzoyl group should be placed. The aromatic substituents (m-OCH -C H , p-F- C H , and so on) were
3
6
4
6
4
selected from the corresponding reference (25, 26). After checking carefully, the total of D-series compounds
amounted to 105. As the same with D-series in the design process, 36 E-series compounds, 72 F-series compounds,
and 30 G-series compounds were prepared according to each procedure. In totally, there were 691 novel molecules
generated as tubulin inhibitors.
Performance of Docking-Based and Binding Energy-based VS
There were three common docking protocols (Glide 5.8 (41), Surflex (42), and CDOCKER (16)) employed in our
lab. However, comparison of molecular docking studies on in-house compounds as B-raf inhibitors using these three
software above, the obtained result by CDOCKER 3.5 exhibited better than the other two docking protocols.
Therefore, we artificially adopted CDOCKER protocol as main docking-based VS platform. Initially all the 692
compounds (691 designed molecules and the positive control colchicine) were docked into the reported X-ray
structures of Tubulin (PDB code: 1SA0) and generated 6920 binding poses. According to the values of
CDOCKER_INTERATION_ENERGY and the possibility of gain for each binding pose, the most optimized
binding conformation of each molecule had been obtained. Compared with the control molecule colchicine (Table
This article is protected by copyright. All rights reserved.

S2), we roughly ranked 322 molecules, of which CDOCKER_INTERATION_ENERGY ranged from -72.02 to -
5.4538 kcal/mol. As expected, the docking energy of the best molecule (D-101) was much lower than that of the
5
control molecule (-55.4257 kcal/mol).
On the basis of molecular docking analyses, the obtained 323 molecules (322 + colchicine) were optimized, and
evaluated by the binding energy term between protein and ligands. All the calculation was performed on the
Discovery Studio suite 3.5, in which the binding energy protocol was involved with two key settings: in situ ligand
minimization and implicit solvent (43). Therefore, all the 323 molecules were firstly simulated in the setting of in
situ ligand minimization using the CHARMm force field, and then added into poisson−boltzmann surface area
(
PB/SA) model for binding energy calculation. The calculated result were displayed in Table S3. Comparison of
these compounds with the positive molecule colchicine, we roughly ranked 215 molecules, of which
CDOCKER_INTERATION_ENERGY ranged from -167.862 to -79.2536 kcal/mol. As expected, the binding
energy of the best molecule (B-94) was much lower than that of the control molecule (-79.2536 kcal/mol).
Performance of 3D-QSAR-Based VS
3
D-QSAR model can correlate compound activities with interaction fields calculated based on protein
crystallography or molecule superimposition. The success of any QSAR model depended on several factors, such as
accuracy of input data, the choice of descriptors and statistical methods for modeling and for validation (44). Firstly,
with regard to the input data, a total of 116 molecules (114 arylthioindoles/thiazoles/benzo[b]furan, plus the two
control molecules colchicine and CA4) as potent tubulin inhibitor was collected from the continuous work of
Silvestri group (19-26) and Kremmidiotis group (27). In this work, the in vitro biological activities of these
compounds for inhibition of tubulin polymerization were converted into the corresponding pIC50 (-log Ki) values
(
Table S4), which were used as dependent variables in the 3D-QSAR model. Due to the discrepancy of the activity
measure used in different papers, we simply modified the IC50 values of compounds 86~97 (25), 98~108 (27),
09~113 (26) by comparing CA4 activity values between in refs (16-21) and refs (22-24) (Table S4). In addition, it
1
was acceptable that the activity value of these compounds ranged from 0.67 to 15 μM. The collected compounds
were randomly divided into a training set (including 93 compounds) and a test set (including 23 compounds) in an
approximate ratio of 4:1. Second , molecular alignment of compounds played a crucial role in the development of
This article is protected by copyright. All rights reserved.

3
D-QSAR models (45). So far, there were two common methods used for molecular alignment of compounds
(
ligand-based superimposition and docking-based superimposition). Owing to structural diversity of the selected 116
compounds, we performed molecular docking study of these molecules, and then selected their optimal pose for
molecular alignment, which can be clearly demonstrated in Figure 4a.
According to the 3D-QSAR procedure in the Discovery Studio suite 3.5, the corresponding model was built easily.
To validate whether a QSAR model was robust for prediction of unknown molecules, three statistical parameters in
2
2
DS 3.5 including especially the cross-validated correlation coefficient (q ), correlation coefficient (r ), and root
mean square error (RMS) should be evaluated. For this 3D-QSAR studies of 116 CSIs, acceptable correlations were
2
observed in the obtained model demonstrated by the toleranted value of q (0.452) and good correlation coefficient
2
r (0.921). As shown in Figure 4b, the correlations between the experimental and the predicted activities for both the
training and test sets for the optimal 3D-QSAR model. Obviously, all the points are rather uniformly distributed
around the regression line, suggesting the satisfactory predictive capability of the models. The detailed values
depicted by the QSAR model were list in Table S4. The 3D-QSAR result was used to be represented as 3D
coefficient contour (Figure 4c and 4d). It showed regions located in the colchicine-binding pocket where variations
of the steric and electrostatic nature of compounds can lead to the increase or decrease in the activity. To aid in
visualization, the positive molecule colchicine was displayed superimposed with the contour maps. For example, to
increase the activity of colchicine, we should add polar groups with positive electrostatic potential in the red area,
and add that with negative electrostatic in the blue area. Figure 4d showed van der Waals (VDW) interactions of
3
D-QSAR. Green indicates positive coefficients; violet indicates negative coefficients. In the same way, to increase
activity, a new molecule should have strong VDW attraction in the green area and weak VDW attraction in the
violet area. As shown in Figure 2, The activities of the 225 molecules for tubulin inhibition was predicted based on
the 3D-QSAR model. Compared with the control molecule colchicine (Table S5), we roughly ranked 51 molecules.
The best compound was D-29 with the predicted pIC value of 6.07586.
5
0
This article is protected by copyright. All rights reserved.

Performance of Pharmacophore-Based VS
Forty-eight molecules was a large number for the result of VS. Therefore, there was a need for reducing the
scope of these novel CSIs. In this section, we selected ten of the most active compounds of each literature (19-27)
discussed above for generating the complex-based pharmacophore models. These molecules were compound 24, 64,
7
2, 81, 82, 85, 94, 100, colchicine and CA4 (Table 1). Meanwhile, we picked up 53 active molecules from the
previous 116 QSAR molecules as the active set of pharmacophore establishment, and prepared 21 inactive in-house
benzimidazole molecules as the inactive set (Table S2). For each complex, all the generated pharmacophore models
could be ranked by the selectivity scores, and the prediction capacity of each pharmacophore model to distinguish
the active inhibitors from the inactive molecules was evaluated by the Accuracy values. One compound could
produce one optimal pharmacophore model. Finally, ten kinds of tubulin colchicine-site pharmacophore models
based on various structure molecules were successfully built (Table 1, pharmacophore a-j). All pharmacophoric
features were derived from the protein−ligand interactions. One other thing to note was that while all the models can
possess several same kinds of features, their coordinates in the active pocket were different.
According to the screening steps provided in Figure 2, the 51 screened molecules were imported into these
pharmacophore models mentioned above for virtual evaluation. To our disappointment, the pharmacophore h model
built based on compound 100 could not work as the others. Therefore, Table 2 just listed the final results screened
by the other nine pharmacophore models. Due to each pharmacophore model possessed different features that can
reveal the structural diversity of these 51 molecules, we ranked all the results by the consensus score that combined
with the fitvalues of nine pharmacophore models. Briefly, if the fitvalue of each compound was in the top 20% of
one pharmacophore’s evaluation system, the consensus number would automatically be added one. Thus the
consensus score was a integer, and the high the value was, the better the compound performed. Analysis of the
results in Table 2 showed most of compounds in front of colchicine could covered several pharmacophore models.
The best molecule A-46 possessed high fitvalues in a series of pharmacophore models (pharmacophore b-e). In
addition, chemical structures of the ultimate 16 molecules was demonstrated in Figure 6. Strikingly, the parent
nucleus of these compounds come from A, B, C, D, and E-series, displaying their structural diversity.
This article is protected by copyright. All rights reserved.

Validation of Availability of The Design Strategy
Among the 16 designed molecules, compound A-132 with in silico moderate activity was synthesized as tubulin
polymerization inhibitor in this manuscript (Scheme 1). The detailed synthesis was described in the Experimental
Section. Figure 7b showed the antiproliferative activities of compound A-132 and the control molecule colchicine,
and obviously compound A-132 can be comparable to colchicine in these three cancer cell assay, particularly
performing better than colchicine in the antiproliferative activities against HepG2 cells. In addition, we examined
the inhibitory activity of compound A-132 for tubulin polymerization assay, of which the IC50 value could reach up
to 1.64 μM, slightly superior to that of colchicine.
Immunofluorescence assay of compound A-132 was performed on HepG2 cells (hepatocelluar carcinoma cell) by
using anti-tubulin antibody and DAPI nuclear staining agent. As shown in Figure 7c, the microtubule network
exhibits normal arrangement and organization in HepG2 in the absence of drug treatment (Figure 7 c1). As expected,
with the increasing concentration of compound A-132 ranging from 0.64 to 6.4 μM after 24 h drug treatment,
HepG2 cell seemed obviously to shrink in size, and cellular microtubule significantly displayed depolymerization in
a concentration-dependent manner(Figure 7 c2 ~ c4). Molecular docking of compound A-132 was performed in the
tubulin active pocket (PDB code: 1SA0, Figure 8). As the same with these known CSIs mentioned in the SAR
analysis section, compound A-132 can be inserted better into the active pocket of tubulin. Although it did not
contain the trimethoxybenzoyl moiety in the subpocket formed around the key residues (CYS241, LEU242, and
ILE378) when compared with colchicine (Figure 8b), this molecule could make full use of the space in the binding
site, particularly the modification of 3-methoxyphenyl and methyl in the C region. Additionally, the predicted
binding mode of compound A-132 was in accordance with the speculated in advance. Thus, insight into deep
analysis of the binding site can contribute to the better design of CSIs.
CONCLUSION
In the present study, we manually curated a small-size CSI database composed of 691 molecules covering several
heterocycle scaffolds in the field of medicinal chemistry, and subsequently performed structure-based virtual
screening against colchicine binding site. In order to improve the success of hits by VS, the design process of these
compounds appeared a bit conservative, of which their structure might be more close to the reported CSIs. For
This article is protected by copyright. All rights reserved.

making a balance between the computational cost and prediction efficiency, we employed a hierarchical strategy by
integrating various computational methods in an increasing order of complexity. Besides, homogeneous
distributions of the number of designed compounds from different series at each stage of VS (Figure 9) validated the
reliability of this tactics or the rationality of various structure compounds. Combining molecular docking, 3D-QSAR
model, and complex-based pharmacophore searching based on tubulin protein structure, sixteen novel compounds
(
Figure 6) with various structural features were obtained, among which compound A-132 has been synthesized, and
preliminary biological assay of A-132 for tubulin polymerization inhibition can provide the accuracy of this kind of
rational drug design. So far, the other fifteen target compounds are being synthesized in our lab, and they can be
served as the starting points for the development of novel tubulin inhibitors.
METHODS AND MATERIALS
Computational Section
The SAR Analysis
All the compounds used in this study were extracted from a series of recently published literature (19-27). The
compounds displayed in Chart 1 were picked up manually, and each molecule exhibited most potent activity in the
corresponding paper. The protein-ligand interaction maps was depicted based on the optimal pose of this compound
from the docking study. Partitions of the modification space were roughly from individual subjective judgment, and
they were not exact. However, this was just to aid the SAR analysis.
Molecular Docking
Molecular docking of compounds into the three dimensional X-ray structure of human tubulin complex (PDB
code: 1SA0) was carried out using the Discovery Studio (version 3.5) as implemented through the graphical user
interface DS-CDOCKER protocol. The 3D structure of tubulin (1SA0) in docking study was downloaded from
Protein Data Bank. The three-dimensional structures of the aforementioned compounds were constructed using
Chem. 3D ultra 12.0 software [Chemical Structure Drawing Standard; Cambridge Soft corporation, USA (2010)],
This article is protected by copyright. All rights reserved.

then they were energetically minimized by using MMFF94 with 5000 iterations and minimum RMS gradient of 0.10.
All bound waters and ligands were eliminated from the protein and the polar hydrogen was added to the proteins.
Each compounds would retain 10 poses, and were ranked by CDOCKER_INTERACTION_ENERGY. The
selection of the conformation of each compound would partly depend on the binding mode in the active pocket,
when compared those known inhibitors, such as CA4, colchicine, and so on.
3D-QSAR
3
D QSAR model provided useful information about the correlation between the molecular fields and the activity.
All the operations were also performed on the the Discovery Studio suite 3.5. 116 compounds were carefully picked
up from the literature mentioned above, and it needed considering the activity range of these compounds. However,
their IC50 values could not span 4 orders of magnitude, mainly due to the less active molecule reported in these
papers had not enough and reliable data. The docking study was carried out to obtain reasonable molecular
alignments for building 3D-QSAR model. In addition, the parameter for cross-validation was set to 10.
Pharmacophore Model
The complex-based pharmacophore model can be obtained directly from ligand-protein co-crystal structure that
can reflect more reliable information about the structural features required for the relating biological potency in the
binding site. The complex-based pharmacophore models for ten tubulin complexes selected from the corresponding
docking study were generated by using the Receptor−Ligand Pharmacophore Generation (RLPG) protocol in DS3.5.
Generally pharmacophoric features include hydrogen bond acceptor/donor (HBA/HBD), hydrophobic (HYD),
negative/positive ionizable (PI/NI), and ring aromatic (RA) features in the DS suite. For each pharmacophore model,
we provided the same set of active molecules(53 compounds selected from the QSAR study) and inactive molecules
(
21 in-house compounds, Chart S2) for the validation of pharmacophore models.
This article is protected by copyright. All rights reserved.

Experimental Section
Chemistry
The synthetic route of compound A-132 is outlined in Scheme 1. Firstly, compound s2 were prepared by
vilsmeier reaction. Second , Compound s3 were prepared according to the procedure reported by Debajyoti Saha et
al. with some modifications. Compound s4 were synthesized from compound s3 reacting with 1-(3-methoxyphenyl)
ethanone at room temperature. The mixture of compound s4 and hydrazine hydrate in CH CH OH was refluxed
3
2
under stirring for 8 h and the precipitate formed was filtered off, washed with cool ethanol to get compound s5.
Finally, compound s5 reacted with 3-methoxybenzoic acid, together with EDC·HCl and HOBt. The mixture was
refluxed under stirring for 24 h. The crude product was crystallized with ethanol to give target compound A-132.
1
Compounds A-132 were fully characterized by H NMR, ESI-MS and elemental analysis.
Chemistry general
All chemicals and reagents of analytical grade used were purchased from Aldrich (USA). Melting points
1
(
uncorrected) were determined on a XT4MP apparatus (Taike Corp, Beijing, China). All the H NMR spectra were
recorded on a Bruker DPX 300 model Spectrometer at 25 with TMS and solvent signals allotted as internal
standards, and chemical shifts were reported in ppm (d). ESI-MS spectra were recorded on a Mariner System 5304
Mass spectrometer. Elemental analyses were performed on a CHN-O-Rapid instrument and were within 0.4% of the
theoretical values. TLC was performed on the glass-backed silica gel sheets (Silica Gel 60 GF254) and visualized in
UV light (254 nm). Column chromatography was performed using silica gel (200-300 mesh) eluting with ethyl
acetate and petroleum ether.
5
-methoxy-1H-indole-3-carbaldehyde (s2) was prepared according to the following procedure. To a solution of
DMF (20 mL) was added POCl (10 mL) in one portion at 0 , and then gradually heated up to room temperature
3
for 2 h. A solution of 5-methoxy-1H-indole (s1, 3.6 g, 24.7 mmol, dissolved in 10 mL DMF) was then added to the
mixture at that temperature. The mixture was further heated at 60 for 1 h. After cooled to room temperature, the
mixture was poured into ice-water slowly and then adjusted pH to 9.0 with NaHCO . The precipitated cake was
3
This article is protected by copyright. All rights reserved.

collected by filtration, after rinsing and drying with water. The crude product was purified by column
chromatography (1 : 2 V (EtOAc)/V (hexanes) to get compound s2, a light yellow powder, yield 80%.
5
-methoxy-1-methyl-1H-indole-3-carbaldehyde (s3) To a suspended solution of NaH (0.60 g, 60% dispersion in
mineral oil, 25 mmol) in THF (5 mL), 5-methoxy-1H-indole-3-carbaldehyde (s2, 1.76 g, 10 mmol) dissolved in THF
5 mL) was added dropwise at 0 . The heterogeneous mixture was stirred at 0 for 15 min and 1 h at room
(
temperature. The mixture was then cooled to 0 with iodomethane (0.84 mL, 13.2 mmol), and allowed from warm
to room temperature. After 30 min, the reaction mixture was cooled to 0 , quenched with saturated NH Cl (20 mL),
4
and extracted with ether (35 mL). The organic layers were combined, washed with brine, dried over anhydrous
Na SO and concentrated in vacuo. The resulting oil was purified by flash column chromatography to provide 5-
2
4
methoxy-1-methyl-1H-indole-3-carbaldehyde, a light yellow oil, yield 67.8%.
(
E)-1-(5-methoxy-1-methyl-1H-indol-3-yl)-3-(3-methoxyphenyl) prop-2-en-1-one (s4) To a stirred solution of 1-
(
(
3-methoxyphenyl) ethanone (1 mmol) and 5-methoxy-1-methyl-1H-indole-3-carbaldehyde (s3, 1 mmol) in ethanol
30 mL), 6 mol KOH (4 mL) was added and the reaction mixture was stirred until the solids fully formed. The
products were filtrated and washed carefully with ice water and cool ethanol, and purified by crystallization from
ethanol in refrigerator to give (E)-1-(5-methoxy-1-methyl-1H-indol-3-yl)-3-(3-methoxyphenyl) prop-2-en-1-one
(
compound s4) a yellow powder, yield 60%.
5
-methoxy-3-(3-(3-methoxyphenyl)-4, 5-dihydro-1H-pyrazol-5-yl)-1-methyl-1H-indole (s5) To a solution of
compound s4 (1 mmol) in isopropanol (5 mL) hydrazine hydrate (0.2 mL, 4 mmol) was added. The mixture was
refluxed under stirring for 8 h, stored at 4~5 for 24 h, and the precipitate formed was filtered off, washed with
cool ethanol. The synthesized compound was purified by crystallization from ethanol in refrigerator to give 5-
methoxy-3-(3-(3-methoxyphenyl)-4, 5-dihydro-1H-pyrazol-5-yl)-1-methyl-1H-indole (s5), yellow powder, yield
4
0%.
(
5-(5-methoxy-1-methyl-1H-indol-3-yl)-3-(3-methoxyphenyl)-4, 5-dihydro-1H-pyrazol-1-yl)(3-methoxyphenyl)
methanone (A-132) To a solution of compound s5 (1 mmol) in dichloromethane (5 mL) and 3-methoxybenzoic acid
1 mmol) was added, together with EDC·HCl (1.2~1.5 mmol) and HOBt (1.2~1.5 mmol). The mixture was refluxed
(
under stirring for 24 h. After completion of the reaction, the contents were cooled, and then evaporated to dryness in
This article is protected by copyright. All rights reserved.

vacuo. Aqueous hydrochloric acid (0.1 M, 30 mL) was added and the mixture extracted with ethyl acetate (3 × 5
mL). The combined ethyl acetate layers were back-extracted with saturated sodium bicarbonate (1 × 5 mL) and
brine (1 × 5 mL), dried over MgSO4, filtered, and evaporated in vacuo. The residue was crystallized from ethanol to
1
give target compound A-132. White crystals, yield 76.7 %. Mp: 128-130 ˚C. H NMR (400 MHz, DMSO, δ ppm):
7
1
3
.45 – 7.34 (m, 5H), 7.30 (t, J = 4.4 Hz, 3H), 7.12 – 7.02 (m, 2H), 6.88 (s, 1H), 6.79 (dd, J = 8.84 Hz, J = 8.84 Hz,
1 2
H), 6.01 (dd, J = 11.64 Hz, J = 11.72 Hz, 1H), 3.88 (dd, J = 18 Hz, J = 17.92 Hz, 1H), 3.78 (d, J = 2.4 Hz, 6H),
1
2
1
2
+
.71 (s, 3H), 3.58 (s, 3H), 3.31 (dd, J = 18 Hz, J = 17.92 Hz, 1H). MS (ESI): 470.53 [M+H] ; Anal. Calcd. for
1
2
C H N O : C, 71.62; H, 5.80; N, 8.95; O, 13.63; Found: C, 71.63; H, 5.79; N, 8.96; O, 13.64.
2
8
27
3
4
Biological assay
Cell proliferation assay (Cell viability was assessed by MTT assay)
We evaluated the antiproliferativivities of compounds A-132 against A549 (carcinomic human alveolar basal
epithelial cell), HepG2 (human hepatocelluar carcinoma cell line), and MCF-7 (breast cancer) cancer cells with
analysis. Cell proliferation was determined using the MTT dye (Beyotime Inst Biotech, China) according to the
3
instructions of manufacture. Briefly, 5 x 10 cells per well were seeded in a 96-well plate, grown at 37 ℃ for 12 h.
Subsequently, cells were treated with compounds (A-132 and Colchicine) at increasing concentrations in the
presence of 10% FBS for 24 h. After 10 μL MTT dye was added to each well, cells were incubated at 37 ℃ for 3~4
h. Then all the solution in the wells was poured out and 150 μL DMSO was added to every well. Plates were read in
a Victor-V multilabel counter (Perkin-Elmer) using the default europium detection protocol. Percent inhibition or
GI50 values of compounds were calculated by comparison with DMSO-treated control wells. The results are shown
in Figure 7b.
Tubulin Polymerization assay
To evaluate the effect of the compounds on tubulin polymerization in vitro, different concentrations of the two
compounds (colchicine and A-132) were preincubated with 10 μM bovine brain tubulin in glutamate buffer at 30 ℃
and then cooled to 0 ℃. After addition of 0.4 mM GTP, the mixtures were transferred to 0 ℃ cuvettes in a recording
This article is protected by copyright. All rights reserved.

spectrophotometer and warmed to 30 ℃. Tubulin assembly was followed turbidimetrically at 350 nm. The IC50 was
defined as the compound concentration that inhibited the extent of polymerization by 50% after 20 min incubation.
Immunocytochemistry
HepG2 cell line plated on coverslips was treated with indicated concentration of test agents for the selected
treatment duration. After treatment, Cells were washed twice in PBS and fixed in 2% paraformaldehyde for 20 min
at room temperature. Cells were permeabilized in PBS with 0.1% Na-Citrate, 0.1% TritonX-100 at room
temperature for 15 min. Nonspecific binding sites were blocked by incubating the cells in 5% BSA. Cells were then
incubated with tubulin antibody (Beyotime, AT108, China) followed by conjugated secondary IgGs (Beyotime,
A0539, China) and DAPI (1 μg/ml). After incubation, cells were washed with PBS. Cellular microtubules were
observed with an Olympus BX50 fluorescence microscope.
Figure/Table/Chart/Scheme Captions
Figure 1 Chemical structures of CSIs in clinical trials currently.
Figure 2 Flowchart of the whole structure-based hierarchical virtual screening strategy.
Figure 3 Design protocol based on different scaffolds. The A, B, C, and D region represented the
modified space in the active pocket of tubulin.
Figure 4 3D-QSAR model. a) The molecular alignment of 116 compounds in the data set based
on the docking study; b) Predicted versus actual pIC for the 3D-QSAR model. The training set
5
0
was marked in black dot and test set marked in red asterisk; c) 3D-QSAR electrostatic contour
plots with the combination of colchicine. Red contours indicate regions where negative charges
increase activity, and blue contours indicate regions where positive charges increase activity. d)
3
D-QSAR steric contour map with the combination of colchicine. Green contours indicate
regions where bulky groups increase activity, and violet contours indicate regions where bulky
groups decrease activity.
Figure 5 Ten complex-based pharmacophore models developed for virtual screening (VS). green:
hydrogen-bonding acceptor; pink: hydrogen-bonding donor; cyan: Hydrophobic; grey: excluded
volume.
Figure 6 Chemical structures of the ultimate sixteen novel compounds.
Figure 7 a) Chemical structures of compound A-132; b) In Vitro inhibition of tubulin
polymerization, growth of three common cancer cells (A549, HepG2, and MCF-7) by compound
This article is protected by copyright. All rights reserved.

A-132 and colchicine; c) Immunofluorescence microscopy detection in HepG2 cells; tubulin (red)
and DAPI nuclear staining (blue); cells were treated for 24 h with compound A-132 (c2~c4).
Figure 8 a) Focused view of predicted binding model of compound A-132 (green, stick) in
context of surrounding residues (cyan, stick); b) Focused view of predicted binding model of
colchicine (green, stick) in context of surrounding residues (cyan, stick).
Figure 9 Distributions of the number of designed compounds belonged to different series for
each stage of virtual screening (VS).
Table 1 Numbers of total features, selectivity score, and accuracy of the pharmacophore models
based on ten potent tubulin inhibitors.
Table 2 The screening results obtained by using the pharmacophore models based on nine active
tubulin inhibitors.
Scheme 1. Synthesis of compounds A-132. Reagents and conditions: (i) DMF, POCl ; (ii) NaH,
3
CH3I, THF; (iii) EtOH, NaOH, r.t.; (vi) NH
2
·NH
2
·H O, EtOH, reflux; (v) EDC·HCl, HOBt,
2
Cl , reflux.
2
H2
This article is protected by copyright. All rights reserved.

Figure 1 Chemical structures of CSIs in clinical trials currently.
This article is protected by copyright. All rights reserved.

Figure 2 Flowchart of the whole structure-based hierarchical virtual screening strategy.
This article is protected by copyright. All rights reserved.

Figure 3 Design protocol based on different scaffolds. The A, B, C, and D region represented the
modified space in the active pocket of tubulin.
This article is protected by copyright. All rights reserved.

Fig. 4 3D-QSAR model. a) The molecular alignment of 116 compounds in the data set based on
the docking study; b) Predicted versus actual pIC for the 3D-QSAR model. The training set
50
was marked in black dot and test set marked in red asterisk; c) 3D-QSAR electrostatic contour
plots with the combination of colchicine. Red contours indicate regions where negative charges
increase activity, and blue contours indicate regions where positive charges increase activity. d)
3
D-QSAR steric contour map with the combination of colchicine. Green contours indicate
regions where bulky groups increase activity, and violet contours indicate regions where bulky
groups decrease activity.
This article is protected by copyright. All rights reserved.

Figure 5 Ten complex-based pharmacophore models developed for virtual screening (VS). green:
hydrogen-bonding acceptor; pink: hydrogen-bonding donor; cyan: Hydrophobic; grey: excluded
volume.
This article is protected by copyright. All rights reserved.

Figure 6 Chemical structures of the ultimate sixteen novel compounds and the corresponding
predicted IC values against tubulin assembly.
5
0
This article is protected by copyright. All rights reserved.

Figure 7 a) Chemical structures of compound A-132; b) In Vitro inhibition of tubulin
polymerization, growth of three common cancer cells (A549, HepG2, and MCF-7) by compound
A-132 and colchicine; c) Immunofluorescence microscopy detection in HepG2 cells; tubulin (red)
and DAPI nuclear staining (blue); cells were treated for 24 h with compound A-132 (c2~c4).
Figure 8 a) Focused view of predicted binding model of compound A-132 (green, stick) in
context of surrounding residues (cyan, stick); b) Focused view of predicted binding model of
colchicine (green, stick) in context of surrounding residues (cyan, stick).
This article is protected by copyright. All rights reserved.

Figure 9 Distributions of the number of designed compounds belonged to different series for
each stage of virtual screening (VS).
This article is protected by copyright. All rights reserved.

Chart 1. Analysis of structure-activity relationship (SAR) for potent tubulin inhibitors
extracted from recently published papers.
This article is protected by copyright. All rights reserved.

This article is protected by copyright. All rights reserved.

Table 1 Numbers of total features, selectivity score, and accuracy of the pharmacophore models
based on ten potent tubulin inhibitors.
Pharmacophore Summary
Pharmacophore
Based-Ligand
Compound 24
Compound 64
Compound 72
Compound 81
Compound 82
Compound 85
Compound 94
Compound 100
Colchicine
Total features
selectivity score
Accuracy
0.768
0.989
0.897
0.966
0.934
0.932
0.965
0.841
0.945
0.866
Pharmacophore_a
Pharmacophore_b
Pharmacophore_c
Pharmacophore_d
Pharmacophore_e
Pharmacophore_f
Pharmacophore_g
Pharmacophore_h
Pharmacophore_i
Pharmacophore_j
DHHHH
2.1117
DHHHHH
HHHHH
3.3300
1.4740
DHHHH
2.1117
DHHHH
2.1117
HHHHHH
HHHHHH
DHHHH
2.6924
2.6924
2.1117
ADHHHH
AHHHH
4.5151
CA4
1.4740
Table 2 The screening results obtained by using the pharmacophore models based on nine active
tubulin inhibitors.
Nine Protein–Ligand-Based Pharmacophores
Cmpd
.
Conse
nsus
Phar
mac
g
Pharm Pharm Pharm Pharm Pharma Pharm
Phar Pharma
ac a
ac b
ac c
ac d
c e
ac f
mac i
c j
0
1
.460 0.633 0.410 0.5076
0.0929
131
A-46
0
0
0
0
9
85
174
0
312
99
0
.455
0.809 0.0002
552 62533
0.0693
767
A-16
A-87
0
0
0
0
0
0
0
0
7
7
6
2
0
0
0
0
.333 0.4010
This article is protected by copyright. All rights reserved.

1
27
0
92
0
0
0
.375
0.922 0.585
582 478
0.5525
1
B-4
C-6
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
7
7
7
6
6
6
6
6
6
6
6
6
6
6
8
46
.816
55
0.055 0.625
0.4150
79
0
0
1
1115
441
0
0
.162 0.0153
55
E-11
A-21
A-33
B-130
B-6
0
0
0
0
0
0
272
0
9
.022 0.131
0.0141
2
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
212
921
0
0
.039 0.112 0.0001
253
0
053
17043
0
0
.605
0.0008
70695
0
0
0
0
0
0
0
0
0
0
74
0
0
.0241
05
0
0
0
0
0
0
5
0
.2762
B-95
C-134
D-29
A-132
D-80
E-16
0
6
6
0
0
.761
57
0
0
1
0
3
.0013
543
0.035 0.0213
2673
0
0
799
.428 0.032
61
0.0012
891
0
0
3
1092
0
.5825
0
0
0
0
0
0
0
0
0
3
0
4
.326
73
0
0
Colchi
cine
0.042 0.086 0.626
6253 2241 284
0.682
723
0
This article is protected by copyright. All rights reserved.