Dynamic and Static Water Molecules Complement the TN16

Dec 15, 2015 - TN16 is one of the most promising inhibitors of α, β dimer of tubulin that occupies the cavity in the β-subunit located at the dimer...
4 downloads 11 Views 9MB Size
Subscriber access provided by UNIV MASSACHUSETTS WORCESTER

Article

Dynamic and Static Water Molecules Complement the TN16 Conformational Heterogeneity Inside the Tubulin Cavity Sarmistha Majumdar, Satyabrata Maiti, and Shubhra Ghosh Dastidar Biochemistry, Just Accepted Manuscript • DOI: 10.1021/acs.biochem.5b00853 • Publication Date (Web): 15 Dec 2015 Downloaded from http://pubs.acs.org on December 21, 2015

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Biochemistry is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Dynamic and Static Water Molecules Complement the TN16 Conformational Heterogeneity Inside the Tubulin Cavity Sarmistha Majumdar, Satyabrata Maiti, Shubhra Ghosh Dastidar* Bioinformatics Centre, Bose Institute, P-1/12 CIT Scheme VII M, Kolkata 700054 *Address for correspondence [email protected] Tel: +91-33-2569-3332 (O) +91-98301-63990 (M)

Funding statement: S Majumdar and S Maiti are supported by funding from Department of Biotechnology of Government of India (DBT, GOI) through the program of ‘Bioinformatics Centre of Excellence’. The high performance computing infrastructure that was available to SGD from DBT, GOI through the project no. BT/PR793/BID/7/370/2011 and has been supportive to this work.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abbreviations: CHARMM: Chemistry at HARvard Macromolecular Mechanics (a software for MD) PDB: Protein Data Bank MD: Molecular Dynamics QM: Quantum Mechanical TBE200I1, TBE200I1, TBE200P, TBUC, TBUS, TN16U: Trajectory names, listed in Table 3 ABT751, TN16, T138067: various ligands of Tubulin GTP/GDP: Guanosine Tri/Di Phosphate CGenFF: CHARMM Generalized Force Field CHARMM-GUI: Graphical User Interface for CHARMM (a webserver) PME: Particle Mesh Ewald ABNR: Adopted Basis Newton Raphson NAMD: NAnoscale Molecular Dynamics (a program for molecular dynamics simulation package) NPT: Normal Temperatute and Pressure NBO: Natural Bond Orbital MMGBSA: Molecular Mechanics Generalized Born and Surface Area GBMV: Generalized Born using Molecular Volume

ACS Paragon Plus Environment

Page 2 of 38

Page 3 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Abstract: TN16 is one of the most promising inhibitors of α, β dimer of tubulin that occupies the cavity in the β-subunit located at the dimeric interface, known as colchicine binding site. The experimentally determined structure of the complex, PDB code 3HKD, presents the conformation and position of the ligand based on the ‘best-fit’, keeping the controversy of other significant binding modes open for further investigation. Computation has already revealed that TN16 experiences fluctuations within the binding pocket, but the insight from that previous report was limited by the shorter windows of sampling and by the approximations on the surrounding environment by implicit solvation. This article reports that in most of the cases straightforward MMGBSA calculations of binding energy revealed gradual loss of stabilization which was inconsistent with the structural observations and thus it indicated the lack of consideration of stabilizing factors with appropriate weightage. Consideration of the structurally packed water molecules in the space between the ligand and receptor successfully got rid of such discrepancies between the structure and stability, serving as the ‘litmus test’ on the importance of explicit consideration of such structurally packed water in the calculations.

Such consideration has further evidenced a quasi-degenerate

character of the different binding modes of TN16 that has rationalized the observed intrinsic fluctuations of TN16 within the pocket; which is likely to be the most critical insight into its entropy dominated binding. QM calculations have revealed a relay of electron density from TN16 to the protein via water molecule in a concerted manner.

Key words: Tubulin, TN16, molecular dynamics, bridging water, rotational relaxation, QM

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Introduction: Microtubules are dynamic protein filaments which play fundamental role during the cell division. They also act as the key regulators in diverse array of cellular functions such as maintaining the cell shape, signalling, motility, transport etc. Repeated association of the hetero-dimer of the α and β tubulin subunits (~50 kDa each) in an alternating manner constitutes the hollow tube like architecture of microtubules.1-3 The microtubules experience “dynamic instability” as the polymerization and depolymerization of the heterodimers is active all the time and depending in the cellular requirement the equilibrium shifts towards the elongation.4-7 Strategic perturbation to this dynamics of the microtubule formation is a promising way to fight against cancer.

1, 8, 9

It could be done in two ways, (i) by disallowing

the microtubule polymerization and thus arresting the cell in pre-mitotic stage, (ii) stabilizing the microtubule and arresting the progress of the course of mitosis in the midway. Compounds like colchicine and vinca are of type (i); while paclitaxel, epothilones etc. are of type (ii).10-12 Since there are several sites for ligand binding on tubulin, the binding sites are identified by one of the promising compounds for each case and primarily those are paclitaxel, peloruside/laulimalide, vinca and colchicines binding sites.13-16 Recently in the year 2014, M. O. Steinmetz and his group reported a new binding site in tubulin for a set of microtubule destabilizing drugs and termed it as ‘maytansine-site’.17 Among all these sites, the colchicine binding site is located at the interface of the dimer and inhibitors suitably occupying it are reported to promote the curved conformation of the complex which is unable to fit into the tubular architecture.18,

19

The popular explanation of the locking of the

conformation into the curved state imposed by the ligand is that the presence of the ligand causes a steric restriction of the T7 loop of the β-subunit. The colchicine binding site is able to accept structurally diverse chemical compounds.18, 20, 21 Audrey Dorléans et al. explained in details from their experimental structures that how the ligands like ABT751, TN16, T138067

ACS Paragon Plus Environment

Page 4 of 38

Page 5 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

compete for and get stabilized within the same colchicine binding site but in a different manner compared with colchicine binding, although they are structurally different from each other. Dorléans et al. also provided information on the difference of the attachment points for all the ligands20 (Figure 1). Later, Molecular dynamics (MD) simulations pointed out the flexible nature of this binding pocket (colchicine domain), which can adjust depending upon the shape of the small molecules and where ligands can also oscillate within the pocket adding to the entropic benefit to the binding.22 Since the interference of the ligand with T7 movement is the widely accepted explanation of change in the conformational preference, the TN16 binding into tubulin is really intriguing as its shape is thinner than bulker colchicine which penetrates deeply into the cavity leaving only a residual portion is available at the α,βinterface to interact with T7. As reported in the previous MD study22 , there is a departure of the TN16 position from crystal structure to move into the core of the β-subunit and remains flexible. The advantage of deeply penetrating binding mode of TN16 relative to that of colchicine binding was understood from the report that TN16 binding is stabilized by enthalpy and entropy both, where colchicine binding is only entropy driven; literature data has been presented in Table 1. Indeed, there was a departure of the ligand’s position from the crystal structure when it went through MD calculations and that must be taken into account as a refinement of the structure because of the following reasons. Dorleans et al. mentioned in their article20 reporting the crystal structure (3HKD) that the placement of the TN16 in the PDB model was based on the Fobs-Fcalc and 2Fobs-Fact omit maps and the optimized polar interactions and underscored the possibility of alternative orientations. A close inspection reveals that the alternative conformations can vary as a function of several factors. For example, in the PDB model the E200 (by default all residues mentioned in the text correspond to β-subunit) side chain of β-subunit is within 3Å of one of the O atoms of the

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

TN16 (Figure 2). In this conformation, it could be reasonably assumed that the carboxylic acid group would be in neutral state to avoid repulsion as well as to get the advantage of Hbonding. Whereas, in an alternate situation it can experience its deprotonated or charged (-1e) state and this can either involve local adjustment to accommodate TN16 or may not bind TN16. Moreover, the size of the cavity, and the space between α and β subunit are sufficient to accommodate water molecules in absence of the ligand and so depending on the shape and size of the ligand, a few residual water molecules could be retained in the complex making significant impact on the stabilization or destabilization. Moreover, the crystallographic data was resolved at 3.7 Å, which is inferior to the quality of resolution (at least < 2 Å) required to see water molecules in the crystal20. So all these factors make the TN16 binding to Tubulin very intriguing and requires rigorous investigation in explicit solvent environment and on a relatively longer timescale that allows the ligand to explore its all possible binding modes. The previous MD study argued two plausible sources of entropic stabilization; one is the movements of the ligands inside the cavity, and also the water molecules which can get released from solvation layer of the ligand upon binding. The role of water molecules in conformational stabilization of protein-ligand association is already known.23-25 Although informative, but one of the major limitations of earlier report on the MD simulation of the earlier work by Chakraborti et al. was implicit nature of the solvent. Although implicit model of solvation was suitable for the objective of that study to complement experimental observations, but it was not able to witness the explicit interactions caused by the water molecules and so offered only a gross qualitative insight into the binding mechanism. This present report aims to understand the possible role of the water that fills up the space created by the movements of the TN16 to obtain much insightful understanding on the thermodynamic characteristics of the Tubulin-TN16 binding and consequences in the ligand design.

ACS Paragon Plus Environment

Page 6 of 38

Page 7 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Methods: Modelling and simulations The modelled systems are listed in Table 2. The structural coordinates of the Tubulin TN16 complex (TBE200I1, TBE200I2) was obtained from Protein Data Bank (PDB) using accession code 3HKD, refined at 3.70Å resolution20. Chain A and B along with the TN16, GTP and GDP were used for modelling. We modelled another independent system of Tubulin TN16 complex (TBE200P) where we neutralized the side chain of E200 by adding a hydrogen atom and kept the remaining coordinates of the complex same as crystal structure. The uncomplexed-curved (TBUC) conformation was taken from Protein Data Bank using PDB code 3HKB refined at 3.65 Å resolution20. The uncomplexed-straight (TBUS) was modelled using the PDB code 1TUB refined at 3.7 Å resolution2. Structurally missing loop regions and truncated segments of crystal structures were modelled using Modeller 9.1026. All the systems were processed using CHARMM package

27, 28

implemented in CHARMM-

GUI webserver29 and CHARMM36 parameter30 set was used to represent them. TN16 was parameterized using CHARMM General Force Field (CGenFF)31. Cubic boxes of TIP3P water model32 was used to solvate each system ensuring minimum ~9 Å thickness of surrounding water everywhere. While solvating the system using a box of water molecules, any water oxygen falling within 2.6 Å of any heavy atom of the solute were deleted by residue. Mg2+ and Cl- ions were added to neutralize the charges as well as to maintain 0.15 M ionic concentration. Each system contained on an average ~175000 atoms. All of the systems were minimized using the adopted basis Newton-Raphson method (ABNR). Thereafter (MD) simulations were run using CHARMM36 force field and NAMD 2.9 simulation package33. The systems were heated to 300 K followed by equilibration. Then NPT simulations at 300K were performed using 1fs integration time step. Short ranged non-bonded interactions were

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

truncated at 12 Å whereas long range electrostatic interactions were computed using particle mesh Ewald (PME)34 technique. The pressure and temperature was controlled using Nosé– Hoover thermostat/barostat algorithm.35, 36 For each system the production run was carried out for 65 ns, as summarized in Table 2; the uncomplexed TN16 was sampled only for 30ns. Movies and figures were prepared using VMD37 and PyMOL38. Throughout the manuscript, the term ‘representative snapshot’ has been meant to be the frames which can represent most of the significant interactions in a single frame, chosen from the last few nanoseconds of a trajectory. Quantum Mechanical (QM) calculations A few representative snapshots from the last 10 ns of the trajectories have been considered form quantum mechanical (QM) calculations to analyse the electronic states. For each frame, only a subset of atoms were used, which were as the following: TN16, E200, Y202, V238, C241, and one or more H2O molecule(s) depending on the frame. All quantum mechanical calculations were carried out by Gaussian 09 package39 . Natural Bond Orbital (NBO) method using Hartree Fock theory with 6-31G* (d,p) basis set [HF/6-31G* (d,p)] was applied to calculate charges and orbital interactions40 on the representative snapshots taken from the different trajectories. Only the TN16 and the surrounding residues bridged by water were included in the calculations to make it computationally feasible. NBO second order perturbation analysis and charge calculation40 was carried out for each case.

Results: Dynamics of TN16: Deprotonated E200, Trajectory 1 (TBE200I1)

ACS Paragon Plus Environment

Page 8 of 38

Page 9 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

The initial models did not accommodate any water molecules in the vicinity of the TN16, even after solvation. The space between the α and β subunit of tubulin were filled with the thin film of water (see supporting information Figure S1) which was very close to the cavity that contains the ligand. Although the system got thermally equilibrated within the initial dynamics 50ps but the structural equilibration took a bit longer time, a few nanosecond. It was observed that at the early stages of the production run the water molecules tend to soak into the core of the protein and eventually two water molecules succeeded to penetrate into the core of the protein by slightly opening of the space between P245, G246 (part of an unstructured loop of tubulin chain B) and C356, D357 (terminal of a β-strand of chain B) shown in Figure 3. One water molecule quickly moved towards the TN16 and bridged the interaction between protein and TN16 (SI movie M1). The other water molecule although initially moved within the binding pocket rapidly, but within less than a nanosecond it gets stabilized by bridging the interaction between protein and TN16. Figure 2c has summarized the water molecules getting trapped between ligand and protein. The most striking event is that the movement of the second water molecules gets concerted with the change of the C-CC-C dihedral angle change of TN16 causing a flip or twist of the ligand as shown in the Figure 2c and SI movie M1. Once it flips, the O atom of the ligand, which was originally close to E200 side chain, moves to the opposite side swapping its position with ‘-NH’ group of the same ring (5 membered, Figure 4) of TN16 (SI movie M1). At this point the two water molecules bridges the ligand with protein from two sides, one forming the interaction between E200 side chain and ‘-NH’ of TN16, another in between the C241 backbone NH and ‘O’ of TN16 (Figure 2c). Throughout the trajectory it has been observed that although 4-5 water molecules visit the binding pocket, but only two of them are of ‘bridging’ type. The flipping of partially polar five membered ring of TN16 occurs at ~2.8 ns which takes a rotation of ~150°-160° and plausibly this event facilitates the entry of other water molecules

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

within the binding pocket (SI movie M1) by creating a space between D199, E200 (part of a β-sheet of chain B) and A256, V257 (part of H8 helix of chain B) shown in Figure 3. The distribution of sampled conformations Figure S2 of complexed and uncomplexed TN16 over the trajectories has revealed that ligand conformations are predominantly confined into one major cluster, but a small but finite population is witnessed in uncompleted state which is similar to the conformations observed in TBE200I1. TN16 experiences two timescales of Bridging-Hbond on its two side; one side the water bridging to V238 and C241 is long lasting as well and the same water molecule is retained throughout, though the Hbond is broken and reformed many times. On the other hand the water molecule bridging the NH of TN16 and E200 side chain is observed to be continuously replaced by another molecule of water but the ‘bridging’ is preserved throughout. The lifetime of these interactions were calculated and discussed in the following section. After the structures attained equilibrium within ~5 ns, the TN16 was observed to be surrounded by the residues which were mostly hydrophobic in nature (Figure 5). Two apolar rings of TN16 remained packed within the vicinity of A354, V318, A316, I378, M259, and I4, L252, V238, L242, L255, F169 etc. of β-tubulin. These hydrophobic networking between the TN16 and protein had an impact on the TN16 binding stability (enthalpic gain) inside the colchicine binding pocket which clue could be extracted from the previous report by Chakraborti et al.22 . Here the simulations have confirmed that the five membered ring of TN16 is comfortable with at least

two different ways of interactions within the binding

pocket. Before the flipping of the TN16 ligand, its two polar groups (-NH and OAB of five member ring, see Figure 4) were positioned in the vicinity of C241 and V238 of B-chain and was coordinated with the protein through a single water (say W1); but after the rotation of the ligand, the space emptied by the -NH group of five member ring got filled up a new incoming

ACS Paragon Plus Environment

Page 10 of 38

Page 11 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

water (say W2) and the water which was already in this zone (W1) moved to the other side of the ligand in a correlated manner with the ligand’s rotation (SI movie M1). W1 now forms new polar interactions between E200, N167 and TN16 and later associates two more incoming water molecules in this coordination (Figure 2c). The W2 remained static for the rest of the time in that position coordinating C241, V238 and CO of TN16.

Thus a

maximization of the polar interactions takes place. Deprotonated E200, Trajectory 2 (TBE200I2) In an independent run of the same system mentioned above, the ‘flip’ of TN16 was not witnessed within 65 ns, and associated the bridging water in a fashion which is much different from TBE200I1 (Figure 2c, 2d). In this trajectory the water molecules also entered into the binding pocket in the early stage of production dynamics (~2.5 ns-3 ns) and retained the two types of bridging interactions with E200 and ‘O’ atom of TN16 as well as with backbone NH of C241 and ‘OAB’ of TN16 as found earlier (Figure 2d). A new interaction is formed between backbone ‘O’ of G237 and ‘OAB’ of TN16 through bridging water. Interestingly these two trajectories independently confirm the path through which the water molecules enter into the binding pocket (Figure 3). Protonated E200 (TBE200P) In the trajectory where E200 side chain was protonated and neutral, the primary difference that was observed from the protonated set was the stability of the conformation of TN16, as it was observed in the crystal structure, i.e. no twist or flip was witnessed. As a consequence, the -COOH group, which is now capable to act as H-bond donor, form H-bond with ‘O’ atom of the TN16 directly (Figure 6a), without the need of any mediator H2O. So the bridging water which finally stabilized only at one site instead of two (as seen for deprotonated cases) as shown in Figure 2b,c,d and Figure 6a,b. In the crystal structure the following polar

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

interactions are observed: (i) between E200 side chain and ‘O’ of TN16 (ii) backbone CO of V238 and ‘–HN’ of TN16 (Figure 2a). Simulation shows that these interactions are preserved and two new interaction is formed: (i) between ‘OAB’ of TN16 (see Figure 4) and backbone NH of C241 (ii) between backbone ‘O’ of G237 and ‘OAB’ of TN16, which are mediated by the water molecule as shown in Figure 2b. Since the rest of the interactions are hydrophobic and provide the ligand an environment analogous to hydrophobic enclosure, the gross features of the interactions remain similar to that of deprotonated state (Figure 5). A ligplot41 shows the surrounding environment in supporting information Figure S3, S4, S5 and S6 for different systems. Water dynamics in the ligand binding pocket: H-Bonds: In Uncomplexed Curved Tubulin (TBUC) and Straight Tubulin (TBUS) The statistics of water occupancy in the binding cavity of the uncomplexed tubulin dimer (TBUC) shows that on an average ~14 water molecules at a time are occupied the cavity, as mentioned in the Table S1 of supporting information. The average lifetime of the proteinwater Hbond was in the order of 10ps for most of the water molecules, except one molecule which was ~126 ps. In case of uncomplexed straight tubulin dimer (TBUS), the average number of water was ~9 because the T7-loop pre-occupied the binding pocket. The H-bond lifetime was in the order of 10 ns for all the water molecules in the cavity except one with ~121 ps. H-bonds: in the Complexes The lifespan of the Hbonds formed by the water molecules bridging the tubulin and TN16 have been measured which gives an estimation of their stabilities compared to bulk water. In TBE200I1, one water molecule bridging between V238 and TN16 breaks only after 278.6 ps

ACS Paragon Plus Environment

Page 12 of 38

Page 13 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

(see Table 3), although the Hbond between the same set of atoms reformed; thus experienced making and breaking of Hbonds for multiple times within 65 ns timescale (Table S1). The water bridging the C241 and TN16 had average lifetime of ~27.8 ps. On the other side, the Hbond between T240 and TN16 had average lifetime ~9 ps which is in the range of several other H-bonds listed in Table 3. In TBE200I2 the longest lifespan of such bridging water molecules were observed to be ~33.5 ps. It seems from this Hbond life time that the flipping of the TN16 in TBE200I1 places two groups appropriate fashion suitable for Hbonds stability, which is absent in TBE200I2. The protonation of E200 sidechain reduces the number of sites occupied by the bridging water molecule (Table 3) anytime during 65 ns; it was maximum in TBE200I1 (total 8). So as a consequence of the flip of the TN16 the number of available site as well as the lifespan of the Hbonds increase. In a following section its consequence on the binding free energy has been discussed. The water molecule adjacent to the sight of G237, V238, T240 and C241 is mainly involved in the long lasting bridging interaction between the ligand and protein. Table 3 presents the list of sites and the corresponding H-bond life time of the bridged water molecules for all the three complexed systems (TBE200I1, TBE200I2 and TBE200P). The H-bond life time was calculated for the bridging water molecules as well as the ‘non-bridging’. Here the term ‘nonbridging refers to those for those water molecules which are not involved in bridge formation but are able to make polar interactions by forming H-bonds with the binding pocket residues over the last 60 ns of the trajectory. During the calculation of H-bond life time, criterion for distance cut-off was set to 2.4 Å and the dihedral angle cut off was kept to 120°. Rotational dynamics: The rotational relaxation of the water molecules, i.e. their reorientational dynamics, were analysed using the following autocorrelation function Γ242-44:

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Γ2 (t ) = P2 (µ (0).µ (t ) )

Page 14 of 38

(1)

Where µ (t) is the unit vector along the water dipole at time t and the angular brackets indicate a time average. P2 is the second-order Legendre polynomials:

P2 ( x) =

(3 x 2 − 1) 2

(2)

The rotational correlation was calculated separately for the water molecule present in four different layers defined as (i) bridging, (ii) interfacial, (iii) hydration and (iv) bulk. The ‘bridging’ has been already defined in the preceding section. The interfacial water is the thin film of water filling up the space between α and β subunits, satisfying the 4Å cut-off distance from both subunits. Hydration was defined as the water molecules within 4Å of tubulin surface excluding the type (i) and (ii). The bulk was mimicked by a box of water simulated for 2 ns separately. As shown in Figure 7, the rotational correlations as a function of time falls rapidly for the bulk water and relatively slower for the case of surface water (hydration water). The ‘glassy’ character of the first hydration water compared to the ‘bulk’ is already known from past studies45,

46

. The noticeable character is the much slower rotational

properties of the interfacial water compared to the hydration water. This is because the interfacial water molecules have the chance to form more than one Hbonds with the two protein surfaces (α and β) at two sides; if not forming Hbonds, their movements still get perturbed by the long range interactions from two sides. The ‘bridging’ water have been observed to be the slowest in the rotational relaxation and its difference with interfacial water is substantial. This is intriguing because it implies that a water molecules gets more stabilization when it is placed between tubulin and TN16 compared to two protein surfaces. Its effect on the overall protein-ligand complex stabilization is reported the following section.

Binding energy calculation using MM-GBSA protocol:

ACS Paragon Plus Environment

Page 15 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

The binding energy was computed using MMGBSA protocol47-49 . To compute the free energy of solvation the GBMV50, 51 method implemented in CHARMM was used. In the first round of calculation, all the water molecules were removed from the trajectory. Then the trajectories of the complex, were split into ligand and receptor to compute the binding energies. The cumulative average of the binding energies (∆E) show (Figure 8) that for all the systems a reasonable convergence of the binding energies is attained within 20ns, although they not necessarily agree with each other by then. The almost horizontal trend of the ∆E of all the systems are observed during 53-65 ns regime, indicating a much tighter convergence. The averaged energies have been listed in the Table 4. The deprotonated systems, TBE200I1 and TBE200I2 have binding energies -24.8 kcal/mol and -20.8 kcal/mol respectively whereas the TBE200P gives the binding energy -33.0 kcal/mol immediately advocating the highest stability of that conformation. The deprotonated situation causes a repulsion with the ligand which is never recovered even after the ligand is flipped to adjust its orientation. Also the other noticeable issue is that the time requires to reach the convergence of the ensemble averaged energy are almost similar in both the trajectories (TBE200I1 and TBE200I2) and independent of the ligand flipping. The second round of the entire calculation was repeated by considering the bridging water molecules explicitly in the MMGBSA binding energy.52 As shown in Table 5, the three systems TBE200I1, TBE200I2 and TBE200P have yielded -44.5 kcal/mol, -44.0 kcal/mol and -41.8 kcal/mol respectively. So it is the explicit presence of the bridging water molecules in the binding energy calculations which reflects a substantial difference in the binding energy.

QM calculations:

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The elementary understanding of chemistry would argue that the tightly packed water molecules would be differently polarized which many not agree with the partial charges on its atoms governed by the TIP3P model of water. Since the polarizable force fields are not yet readily available for MD simulation, to understand the extent of polarization and the alteration in the electronic levels involved in this complexation the NBO second order perturbation analysis and charge calculation40 were carried out on a few representative snapshots taken from the TBE200I1 trajectory. The snapshot contained two water molecules at a time, each of which was considered in the QM in separate runs, say QMW1, QMW2. In QMW1, where the water molecule bridges the TN16 with C241 and V238 (similar to the situation shown Figure 2c), the calculations have revealed that the carbonyl moiety of TN16 holds this water involving its lone pair on oxygen as donor and O-H σ* bond of water as acceptor ( En→σ ∗ = 7.8 kcal/mol) at O…H…O angle 153.3°. Similar type of interaction has been perceived involving lone pair on oxygen of acid group of V238 as donor and σ O∗ −H bond of water as acceptor (For lone pair 1: En→σ ∗ = 6.30 kcal/mol, for loan pair 2: En→σ ∗ = 7.2 kcal/mol) at O…H….O angle 173.1°; both have been summarized in Figure 9a. While above stated n → σ ∗ from TN16 to water is active, another strong n → σ ∗ interaction is arising from lone pair on oxygen of that interacting water as donor and σ N∗ − H of C241 as acceptor ( En→σ ∗ = 16.3 kcal/mol), as shown in Figure 9b. As a result partial charge on N atom of amino group of C241 has decreased by 0.04e with respect to the same in absence of water molecule. It also indicates that above type of interaction has increased the electron density on N atom of amino group of C241. Two simultaneous electron delocalization events create a non-bonded bridging interaction in which water plays both acceptor and donor. It could be proposed that this static water has some pharmacological importance as they can determine or discriminate between the different binding modes of TN16.

ACS Paragon Plus Environment

Page 16 of 38

Page 17 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

In QMW2, the bridging water molecule is involved in n → σ ∗ interaction (Figure 10a) with surrounding amino acid and TN16. Carbonyl moiety of pyrrolidine of TN16 interacts with this water molecule involving its two lone pairs as donor and σ O∗ − H of water as acceptor ( En→σ ∗ * = 2.3 kcal/mol for lone pair 1, En→σ ∗ = 2.70 kcal/mol for lone pair 2) at the geometry having O…H….O angle is 159.5°. Here the water molecule is involved in above type interaction arising from two lone pairs on oxygen of the E200 side chain as donor and σ O∗ − H of water as acceptor ( En→σ ∗ = 9.1 kcal/mol for lone pair 1, En→σ ∗ = 8.0 kcal/mol for lone pair 2); at the geometry having O…H…O angle 164.1°. Here we can propose that water molecule is clasped by carbonyl moiety of pyrrolidine of TN16 and side chain acid group of E200 due to non-bonded interaction. From NBO analysis of the representative snapshot from TBE200P it has been observed that carbonyl moiety of pyrrolidine of TN16 interacts with one water involving the two lone pairs on the carbonyl oxygen as donor and σ O∗ − H as acceptor ( En→σ ∗ = 9.5 kcal/mol for lone pair 1,

En→σ ∗ = 1.5 kcal/mol for lone pair 2); at the geometry having O…H….O angle is 143.7°. Simultaneously lone pair on oxygen of water donates electron pair to the σ N∗ − H of cystine ( En→σ ∗ = 11.6 kcal/mol). Besides, σ O − H of water weakly interacts (Figure 10b) with σ N∗ −H of C241

( Eσ →σ ∗ = 0.7 kcal/mol). So water acts like a bridge in this electron shift relay

between TN16 and C241. This gives stability to the ligand binding mode. Although this effect of polarization is not reflected in the CHARMM force field driven dynamics in the present study, but it could be predicted that in the real system the molecular complex would gain form this additional benefit.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Discussion: As mentioned in the introduction, set of the ligands which can occupy the colchicine binding site vary widely in structures. The experimentally determined structures of some of them show (Figure 1) that the mode of binding and the interaction sites also differ from ligand to ligand, indicating the presence of redundant opportunities of interactions; where each ligand are able to exploit only a subset of it. As mentioned in the introduction (and in the article20 reporting 3HKD) the TN16 was placed in the PDB model by optimizing several parameters, i.e. the positions of the ligands had not been possible to be resolved with desired resolution. The simulations have revealed at least two possible distinct orientations, which vary as a function of two other factors, the first one is the side chain ionization state of E200 of chain B and the other one is the number of the explicit water molecules which directly take part in the binding of tubulin and TN16. Since the interactions occur in a complex in many body fashion, the dependence on each factor cannot be distilled out into independent components. The TBE200I1 trajectory started with the crystallographic orientation of TN16 and the ionized side chains of E200. The ligand quickly re-adjusted the repulsion between the E200 and the partially negatively charged ‘O’ of the carbonyl group of TN16 by swapping the interactions by a flip placing the H-N group near the COO- and moving the CO on the opposite side. As observed from the simulation, the swapping of interactions takes place only after the two water molecules enter into the cavity, by breaking some weaker interactions on the surface, as mentioned in the result section. The stabilization of the ligand’s orientation and the coordination of the interactions by the water molecules take place almost simultaneously. In the trajectory where the E200 side chain charge has been neutralized by protonation (TBE200P), does not show any flip of the ligand. The penetration of the water molecules through the surface pore was observed to be same as TBE200I1 and TBE200I2, but the difference was in the number of water molecules that form the bridging interaction;

ACS Paragon Plus Environment

Page 18 of 38

Page 19 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

TBE200P had 1 water molecule whereas the rest had 2. These waters plausibly are also not detected in the experiments due to their positional exchanges due to ligand’s fluctuations and the relatively lower resolution (3.7 Å) at which the structure was resolved. An independent run with ionized E200 side chain (TBE200I2) did not show any change in the ligand’s orientations, and also associate lesser number of bridging water, which has not been found to be energetically most preferred situation. The next issue that was investigated was the extent of damping the rotational and transitional motions of the waters around the ligand and/or within the binding pocket. For comparison, the water molecules in different layers around the protein was also investigated. The rotational correlation of the water dipole decays rapidly for the bulk water, from 1.0 goes below 0.2 within ~4.5 ps and decreases to 0.05 in 10 ps. Whereas the same for the water molecules in the first hydration layer (within 4 Å on the surface) of the protein for all the systems decays from 1 to 0.2 in 10ps, i.e. significantly slower than the bulk water. The water at the α, β interface of the dimer is even slower than the surface water, as they are sandwiched between the two domains of proteins (in contrary to surface, where one side of the water is the ‘bulk’) and have much stronger perturbation on the movements through Hbond and other interactions. The bridging water molecules have been slowest in rotation; in fact they have been almost restricted from rotating at all, as the correlation doesn’t go below 0.5 in 10 ps (Figure 7). This is in agreement with the fact that these water molecules brings extra stabilization to the protein-ligand interactions. The strength of the H-bond depends on the H-bond angle as well as the donor acceptor distance. So although the bridging water are rotationally frozen, whether the distance criteria is maintained or not to provide stability is a critical information about the system. The H-bond lifetime was measured using the dihedral angle cut-off to 120° and 2.4 Å distance cut-off, which revealed the long lasting nature of the bridging H-bonds particularly in TBE200I1

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

So it is clear that without having the insight of this water mediated interactions the understanding on the mechanism of TN16 binding into tubulin remains substantially incomplete. The dynamics of these water molecules are also significant to compare the alteration of the micro details of the complexed and uncomplexed situations. The ensemble of the uncomplexed states show that the ligand binding cavity is filled up with at least 9-14 water molecules (varying due to the T7 orientations in straight and curved conformations; Figure S7); their average lifetime of H-bond is either in tens of picosecond or in the range of 120 ps to185 ps. Where in TBE200I1, the number of water in the cavity become ~4 and at least one H-bond became highly stable (~744 ps, see Table S1), in the order of nanosecond. So the ligand binding leads to the removal of at least 4-9 water molecules shifting the H-bond life time in the nanosecond regime. The removal of water shall add to the entropic gain, which is in agreement with the experiment (Table 1) and the mechanism hypothesized in the previous report with implicit water simulation22; but the increase of lifetime of the H-bond would certainly diminishes the entropic benefit but the experimental data evidences that the resultant effect is still have the entropic favour (possibly due to overwhelming number of released water compared to the number of water that become frozen upon the ligand binding). The computed binding energies using MMGBSA protocol rationalizes the structural observations and the dynamic events witnessed. As presented in the results, the two systems, protonated (TBE200P) and deprotonated E200 (TBE200I1 and TBE200I2), differ widely in the binding energy when the solvation energy of the systems are computed implicitly. In a straightforward MMGBSA calculation the TBE200P appears to be the most stable conformation of the complex. But a careful look at the cumulative average (Figure 8) shows that there is a gradual increase of the binding energy (i.e. decrease in stability) in case of TBE200I2 and TBE200P. This is in contrary to the understanding that as the simulations progress the complex should adopt to a more stable conformations by maximizing the

ACS Paragon Plus Environment

Page 20 of 38

Page 21 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

interactions at the binning interface. But the rise of binding energy for a stable system is intriguing (for clarity, the first 5 ns of the simulations have been shown in Figure S8) and suggests that possibly some component, that might be strongly responsible for the stability, is not being considered in the total binding energy with appropriate weightage. Apparently, the only structural variation with time has been filling up the space in the cavity by bridging water and so plausibly their tight packing into the cavity requires more sophisticated approach (compared to implicit solvation) to include the effect. Indeed, the re-calculation of the binding energy dramatically changes the binding energies, and for both systems it shows a growth of stabilization, which is generally expected from such sampling techniques at room temperature. For the system TBE200I1, there is a rapid decrease of binding energy (i.e. increased stabilization) which agrees with the events of ligand flip and the capture of two bridging water in the binding. The binding energy obtained from TBE200I2 was -44.0 kcal/mol which is almost equal to the same for TBE200I1 (-44.5 kcal/mol), which reflects that the structural differences introduced by the ligand flip is compensated by the variation in the bridging water association. The binding energy from the TBE200P was -41.8 kcal/mol which is slightly weaker and doesn’t reflect any benefit of protonating the E200 side chain to avoid repulsion. These three trajectories yield a range of binding energies which exist within a band of ~3 kcal/mol. Given the degree of accuracy of standard force fields in general and approximations involved in the MMGBSA methods one can argue that the ~3 kcal/mol range of window is not a sufficiently broad range of energy to confidently claim that these three trajectories are energetically well separated. Instead, it gives a clue that these ensembles of conformations energetically overalp and thus might experience interconversion between conformations under thermal fluctuations. So the calculations clearly revel that the tubulin can accommodate TN16 in multiple ways and TN16 can have at least two distinct orientations; all being interconvertible without much expense of energy. This is surely

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

contributing to the increase of conformational heterogeneity of the ensemble leading to entropic benefit of the complex, which is in full agreement with the experimentally observed trend (Table 1). Investigations on such events of molecular level would require truly polarizable force field along with appropriate tuning with respect to thermodynamic data. Such a force field is not yet available with suitable benchmark for regular use. So to investigate whether the bridging water molecules are really inducing polarization on the ligands through the strong H-bond, quantum mechanical calculations were carried out. The NBO calculations show a concerted interactions between the protein and the ligand via water molecules, as n → σ ∗ (V238 to water) and then n → σ ∗ (water to TN16) and n → σ ∗ (TN16 to water) and then n → σ ∗ (water to C241), giving substantial stabilization to the system. Although, such stabilization is not extractable from the MMGBSA methods (using standard force field) but it gives an indication such extra sources of energies can further contribute to the conformational stabilizations as well as conformational interconversions. Several other interactions of the similar type have been noticed. So it is clear that the fewer water molecules, which becomes static to be a constitutional unit of the interactions between protein and ligand, experience molecular-orbital overlaps to bring extra stabilization and therefore they have a much stronger pharmacological importance.

Conclusions: The binding mechanism of the TN16 inside in the colchicine binding site of the tubulin has been investigated using modelling and simulations. The topology of TN16-tubulin interactions provided in the PDB model (3HKD) was not fully convincing and the same concern was expressed by Dorleans et al.1 expressing their opinion about the possibility of other possible interactions. The present report provides several major insights into the mechanism to overcome the limitations of the understanding of past: (i) the ligand can

ACS Paragon Plus Environment

Page 22 of 38

Page 23 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

experience multiple ways of binding by rotating itself within the deeply buried shallow binding cavity which feature cannot be captured in a single model (ii) the resolution (3.7 Å) of 3HKD was much inferior to report any structural water, whereas the fluctuations of the ligand is highly coordinated with the intake of ‘bridging water’ molecules (iii) the multiple modes of binding lie on with a short band of ~3 kcal/mol and can experience interconversion, reflecting a quasi-degenerate situation (iv) QM calculations a relay of electron density between TN16 and protein via structural water, making it an integral part of stabilization (v) since there is no experimental benchmark available to compare the confidence level of the binding energy calculation with and without explicit water, a ‘reliability check’ strategy (for instance, the cumulative average adopted for the present report) could be formulated by comparing the structure and stabilization as a function of time, which could be tested on other systems in general. The insight into the quasi-degeneracy of the conformational states of the TN16 and the structural water appear to be have very high pharmacologically significance.

Acknowledgement:

SGD thanks B Bhattacharyya of Bose Institute for insightful

discussion.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 38

Tables:

Table 1: Thermodynamic data on the Tubulin-TN16 binding, as available from literature. The data for Colchicine has also been shown as the binding site is commonly referred as the ‘colchicine binding site’.

Thermodynamic Parameters

TN1622

Colchicine53, 54

∆H(Kcal/mol)

-1.7 ± 0.58

10.0

T∆S(Kcal/mol)

5.7

19.6

∆G(Kcal/mol)

-7.4

-9.6

Table 2: List of trajectories Serial No.

Trajectory name

Description of the systems

1

TBE200I1

Tubulin E200 ionized_TN16 complex (curved, HKD)- Trajectory 1

65

2

TBE200I2

Tubulin E200 ionized_TN16 complex (curved, HKD)- Trajectory 2

65

3

TBE200P

Tubulin_TN16 complex (curved, HKD)

65

4

TBUC

Tubulin Uncomplexed (curved, 3HKB)

65

5

TBUS

Tubulin Uncomplexed (straight, 1TUB)

65

6

TN16U

TN16 Uncomplexed

30

Total

ACS Paragon Plus Environment

Length (ns)

355

Page 25 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Table 3: Average lifetime of the protein-water-ligand bridge. Two alternative ways of forming the H-bond and their lifetime are separated by front slash ‘/’. The lifetime does not necessarily mean the replacement of that particular water molecule, as the same water molecule can either reform the bridge or can get replaced by an incoming water, as elaborated in the text. System name

TBE200I1

TBE200I2

TBE200P

Amino acid resid and atom name

TN16-atom name

Protein-ligandwater bridging Hbond life-time(ps)

V238-O

N16-O

278.6

C241-HN/HG1

N16-O

27.8/7

E200-OE1

N16-OAB/HN

10.4/8.9

E200-OE2

N16-OAB/HN

11.0/8.5

T240-HG1/HN

N16-O

9.0/5.5

Y202-OH

N16-OAB/HN

5.7/5.7

Y202-HH

N16-HN

7.5

N167-OD1/HD21

N16-HN

7.9/5.0

C241-HN/HG1

N16-OAB

33.5/8.1

G237-O

N16-OAB

22.1

Y202-OH/HH

N16-O

8.1/7.4

T240-HN/OG1/HG1

N16-OAB

5.0/6.0/5.2

E200-OE1/OE2

N16-O

12.4/13.0

V238-O

N16-OAB

5.5

C241-HN/HG1

N16-OAB

54.0/6.7

G237-O

N16-OAB

17.3

V238-O

N16-OAB

6.6

T240-HN/HG1

N16-OAB

6.2/5.0

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 38

Table 4: Components of binding energy (in kcal/mol) in different trajectories representing the solvation only by implicit model of water. Standard deviations are indicated within the parentheses.

Components of energy

TBE200I1

TBE200I2

TBE200P

∆Eelec

-36.2(5.0)

-13.8(2.9)

-28.3(3.5)

∆Evdw

-42.0(2.4)

-38.7(2.1)

-40.7(2.5)

∆Egb(polar)

57.7(8.5)

35.9(6.7)

40.0(5.8)

∆Egb(nonpolar)

-4.0(0.2)

-4.2(0.2)

-3.9(0.2)

∆Egb(polar)+elec

21.5(8.5)

22.1(6.6)

11.7(6.1)

-24.4(8.7)

-20.8(6.7)

-33.0(6.0)

∆Etotal

Table 5: Components of binding energy (in kcal/mol) in different trajectories representing the solvation by implicit model of water and including the structural water explicitly, as mentioned in the text. Standard deviations are indicated within the parentheses. Components of energy

TBE200I1

TBE200I2

TBE200P

∆Eelec

-50.2(5.3)

-34.4(5.8)

-33.6(4.2)

∆Evdw

-40.8(3.0)

-39.4(2.8)

-41.9(2.8)

∆Egb(polar)

50.4(9.6)

33.6(9.5)

37.5(6.6)

∆Egb(nonpolar)

-3.9(0.2)

-3.8(0.2)

-3.7(0.2)

∆Egb(polar)+elec

0.2(10.1)

-0.7(10.5)

3.9(7.1)

∆Etotal

-44.5(10)

-44.0(10.4)

-41.8(7.1)

4

9

6

Average no. of water

ACS Paragon Plus Environment

Page 27 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Figures:

Figure 1: Four examples of different ligands (magenta) occupying the colchicine binding site of tubulin revealing the diversity of the interactions and heterogeneity of the ligands’ orientations, as observed from their crystal structures (a) 1SA0, (b) 3HKC, (c) 3HKD, and (d) 3HKE. The bond line formula of the ligands are shown at the bottom of each structure. The protein has been shown in (green) cartoon representation. A few sidechains are shown in sticks.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2: (a) Polar interactions witnessed from the crystal structure (PDBID: 3HKD) of the Tubulin-TN16 complex. Representative snapshots showing the water mediated interactions in the following trajectories: (b) TBE200P (c) TBE200I1 (d) TBE200I2. The colouring scheme is same as Figure 1.

ACS Paragon Plus Environment

Page 28 of 38

Page 29 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Figure 3: A schematic representation of the path for the entry of water molecules into the TN16 binding pocket. The cyan stick refers the TN16 and the protein has been shown in green cartoon along with important side chains are show in sticks. The ball and stick model has been used or the water.

Figure 4: The TN16 molecule. The labels (atom names) have been shown for convenience as some of them have been referred in the text.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 5: Amino acid residues within the 4Å of TN16: (a) crystal structure (PDBID: 3HKD) (b) a representative snapshot taken around 30 ns of TBE200I1.

Figure 6: Representative snapshots from the trajectories having (a) neutral E200 sidechain and (b) ionized side chain of E200. The difference in the orientations of the 5 membered ring is noticeable, which reveals two distinctly different topologies of interactions of TN16 with the protein and they also accommodate different numbers of bridging water molecules.

ACS Paragon Plus Environment

Page 30 of 38

Page 31 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Figure 7: Rotational relaxation of water molecules in different layers around protein and TN16

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 8: Cumulative average of the binding energy (∆E) calculated using MMGBSA method, as a function of simulation time, for the systems with E200 protonated (TBE200P, green lines) and ionized (TBE200I1 and TBE200I2, red and blue lines respectively). The ‘IMPONLY’ refers the standard MMGBSA calculations without any explicit water, which have been shown with dashed lines. The ‘IMP+EXP’ refers to the calculations involving structural water molecules (elaborated in the text) which have been shown with smooth lines. Same colour for the smooth and dashed lines indicates same trajectory but two different calculations. A vertical line has been shown at 53ns to highlight better convergence of data in the 53-65ns regime.

ACS Paragon Plus Environment

Page 32 of 38

Page 33 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

Figure 9: Qualitative representative of the orbital interactions revealed by the NBO calculations: (a) electron transitions ( n → σ ∗ ) between TN16-H2Oand between V238-H2O.

(b) Concurrent electron transition ( n → σ ∗ ) between TN16-H2O and H2O-C241.

Figure 10: Qualitative representative of the orbital interactions revealed by the NBO calculations: (a) Concurrent electron transition ( n → σ ∗ ) between E200-H2O and TN16-H2O. (b) Electron transition ( σ → σ ∗ ) between water and C241.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Supporting Information Supplementary text and figures giving further details on the structure and dynamics, HBonds life-time, statistics of the TN16 conformations, ligplots showing proteinwater-ligand interactions. Movie showing water entry and ligand’s dynamics.

References: 1. 2. 3. 4.

5. 6. 7. 8. 9.

10. 11. 12. 13.

14.

15.

16.

Jordan, M. A., and Wilson, L. (2004) Microtubules as a target for anticancer drugs, Nature reviews. Cancer 4, 253-265. Nogales, E., Wolf, S. G., and Downing, K. H. (1998) Structure of the alpha beta tubulin dimer by electron crystallography, Nature 391, 199-203. Verhey, K. J., and Gaertig, J. (2007) The tubulin code, Cell cycle 6, 2152-2160. Alberts, B., Johnson, A.,Lewis, J., et al. (2002) Molecular Biology of the Cell: The SelfAssembly and Dynamic Structure of Cytoskeletal Filaments, 4th ed., pp 908-927, Garland Science. van der Vaart, B., Akhmanova, A., and Straube, A. (2009) Regulation of microtubule dynamic instability, Biochemical Society transactions 37, 1007-1013. Janson, M. E., de Dood, M. E., and Dogterom, M. (2003) Dynamic instability of microtubules is regulated by force, The Journal of cell biology 161, 1029-1034. Mitchison, T., and Kirschner, M. (1984) Dynamic instability of microtubule growth, Nature 312, 237-242. Stanton, R. A., Gernert, K. M., Nettles, J. H., and Aneja, R. (2011) Drugs that target dynamic microtubules: a new molecular perspective, Medicinal research reviews 31, 443-481. Gupta, K., and Panda, D. (2002) Perturbation of microtubule polymerization by quercetin through tubulin binding: a novel mechanism of its antiproliferative activity, Biochemistry 41, 13029-13038. Chen, S. M., Meng, L. H., and Ding, J. (2010) New microtubule-inhibiting anticancer agents, Expert opinion on investigational drugs 19, 329-343. Dumontet, C., and Jordan, M. A. (2010) Microtubule-binding agents: a dynamic field of cancer therapeutics, Nature reviews. Drug discovery 9, 790-803. Guha, S., and Bhattacharyya, B. (1997) The colchicine-tubulin interaction: A review, Curr Sci India 73, 351-358. Begaye, A., Trostel, S., Zhao, Z., Taylor, R. E., Schriemer, D. C., and Sackett, D. L. (2011) Mutations in the beta-tubulin binding site for peloruside A confer resistance by targeting a cleft significant in side chain binding, Cell cycle 10, 3387-3396. Mooberry, S. L., Taoka, C. R., and Busquets, L. (1996) Cryptophycin 1 binds to tubulin at a site distinct from the colchicine binding site and at a site that may overlap the vinca binding site, Cancer letters 107, 53-57. Field, J. J., Waight, A. B., and Senter, P. D. (2014) A previously undescribed tubulin binder, Proceedings of the National Academy of Sciences of the United States of America 111, 13684-13685. Kanakkanthara, A., Wilmes, A., O'Brate, A., Escuin, D., Chan, A., Gjyrezi, A., Crawford, J., Rawson, P., Kivell, B., Northcote, P. T., Hamel, E., Giannakakou, P., and Miller, J. H. (2011) Peloruside- and laulimalide-resistant human ovarian carcinoma cells have betaI-tubulin

ACS Paragon Plus Environment

Page 34 of 38

Page 35 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

17.

18.

19.

20.

21.

22.

23.

24. 25. 26.

27.

28.

29. 30.

mutations and altered expression of betaII- and betaIII-tubulin isotypes, Molecular cancer therapeutics 10, 1419-1429. Prota, A. E., Bargsten, K., Diaz, J. F., Marsh, M., Cuevas, C., Liniger, M., Neuhaus, C., Andreu, J. M., Altmann, K. H., and Steinmetz, M. O. (2014) A new tubulin-binding site and pharmacophore for microtubule-destabilizing anticancer drugs, Proceedings of the National Academy of Sciences of the United States of America 111, 13817-13821. Ravelli, R. B. G., Gigant, B., Curmi, P. A., Jourdain, I., Lachkar, S., Sobel, A., and Knossow, M. (2004) Insight into tubulin regulation from a complex with colchicine and a stathmin-like domain, Nature 428, 198-202. Barbier, P., Dorleans, A., Devred, F., Sanz, L., Allegro, D., Alfonso, C., Knossow, M., Peyrot, V., and Andreu, J. M. (2010) Stathmin and interfacial microtubule inhibitors recognize a naturally curved conformation of tubulin dimers, The Journal of biological chemistry 285, 31672-31681. Dorleans, A., Gigant, B., Ravelli, R. B., Mailliet, P., Mikol, V., and Knossow, M. (2009) Variations in the colchicine-binding domain provide insight into the structural switch of tubulin, Proceedings of the National Academy of Sciences of the United States of America 106, 13775-13779. Bhattacharyya, B., Panda, D., Gupta, S., and Banerjee, M. (2008) Anti-mitotic activity of colchicine and the structural basis for its interaction with tubulin, Medicinal research reviews 28, 155-183. Chakraborti, S., Chakravarty, D., Gupta, S., Chatterji, B. P., Dhar, G., Poddar, A., Panda, D., Chakrabarti, P., Ghosh Dastidar, S., and Bhattacharyya, B. (2012) Discrimination of ligands with different flexibilities resulting from the plasticity of the binding site in tubulin, Biochemistry 51, 7138-7148. Bhat, T. N., Bentley, G. A., Boulot, G., Greene, M. I., Tello, D., Dall'Acqua, W., Souchon, H., Schwarz, F. P., Mariuzza, R. A., and Poljak, R. J. (1994) Bound water molecules and conformational stabilization help mediate an antigen-antibody association, Proceedings of the National Academy of Sciences of the United States of America 91, 1089-1093. Garcia, A. E., and Hummer, G. (2000) Water penetration and escape in proteins, Proteins 38, 261-272. Ladbury, J. E. (1996) Just add water! The effect of water on the specificity of protein-ligand binding sites and its potential application to drug design, Chemistry & biology 3, 973-980. Eswar, N., Webb, B., Marti-Renom, M. A., Madhusudhan, M. S., Eramian, D., Shen, M. Y., Pieper, U., and Sali, A. (2007) Comparative protein structure modeling using MODELLER, Current protocols in protein science / editorial board, John E. Coligan ... [et al.] Chapter 2, Unit 2 9. Brooks, B. R., Brooks, C. L., Mackerell, A. D., Nilsson, L., Petrella, R. J., Roux, B., Won, Y., Archontis, G., Bartels, C., Boresch, S., Caflisch, A., Caves, L., Cui, Q., Dinner, A. R., Feig, M., Fischer, S., Gao, J., Hodoscek, M., Im, W., Kuczera, K., Lazaridis, T., Ma, J., Ovchinnikov, V., Paci, E., Pastor, R. W., Post, C. B., Pu, J. Z., Schaefer, M., Tidor, B., Venable, R. M., Woodcock, H. L., Wu, X., Yang, W., York, D. M., and Karplus, M. (2009) CHARMM: The Biomolecular Simulation Program, Journal of computational chemistry 30, 1545-1614. Brooks, B. R., Bruccoleri, R. E., Olafson, B. D., States, D. J., Swaminathan, S., and Karplus, M. (1983) Charmm - a Program for Macromolecular Energy, Minimization, and Dynamics Calculations, Journal of computational chemistry 4, 187-217. Jo, S., Kim, T., Iyer, V. G., and Im, W. (2008) CHARMM-GUI: a web-based graphical user interface for CHARMM, Journal of computational chemistry 29, 1859-1865. Huang, J., and MacKerell, A. D., Jr. (2013) CHARMM36 all-atom additive protein force field: validation based on comparison to NMR data, Journal of computational chemistry 34, 21352145.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

31.

32.

33.

34.

35. 36. 37. 38. 39.

40. 41. 42.

43. 44. 45.

46. 47.

Vanommeslaeghe, K., and MacKerell, A. D., Jr. (2012) Automation of the CHARMM General Force Field (CGenFF) I: bond perception and atom typing, Journal of chemical information and modeling 52, 3144-3154. Jorgensen, W. L., Chandrasekhar, J., Madura, J. D., Impey, R. W., and Klein, M. L. (1983) Comparison of Simple Potential Functions for Simulating Liquid Water, J Chem Phys 79, 926935. Phillips, J. C., Braun, R., Wang, W., Gumbart, J., Tajkhorshid, E., Villa, E., Chipot, C., Skeel, R. D., Kale, L., and Schulten, K. (2005) Scalable molecular dynamics with NAMD, Journal of computational chemistry 26, 1781-1802. Darden, T., Perera, L., Li, L., and Pedersen, L. (1999) New tricks for modelers from the crystallography toolkit: the particle mesh Ewald algorithm and its use in nucleic acid simulations, Structure 7, R55-60. Nose, S. (1984) A Unified Formulation of the Constant Temperature Molecular-Dynamics Methods, J Chem Phys 81, 511-519. Hoover, W. G. (1985) Canonical dynamics: Equilibrium phase-space distributions, Physical review. A 31, 1695-1697. Humphrey, W., Dalke, A., and Schulten, K. (1996) VMD: visual molecular dynamics, Journal of molecular graphics 14, 33-38, 27-38. Delano, W. L. (2002) The PyMOL Molecular Graphics System. Frisch, M. J. T., G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ö.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. (2013) Gaussian 09, Revision D.01, Wallingford CT, . Glendening, E. D. R., J. E.; Carpenter, J. E.; Weinhold, F. Natural Bond Orbital (NBO),version 3.1. Laskowski, R. A., and Swindells, M. B. (2011) LigPlot+: multiple ligand-protein interaction diagrams for drug discovery, Journal of chemical information and modeling 51, 2778-2786. Rocchi, C., Bizzarri, A. R., and Cannistraro, S. (1998) Water dynamical anomalies evidenced by molecular-dynamics simulations at the solvent-protein interface, Phys Rev E 57, 33153325. Zichi, D. A., and Rossky, P. J. (1986) Solvent Molecular-Dynamics in Regions of Hydrophobic Hydration, J Chem Phys 84, 2814-2822. Abseher, R., Schreiber, H., and Steinhauser, O. (1996) The influence of a protein on water dynamics in its vicinity investigated by molecular dynamics simulation, Proteins 25, 366-378. Dastidar, S. G., and Mukhopadhyay, C. (2003) Structure, dynamics, and energetics of water at the surface of a small globular protein: a molecular dynamics simulation, Physical review. E, Statistical, nonlinear, and soft matter physics 68, 021921. Dastidar, S. G., and Mukhopadhyay, C. (2004) Anomalous behavior of water around sodium dodecyl sulphate micelles, Phys Rev E 70. Gilson, M. K., Given, J. A., Bush, B. L., and McCammon, J. A. (1997) The statisticalthermodynamic basis for computation of binding affinities: a critical review, Biophysical journal 72, 1047-1069.

ACS Paragon Plus Environment

Page 36 of 38

Page 37 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Biochemistry

48.

49.

50. 51.

52.

53.

54.

Swanson, J. M., Henchman, R. H., and McCammon, J. A. (2004) Revisiting free energy calculations: a theoretical connection to MM/PBSA and direct calculation of the association free energy, Biophysical journal 86, 67-74. Kollman, P. A., Massova, I., Reyes, C., Kuhn, B., Huo, S., Chong, L., Lee, M., Lee, T., Duan, Y., Wang, W., Donini, O., Cieplak, P., Srinivasan, J., Case, D. A., and Cheatham, T. E., 3rd. (2000) Calculating structures and free energies of complex molecules: combining molecular mechanics and continuum models, Accounts of chemical research d33, 889-897. Lee, M. S., Salsbury, F. R., and Brooks, C. L. (2002) Novel generalized Born methods, J Chem Phys 116, 10606-10614. Lee, M. S., Feig, M., Salsbury, F. R., Jr., and Brooks, C. L., 3rd. (2003) New analytic approximation to the standard molecular volume definition and its application to generalized Born calculations, Journal of computational chemistry 24, 1348-1356. Zhu, Y. L., Beroza, P., and Artis, D. R. (2014) Including explicit water molecules as part of the protein structure in MM/PBSA calculations, Journal of chemical information and modeling 54, 462-469. Bane, S., Puett, D., Macdonald, T. L., and Williams, R. C., Jr. (1984) Binding to tubulin of the colchicine analog 2-methoxy-5-(2', 3', 4'-trimethoxyphenyl)tropone. Thermodynamic and kinetic aspects, The Journal of biological chemistry 259, 7391-7398. Chakrabarti, G., Sengupta, S., and Bhattacharyya, B. (1996) Thermodynamics of colchicinoidtubulin interactions. Rrol of B-ring and C-7 substituent, The Journal of biological chemistry 271, 2897-2901.

ACS Paragon Plus Environment

Biochemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

For Table of Contents Use Only

ACS Paragon Plus Environment

Page 38 of 38