Article Information

• PDF (829 kb)
• Full Text with Large Figures
• Cited by in Scopus (1)

PubMed

• Articles by Tyler Luchko
• Articles by Jack Tuszynski

Related Articles

• An Analysis of Core Deformations in Protein Superfamilies

  An Analysis of Core Deformations in Protein Superfamilies Biophysical Journal, Volume 88, Issue 2, 1 February 2005, Pages 1291-1299 Alejandra Leo-Macias, Pedro Lopez-Romero, Dmitry Lupyan, Daniel Zerbino and Angel R. Ortiz Abstract An analysis is presented on how structural cores modify their shape across homologous proteins, and whether or not a relationship exists between these structural changes and the vibrational normal modes that proteins experience as a result of the topological constraints imposed by the fold. A set of 35 representative, well-populated protein families is studied. The evolutionary directions of deformation are obtained by using multiple structural alignments to superimpose the structures and extract a conserved core, together with principal components analysis to extract the main deformation modes from the three-dimensional superimposition. In parallel, a low-resolution normal mode analysis technique is employed to study the properties of the mechanical core plasticity of these same families. We show that the evolutionary deformations span a low dimensional space of 4–5 dimensions on average. A statistically significant correspondence exists between these principal deformations and the ∼20 slowest vibrational modes accessible to a particular topology. We conclude that, to a significant extent, the structural response of a protein topology to sequence changes takes place by means of collective deformations along combinations of a small number of low-frequency modes. The findings have implications in structure prediction by homology modeling. Abstract | Full Text | PDF (369 kb)

• A Molecular Dynamics Study of Slow Base Flipping in DNA usin...

  A Molecular Dynamics Study of Slow Base Flipping in DNA using Conformational Flooding Biophysical Journal, Volume 93, Issue 3, 1 August 2007, Pages 770-786 Benjamin Bouvier and Helmut Grubmüller Abstract Individual DNA bases are known to be able to flip out of the helical stack, providing enzymes with access to the genetic information otherwise hidden inside the helix. Consequently, base flipping is a necessary first step to many more complex biological processes such as DNA transcription or replication. Much remains unknown about this elementary step, despite a wealth of experimental and theoretical studies. From the theoretical point of view, the involved timescale of milliseconds or longer requires the use of enhanced sampling techniques. In contrast to previous theoretical studies employing umbrella sampling along a predefined flipping coordinate, this study attempts to induce flipping without prior knowledge of the pathway, using information from a molecular dynamics simulation of a B-DNA fragment and the conformational flooding method. The relevance to base flipping of the principal components of the simulation is assayed, and a combination of modes optimally related to the flipping of the base through either helical groove is derived for each of the two bases of the central guanine-cytosine basepair. By applying an artificial flooding potential along these collective coordinates, the flipping mechanism is accelerated to within the scope of molecular dynamics simulations. The associated free energy surface is found to feature local minima corresponding to partially flipped states, particularly relevant to flipping in isolated DNA; further transitions from these minima to the fully flipped conformation are accelerated by additional flooding potentials. The associated free energy profiles feature similar barrier heights for both bases and pathways; the flipped state beyond is a broad and rugged attraction basin, only a few kcal/mol higher in energy than the closed conformation. This result diverges from previous works but echoes some aspects of recent experimental findings, justifying the need for novel approaches to this difficult problem: this contribution represents a first step in this direction. Important structural factors involved in flipping, both local (sugar-phosphate backbone dihedral angles) and global (helical axis bend), are also identified. Abstract | Full Text | PDF (469 kb)

• The Cardiac Ca-Sensitive Regulatory Switch, a System in Dyna...

  The Cardiac Ca-Sensitive Regulatory Switch, a System in Dynamic Equilibrium Biophysical Journal, Volume 95, Issue 10, 15 November 2008, Pages 4772-4789 John M. Robinson, Herbert C. Cheung and Wenji Dong Abstract The Ca-sensitive regulatory switch of cardiac muscle is a paradigmatic example of protein assemblies that communicate ligand binding through allosteric change. The switch is a dimeric complex of troponin C (TnC), an allosteric sensor for Ca, and troponin I (TnI), an allosteric reporter. Time-resolved equilibrium Förster resonance energy transfer (FRET) measurements suggest that the switch activates in two steps: a TnI-independent Ca-priming step followed by TnI-dependent opening. To resolve the mechanistic role of TnI in activation we performed stopped-flow FRET measurements of activation after rapid addition of a lacking component (Ca or TnI) and deactivation after rapid chelation of Ca. Time-resolved measurements, stopped-flow measurements, and Ca-titration measurements were globally analyzed in terms of a new quantitative dynamic model of TnC-TnI allostery. The analysis provided a mesoscopic parameterization of distance changes, free energy changes, and transition rates among the accessible coarse-grained states of the system. The results reveal that 1), the Ca-induced priming step, which precedes opening, is the rate-limiting step in activation; 2), closing is the rate-limiting step in de-activation; 3), TnI induces opening; 4), there is an incompletely deactivated population when regulatory Ca is not bound, which generates an accessory pathway of activation; and 5), there is incomplete activation by Ca—when regulatory Ca is bound, a 3:2 mixture of dynamically interconverting open (active) and primed-closed (partially active) conformers is observed (15°C). Temperature-dependent stopped-flow FRET experiments provide a near complete thermokinetic parameterization of opening: the enthalpy change (Δ=−33.4kJ/mol), entropy change (Δ=−0.110kJ/mol/K), heat capacity change (Δ=−7.6kJ/mol/K), the enthalpy of activation (=10.6kJ/mol) and the effective barrier crossing attempt frequency (=1.8×10 s). Abstract | Full Text | PDF (616 kb)

• …more

Copyright © 2008 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 94, Issue 6, 1971-1982, 15 March 2008

doi:10.1529/biophysj.107.115113

Biophysical Theory and Modeling

Conformational Analysis of the Carboxy-Terminal Tails of Human β-Tubulin Isotypes

Tyler Luchko^*, †, J. Torin Huzil^‡, Maria Stepanova^† and Jack Tuszynski^*, ‡, ,  

^* Department of Physics, University of Alberta, Edmonton, Alberta T6G 2G7, Canada
^† National Institute for Nanotechnology NRC, Edmonton, Alberta T6G 2M9, Canada
^‡ Department of Oncology, University of Alberta, Edmonton, Alberta T6G 1Z2, Canada

Address reprint requests to Jack Tuszynski, Room 3336, Cross Cancer Institute, 11560 University Avenue, Edmonton AB T6G 1Z2, Canada. Tel.: 780-432-8906; Fax: 780-432-8892.

Microtubules (MTs) are hollow cylinders constructed from linear chains of the protein tubulin. Tubulin is highly conserved throughout the entire eukaryotic kingdom, and the two main classes, α/β, are expressed from multiple genes, producing several isotypes with seemingly identical functions 1,2,3. At the cellular level, it has been hypothesized that MT stability is regulated by subtle variations observed between α and β isotypes 4,5,6. There is ∼80%–95% sequence identity between isotypes; however, the extreme carboxy-terminal tails (CTTs) exhibit considerable differences, having only 50%–60% identity in this region 7,8. Early phylogenetic comparisons of the vertebrate β-tubulin families identified the CTTs, along with an internal variable domain, as the primary isotype defining features of the β-tubulin protein 9.

The importance of the CTTs has been demonstrated through their removal using limited proteolytic cleavage by substilisin. After cleavage, the critical tubulin concentration required for polymerization was found to be ∼50 times lower than that for intact tubulin 10. This rate was shown to decrease further through the addition of MT-stabilizing compounds such as paclitaxel 11. Proteolyzed tubulin also exhibits altered protofilament-bending, resulting in the formation of sheets, bundles of twisted filaments, rings, unstructured aggregates, or MTs with reversed polarity 10,12. Finally, the presence of the highly charged CTTs is thought to obstruct tubulin/tubulin interactions, regulating the rate and conformation of MT assembly through unfavorable electrostatic interactions 13. This hypothesis is supported by the observation that divalent cations, or the substitution of acidic residues, results in an increased rate of MT assembly and a decreased rate of MT disassembly 14,15.

In addition to tubulin interactions within the MT itself, the functional stability of MTs in vivo has been shown to depend on interactions with MT-associated proteins (MAPs) through interactions with the CTTs. For example, the correct assembly of the kinetochore and the integrity of the mitotic spindle is dependent on the ability of the Dam1 complex to bind the CTTs of β-tubulin 16. Interactions between the CTTs and MAP2 or MAP τ have also been shown to result in the stabilization of MT structures in vivo 17,18. Interactions between kinesin and MTs may also be significantly modulated through direct contacts with the CTT 19 or through the indirect modulation of kinesin's ability to bind ATP 20. Recent evidence also suggests that the CTTs may play a role in apoptosis, as interactions between MTs and the proapoptotic and antiapoptotic members of the Bcl-2 family are also affected by the presence of the CTT regions of α- and β-tubulin 21, an observation that may offer a plausible explanation for the previously observed interactions between Bcl-2 and paclitaxel 22.

Although the presence of intact CTTs is seemingly essential for proper MT assembly and function, little is known about their structure and how this may affect interactions with other proteins. Several early studies have examined possible conformations of the CTTs using NMR and circular dichroism (CD) spectroscopy and demonstrated that there was inherent disorder within the CTTs 23,24. Although a region of increased helicity toward the amino terminus was identified, a finding that was subsequently confirmed by electron crystallographic analysis, this has provided little additional information about CTT structure, as the last 10 residues of α-tubulin and the last 18 residues of β-tubulin were not visualized due to lack of density 25. A possible explanation for this lack of density is the nonhomogenous presence of tubulin isotypes in the MT preparations used in the crystallographic analysis. However, homogenous samples of αβII and αβIII tubulin failed to improve the quality of density maps, and Nogales et al. 25 concluded that their lack of resolution was indeed due to disorder in this region.

Fortunately, molecular modeling provides us with the unique ability to examine conformations of the CTTs at a level of detail experimental analysis is unable to yet provide. However, even modeling has its limitations in this regard, where a persistent problem with low temperature protein-folding simulations is that of obtaining adequate sampling. This problem exists because simulations generally become trapped in one of a large number of local energy minima. Several generalized ensemble algorithms, based on non-Boltzmann probability weight factors, are capable of overcoming this problem by introducing a random walk in energy space 26. However, it is often not a trivial matter to determine the non-Boltzmann weight factors, and for this reason we chose replica exchange molecular dynamics (REMD) as a method to examine possible conformations of the tubulin CTTs 27,28. REMD utilizes a large number of parallel simulations at different temperatures, with exchanges between trajectories attempted periodically using Metropolis criteria. All replicas of the system then perform a random walk through temperature space and, as a consequence, also through energy space. A replica may therefore overcome an energy barrier through exchange with replicas at higher temperatures. This allows configurations to be sampled at a given temperature on timescales not otherwise possible while still maintaining a thermodynamically consistent ensemble of configurations. Here, we discuss the creation of nine models of human β-tubulin CTT peptides using REMD and their analysis for relative conformational flexibility within the ensemble.

Based on our previous examination of β-tubulin isotypes, we chose a set of nine peptides corresponding to the consensus CTT sequences of βI (GI:34222261), βII (GI:68299771 and GI:42476191), βIII (GI:50592995), βIVa (GI:21361321), βIVb (GI:68051719), βV (GI:14210535), βVI (GI:41152077), βVII (TUBB4Q) (GI:55770867), and βVIII (TUBB8) (GI:42558278) (Table 1) 8. Accession numbers correspond to Entrez Nucleotide expressed mRNA sequence IDs contained only within the human genome. Using the crystallographic structure of tubulin (Protein Data Bank (PDB) ID: 1JFF) as a guide, we selected the conserved DA[TK]A motif that is positioned at the extreme end of helix H12 as the initiation point for the construction of all the peptides 29. This position marks the beginning of the structurally undefined region of the CTT in the structure and the domain examined using NMR and CD spectroscopy 23,24. All nine CTTs were prepared identically in an extended conformation, having the N-terminal charge neutralized by capping with an acetyl group using the PyMol v0.99 residue and fragment builder facility 30. The conventional protonation state at pH 7.0 was used for all residues, with a focus on His, which was protonated on the ɛ nitrogen throughout.

All calculations utilized the AMBER99 force field 31,32, which was applied using PDB2GMX in GROMACS 3.3 33. Each peptide was placed in a rhombic-dodecahedron unit cell consisting of ∼2800 TIP3P waters 34. Sodium and chlorine counterions were added to the most energetically favorable locations (as determined using GENION) such that the net charge of the system was neutralized and a final ion concentration of 100mM was established. For all molecular dynamics (MD) and REMD calculations, rigid bonds were maintained using the LINCS algorithm for the peptide and SETTLE for waters. Particle Mesh Ewald (PME) was used for electrostatic interactions with a cutoff of 0.8nm and a Lennard-Jones cutoff of 1.0nm. Constant pressure was maintained by allowing the box size to fluctuate isotropically. Each system was minimized using a low-memory Broyden-Fletcher-Goldfarb-Shanno minimizer in GMXRUN until machine precision was achieved. This was followed by 20ps of heating and 1.1ns of equilibration.

To achieve a transition probability of 0.3, 43 target temperatures between 273 and 382K were selected. Each replica was then heated/cooled to its target temperature over 50ps and simulated without exchange for 500ps, followed by 500ps of REMD. Production REMD runs consisted of 10ns of dynamics for each replica with exchanges attempted every 5ps. Energy and conformational snapshots were saved every 1 and 10ps, respectively. For each isotype, this produced an ensemble of 43,000 conformations and an aggregate simulation time of 430ns over all temperatures.

A cluster analysis of the resulting REMD conformations was used to determine preferred conformations and relative populations of each peptide. The GROMOS clustering algorithm 39, implemented in G_CLUSTER, was used for this purpose. All Cα atoms were root mean-square deviation (RMSD) fit to the extended starting structure, and an RMSD cutoff was set for the Cα atoms; for each structure from the ensemble, all structures within this cutoff were assigned as neighbors. The structure with the most neighbors was then identified as the center of a cluster and was removed from the pool with all its neighbors. This process was repeated until all structures had been assigned to a cluster. The RMSD cutoff was chosen such that the largest cluster contained 50% plus 1% of the structures.

As MD simulations tend to produce immense quantities of data, principal component analysis (PCA) is a powerful mathematical tool used to detect correlations in MD trajectories 40,41. To perform PCA, all of the conformations from the REMD ensemble were RMSD least squares fit to the reference structures, effectively removing all rotations and translations. Then, for non-mass-weighted PCA, the covariance matrix can be calculated as

where r₁…r_3N are the Cartesian coordinates of the N atoms used in the analysis and 〈…〉 denotes the ensemble average. The resulting matrix can then be diagonalized, and the resulting 3N-dimensional eigenvectors, are organized in descending order of eigenvalues, . Eigenvalues represent the variance along their associated eigenvector, and the larger the eigenvalue the more significant the correlated motion.

A PCA of all nine isotypes of β-tubulin was performed at 311K using the positions of all Cα atoms. To mimic the CTTs bound at their N-terminus, the backbone atoms of the first three amino acids of the ensemble were RMSD fit to the reference structure for the isotype at 311K, producing an average structure (Fig. 1). This was essential for the physiological relevance of the PCA calculation and the considerable motion of the CTTs. Eigenvalues of different isotypes cannot be directly compared, as different numbers of atoms were used in the covariant analysis. This can be overcome by normalizing the eigenvalues by the number of atoms used in the analysis:

Display large version of this figure

Figure 1

Human β-tubulin CTT sequence alignment and consensus. Illustration of reference and average structures used for PCA. The N, C, and Cα atoms from the three N-terminal residues of each conformation from the βI ensemble were RMSD fit to a representative structure backbone. The backbone atoms of 50 randomly selected conformers of the ensemble are shown in black, gray, and white to delineate conformers in no specific order. The first four residues of each peptide are shown in yellow. The average structure, used for subsequent PCA calculations, is shown in light blue. Image produced using VMD 63.

where σ_i is the root normalized eigenvalue representing the standard deviation along the ith eigenvector.

Sequence alignments of CTTs were performed with the default values using the European Bioinformatics Institute Clustal server (http://www.ebi.ac.uk/clustalw/), with the exception that the extended gap penalty was increased to a value of 0.5 (Table 1).

Structural motifs within the CTT ensemble were identified by performing pairwise alignment of each CTT sequence using the ClustalW alignment facility in MacVector (MacVector, Cary, NC). Motifs were identified as those sequences containing at least four consecutive identical residues. A structural similarity search was then performed using only sequence runs of four or more identical or similar residues. Structural similarity was determined by calculating the ρ_sc between all pairs of conformations from the ensembles of each isotype to create a 1001×1001 matrix. Here, ρ is defined as

The absolute length of each human β-tubulin CTT ranges from 9 to 25 residues; however, their typical length falls between 18 and 20 residues (Table 1). The exceptions are βVII, which contains only nine residues, and βIII/βVI, which have 24 and 25 residues, respectively. Although the overall amino acid composition of each CTT is quite similar, there are some notable factors to consider. When we examine all nine human CTT sequences in aggregate, the residue that occurs most frequently is Glu (42%). The occurrence of all other residues drops significantly from this value: Ala (18%), Asp (12%), Gly (8%), Thr (5%), Val (4%), and Phe (3%). The βVIII CTT contains no Gly, the βIVa/b CTTs have a proportionally low Asp content, and, finally, the βVII CTT has a reduced Glu content, which may simply be a result of its length. Interestingly, the remainder of the amino acids occurs exclusively in the βIII (Lys, Met, Pro, Gln, Ser, and Tyr), βV (Ile, Asn), and βVI (His, Lys, Leu, and Pro) CTTs at frequencies of 1% or less. Finally, Cys, Trp, and Arg are absent from all the CTT sequences.

A notable exclusion from these simulations is the tubulin protein itself, a factor that will undoubtedly influence the results as they have been explained here. The decision to exclude the bulk tubulin was made due to its large system size and the computational resources that would be required to explore the conformational space of an entire tubulin dimer. Therefore, simulations of only the CTTs become a problem of protein folding, and adequate sampling is critical to the proper interpretation of the results. To increase sampling efficiency at different temperatures, we used replica exchange MD (REMD), which does not provide information on the dynamics of the system. However, as our ultimate interest is in the overall conformations of the CTTs, the loss of dynamics data was an acceptable compromise to gain increased sampling.

The completeness of sampling can be determined by calculating the normalized overlaps between two different parts of an MD trajectory. Such overlaps can indicate whether or not both parts are spanning over the conformational space equally and not diffusing to new parts 43,44,45. Subspace overlap between two sets of n orthonormal vectors and is defined as

When the overlap has a value of 1, the sets v and w can be considered to span the same subspace. This measure, however, does not account for the magnitudes of the eigenvalues, meaning that differences between all eigenvectors contribute equally. Furthermore, when two or more eigenvalues are equal, the corresponding eigenvectors are random, causing a random variation in the subspace overlap number. An alternative method, suggested by Hess 45, is the normalized overlap, defined as

where d is the difference between the covariant matrices, C₁ and C₂, and tr is the trace. Here, if the overlap is 0, then the two sets are considered to be orthogonal, whereas an overlap of 1 indicates that the matrices are identical. Covariant analysis of the trajectories from each isotype at 311K, chronologically divided into thirds, was performed using the same procedure used for the PCA. The subspace and normalized overlaps calculated between each of these thirds are reported in Table 2. The high overlap between the thirds indicates that each part of the simulation is sampling approximately the same conformational space, and it is unlikely that there are unexplored regions missed earlier in the run. Although not a guarantee of complete equilibrium sampling, we concluded that the overlap using both of the above mentioned methods is acceptable and that adequate sampling within all of these systems has been obtained.

Clustering each CTT peptide can provide a representation of probable folded conformations; however, a bell-shaped ρ_sc distribution about a representative structure demonstrates that there is actually no native folded conformation for most of the isotypes (Fig. 2). The only exception was the βVII CTT, which exhibited a significant population of folded structures with a ρ_sc<0.5 and a second population of unfolded structures with a ρ_sc>0.5. These results recapitulate previous observations that the CTTs are extremely flexible, and any structures within them are transient in nature and cannot be captured by the clustering methodology used in this study. Therefore, we felt that a more appropriate analysis was to consider the ensemble average of secondary structures as calculated from the entire 10ns of ensemble conformations at 311K by STRIDE 46 (Fig. 3). These results demonstrated that, whereas the ρ_sc distribution showed no native folded conformations, many of the CTTs contain regions that are either α- or 3–10-helical at least 40% of the time.

Display large version of this figure

Figure 2

ρ_sc distributions of CTTs about representative structures. The large spread in the histograms of the ρ_sc distributions for each CTT indicates the lack of a single folded conformation. A distinct population of 67% conformers have a ρ_sc<0.5 for the βVII isotype, suggesting a distinct folded conformation for this CTT.

Display large version of this figure

Figure 3

Time average of each type of secondary structure for each CTT. The total fractional secondary structure content was determined as an average over the entire 10ns MD simulation at 311K. The fractional time averages of each type of secondary structure are stacked to sum up to one. Although experimental studies of tubulin suggest that the CTTs are unstructured, these results suggest that many of the CTTs contain a significant amount of transient helix.

Having established the presence of transient secondary structures within the ensembles, we performed pairwise alignments of all the CTT sequences to identify potential sequence motifs with which to correlate our modeling results. A ρ_sc matrix was then calculated for all the conformations across each isotype at 311K. Only those motifs having a high fraction of low ρ_sc were determined to be significantly similar (ρ_sc of 0.3 or less) (Table 3). Through this analysis, we identified two motifs that showed both sequence and conformational similarities. To visually compare the motifs, the conformations of the aligned subsequences were clustered as described in Materials and Methods. Illustrated in Fig. 4 are RMSD alignments of representative structures (those with the most neighbors) from each isotype containing the respective motif. The first of these motifs was identified as a probable casein kinase-2 (CK-2) binding motif at the N-terminal end of the peptide (Figure 4A). This motif was also independently identified using a Prosite search for motifs within each of the CTT sequences 47 (not shown). The second motif was determined to correspond to a previously identified MAP2 binding motif found within α-tubulin (Figure 4B) 48. We should note that although a common conformation for each of the motifs across isotypes has been identified, these are not stable folds; and depending on the isotype and residues included in the search, anywhere between 1% and 69% of the ensemble structures have a ρ_sc<0.3 when compared to the motif conformation.

Display large version of this figure

Figure 4

Sequence motif identification. High scoring sequence pairs from Table 3 were RMSD fit to illustrate the overall fold and position of side chains in the structural motifs. In both images, red residues are acidic, white residues are nonpolar, green residues are polar, and all side chains are blue. Image created using VMD 63.

Not only is the secondary structure of each CTT significant, but their flexibility is also of critical importance when attempting to understand how they are able to conform when interacting with MAPs. Because of the complex interplay of many factors, sophisticated statistical analysis of conformational ensembles, such as PCA, appears to be the most appropriate computational tool to characterize the flexibility. PCA was performed on the ensemble at 311K, with the backbone atoms of the first three residues RMSD fit to a reference structure (Fig. 1). Resulting eigenvectors were ordered by descending eigenvalues, which represent the variance of the motion along the principal components. In Fig. 5, the six largest root normalized eigenvalues are shown for each isotype. Except for βVII and βIVa, to a lesser extent, three components can be identified with prominent eigenvalues of comparable magnitude, which is significantly higher than the magnitude of other eigenvalues. The three components with the largest eigenvalues represent correlated motions of the peptide fragments with the most significant standard deviations of the motion along the corresponding orthogonal directions. Since the standard deviations shown in Fig. 5 were normalized for the number of residues, they can be employed as a universal measure for comparison of conformational motions in different CTT peptides.

Display large version of this figure

Figure 5

Root mean-square normalized eigenvalues for each isotype. The magnitudes of the first six root normalized eigenvalues (white to black) are shown for each isotype and represent the standard deviation of the motion along the eigenvector normalized for the overall length of the peptide. The relative size of the root-normalized eigenvalues indicates the relative flexibilities of the individual peptides.

A more detailed representation of the conformational motions in peptides is provided by projections of the ensemble of conformations onto the planes spanned by the most important principal components (Fig. 6). Here, the spatial distributions of occupancies of the various conformational states are shown over the planes spanned by the first and second and by the first and third principal components. A representative isotype, βV, illustrates that both distributions are smooth and without evidence of considerable clustering. Similar behavior exists in all isoforms considered except for βIVa and βVII, which show significant clustering of conformational occupancy, with maxima at (0.1,0) and (−0.2,0), respectively. The widths of the distributions shown in Fig. 6 are not directly dependent on the length of the peptide, as the eigenvalues have been normalized by the number of residues. The reduced width of βVII, therefore, indicates its reduced mobility. The bin size of the histogram has been chosen in proportion to the width of the distribution; so the maxima seen for βIVa and βVII are not artifacts of the sampling but represent regions of increased occupancy. Since the occupancy of a conformational state is inversely proportional to the free energy of that state 49,50, the lack of clustering in the other isotypes suggests that the global minima are broad or that there are generally no significant local minima. In contrast, βVII has a significant basin of attraction, indicative of a folded conformation. This observation is consistent with the preceding cluster analysis, which identified a well-defined folded state in the βVII isoform only. Although less defined, the occupancy distribution for βIVa also exhibits a preferred folded conformation, demonstrating that the PCA methodology is more sensitive to transient structural motifs than the clustering analysis. One more noteworthy feature is that the two-dimensional distributions over the first and second components are largely similar to those over the first and third components, suggesting a significant level of symmetry among principal components.

Display large version of this figure

Figure 6

Projections of the ensemble of conformations onto the planes of the three most important principal components. The first and second principal components (upper row) and the first and third principal components (lower row) are plotted on the x and y axes, respectively, for βV, βIVa, and βVII. The histograms represent the occupancies of the corresponding conformation states, with lighter colors indicating more frequently visited areas. The same coloring scheme is used for all isotypes and is capped at 20 counts. βV is typical of the remaining isotypes in the width of the distribution and lack of a significant maxima.

In addition to information regarding the motion of each CTT peptide, PCA provides the ability to compare relative flexibilities. As we studied only the nine human β-tubulin isotypes, it is not possible to comment on the flexibility of the CTT peptides in an absolute sense. However, there have been a few studies discussing the flexibility of small peptides. Although they use a different approach to calculating the flexibility, Ma et al. made a survey of 28 short peptides and determined that the native helical structures of the peptides were more flexible than random or disordered conformations, in agreement with our observations (Figure 7B) 51. Unfortunately, their methodology, calculating the vibrational free energy, does not allow comparison of flexibility between peptides. Here, as a measure of the relative flexibility of each peptide, the eigenvalues of principal components were compared. To characterize the flexibility, we employed the value where σ₁ and σ₂ are the normalized eigenvalues for the first and the second principal components, respectively. The magnitude of the PCA eigenvalue parameter σ₁₂ can be viewed as the mean square width of the two-dimensional REMD configuration ensembles shown in Fig. 6, and thus σ₁₂ can be directly interpreted as the flexibility of the peptide associated with its major correlated motions. Although the CTTs are characterized by three principal components (Fig. 5), the considerable symmetry between the second and third component demonstrated in Fig. 6 allows the use of only two components out of three to quantify the flexibility of CTTs. The correlations between the normalized distance (the distance between the N- and C-termini Cαs, divided by the number of residues), the time-averaged helical content, and the average clustering parameter 〈ρ_sc〉 with the value σ₁₂, which we use as the major measure of flexibility, can be seen in Fig. 7.

Display large version of this figure

Figure 7

Peptide's flexibility indicators. The end-to-end distance normalized by the number of residues (A), time average helicity (B), and (C) of each CTT peptide at 311K plotted as functions of The plots demonstrate the correlation between each of the above values and the major flexibility indicator, In b, βVI can be excluded from the trend due to the presence of a Pro, which disrupts the secondary structure but contributes to the flexibility of the peptide. The error in all three plots is the standard error calculated using box averaging 64.

The CTT's involvement in MT regulation is critical; it is essential to understand the conformational differences adopted by different tubulin isotypes. We suggest that the sequence variability within the CTTs may have arisen as a mechanism to conserve the overall tubulin structure to maintain proper MT assembly, while still providing a flexible, solvent-exposed region of increased variability for interactions with proteins that affect MT function (see alignment in Table 1). This is both an attractive and reasonable hypothesis, as the role of intrinsically disordered proteins in protein interactions has been implicated as a mechanism to enhance specificity 52. The following discussion addresses three main points about the CTTs and their interactions with MAPs: first, the relative flexibilities of each of the CTT peptides; second, the presence of transient structure within these peptides; and finally, how this may influence protein interactions. We will discuss this possibility further when examining CTT flexibility in the context of MAP interactions, in particular interactions with the CK-2 and MAP2 binding motifs that we identified in this work.

Accurate flexibility and secondary structure measurements are critical when the goal is to propose mechanisms for MAP interactions with the large number of possible CTT conformations on the MT surface. The overall distribution of amino acids within the CTT will obviously have an effect on their flexibility and result in secondary structure. For example, the overall helical propensity of βII, βIII, βIVa, βIVb, and βV is interrupted by the presence of Gly (Fig. 3). Additionally, regions with uninterrupted stretches of Glu also tend to produce regions with greater helical propensity, providing an argument for the prevalence of this residue within the CTT sequences. This observation can be compared with the results of Roe et al. for simulations performed on decaalanine 53. Although shorter than all but the βVII peptide, the overall α- and 3–10-helical content within decaalanine was similar to the helical content observed here. However, through the analysis of the distribution of secondary structure along the CTT sequences, we note that our results differ from that of Roe et al. in that the helical content of the CTTs drops near the middle of the sequence for most isotypes whereas for decaalanine the central residues have a maximal amount of helix.

Three characteristics emerge which are representative of the behaviors of CTTs: the end-to-end distance, normalized by the number of residues, the PCA-based eigenvalues, and the averaged helical content (Figure 7AB). To test the applicability of PCA as a measure of relative flexibility, we compared PCA-based flexibilities with the average parameter of structural similarity 〈ρ_sc〉, which we computed from the distributions in Fig. 2 (Figure 7C). Although less detailed than PCA, the distributions of conformations over the parameter ρ_sc provide an alternative statistical characterization of the CTTs configurations, and thus they can be expected to correlate with the PCA-based results. The correlation between the PCA-based flexibilities and 〈ρ_sc〉 is more pronounced than any other correlation that we investigated. Most of the CTTs show a clear proportionality between σ₁₂ and 〈ρ_sc〉, and the variability of 〈ρ_sc〉 that corresponds to similar σ₁₂ is <15%. This similarity of results from two fundamentally different statistical tools demonstrates that these statistical methodologies are the best suited to characterize CTTs’ flexibility. Interestingly, rather than the expected exponential decay, the magnitudes of the first three eigenvalues are of similar magnitude and significantly larger than the remaining eigenvalues, indicating an isotropy of occupancies in the three-dimensional space (Fig. 5).

More simply, the CTTs are flexible enough for the ensemble of their configurations to be spherically symmetric with respect to an immobilized base of amino acids (Fig. 1). Since these three eigenvalues are significantly larger than all the others, three-dimensional symmetry is likely reached through a few highly flexible bonds within the peptides that allow rotations of large angles, thereby generating nearly spherically symmetric distributions of occupancies. In Fig. 6 this is illustrated by a comparison of the histogram for βV, which is representative of a significant level of symmetry typical for most isotypes, with the less symmetric ones for βIVa and βVII. As the magnitudes of the eigenvalues shift to an exponential decay, spherical symmetry is lost and correlated motions begin to appear. At the same time, an individual CTT may preserve relatively stable motifs contained within some secondary structure. In terms of the PCA results presented in Fig. 5, these motifs are represented by the fourth, fifth, and higher eigenvalues. The fact that these components, although significantly less pronounced than the first three, still have considerable nonzero eigenvalues indicates that the conformational motifs within the CTTs are of transient nature and subject to variability over the trajectory. This is consistent with our observation of a significant amount of transient helicity (Fig. 3), indicated by a substantial amount of 3–10-helix.

Considering the PCA eigenvalue parameter σ₁₂ as the measure of flexibility, it is evident that the βIII and βVI CTTs are the most flexible of all CTTs. These are the longest two sequences, which may account for some of their flexibility, but interestingly they also have the shortest end-to-end distance per residue. The cause of these shared structural properties differs in the two cases: whereas the βIII CTT has among the highest degree of secondary structure, the βVI CTT has among the lowest. The βVI CTT differs from all the other fragments, as it contains a Pro near the C-terminus of the peptide. This Pro breaks the helical structure of the fragment and, at the same time, does not facilitate an extended conformation. This accounts for the lack of secondary structure but is in contrast to βIII, which lies predominantly in the α-helical region of the Φ/Ψ distribution. When accounting for the lone Pro and omitting the data for βVI, the correlation between helicity and flexibility becomes greater. The βII, βIVb, and βV CTTs show the same spring-like flexibility of βIII, exhibiting a similar helical content and moderate end-to-end distance per residue. These CTTs also contain Gly at the same position that breaks the helical content of the fragment. In contrast, the βIVa and βVII CTTs displayed the least amount of flexibility. In the case of βIVa, one reason for this reduction in flexibility may be the early occurrence of Gly in its sequence, which stiffens this region of the peptide by decreasing the helicity. This decrease in helicity in the first seven residues, compared to βIVb, can be seen in Fig. 3. βIVa also contained the most stable DA[TK]AEE motif (data not shown), which resulted in the largest normalized distance of all the common isotypes.

The observation that peptides containing a significant amount of transient helical conformations are more flexible than an extended one may seem counterintuitive. The correlation between the normalized distance and flexibility in Figure 7A is most easily understood in terms of entropy. Take the following example: the freely jointed, or ideal chain, polymer model states that compact forms of the polymer contain more states that are accessible (i.e., have greater entropy) than do extended states 54. Thus, an entropic force drives the polymer or polypeptide into compact conformations. A greater number of accessible states suggests that the chain is more flexible and that extended states are more rigid. Although the ideal chain has no internal energy and is free to collapse, a peptide fragment does contain internal energy that can stiffen the peptide. The source of the force that causes the peptides to be in a rigid, extended state is the peptide backbone. Individual residues that induce a compact conformation through helical propensity or residues that induce a turn (such as proline) provide flexibility to the structure of the peptide while reducing the normalized distance. However, there is not a complete one-to-one correspondence between the flexibility and the normalized distance. For example, βI and βVIII have a normalized distance comparable to βIVa but a σ₁₂ closer to the more compact βV. Although the normalized distance of the peptide is a strong indicator of flexibility, the details of the compact structure, such as the degree of helicity, also contribute. We are not the first to observe this phenomenon in peptides; Ma et al. calculated the vibrational free energy of 28 short peptides in native helical and random extend states and determined that the native helical structures of the peptides were more flexible in all cases 51.

It is becoming clear that each β-tubulin isotype has a unique pattern of expression ranging from specific for βII, βIII, βIVa, and βVI to constitutive for βI and βIV 56,57. The βI isotype is the most commonly expressed in humans and, as such, is also the most common isotype found in cancer cells 58. Our observations also suggest that the CTT of this isotype shows “typical” behavior in our analysis. It is difficult to correlate our results to biological function; however, a clear pattern illustrating βIII as a statistical outlier has emerged. The βIII CTT was one of the most flexible peptides and contained a large amount of transient secondary structure. This is interesting, as βIII-tubulin has been observed at increased levels in human tumors and implicated in the development of drug resistance to standard chemotherapy treatments 59,60,61. In addition to the altered expression of βIII, the alteration of tubulin regulatory proteins such as MAP4 has also been implicated in changes to MT dynamics and the development of drug resistance 62. Because several MAPs have now been shown to interact directly with the CTTs 16,17,18,19,21, the differences that we observed here may have a significant impact on their affinity to the surface of an MT.

Using the conformational similarities for each CTT peptide, in combination with sequence comparisons, we identified two distinct transient structural motifs (Table 3). Several similar regions were observed as significant conformations throughout the trajectory files, indicating that they may be common motifs that could potentially be recognized by MAPs. The first motif was identified as a putative CK-2 phosphorylation site for which the consensus motif has been identified as [ST]XX[DE] 47. All of the CTTs, with the exception of βVI and βVII, were shown to contain this motif within the conserved DA[TK]A sequence. The second motif identified was the MAP2 binding site, for which the consensus sequence has previously been characterized as EEAEEE 48. Only the CTT of βI, βIVa, and βIVb strictly contain this motif; however, most of the other CTTs also contain a similar motif (Table 3). The conformations identified for βI, βIVa, βIVb, and βVI most commonly form a transient helix that is between four and five residues in length. In βVIII, there seems to be an increasing amount of unstructured loop associated with this region, which is most likely a result of the displacement of the amino end due to the presence of a bulky Tyr residue. Both the CK-2 and MAP2 motifs display a significantly higher degree of helical propensity when compared to the entire ensemble of CTT conformations (Fig. 3). Interestingly, with the exception of βVI and βVII, the increased helical propensity of the proposed CK-2 domain suggests that the H12 helix may extend several additional residues farther than those that have already been observed in the crystal structures of tubulin.

Finally, the motifs identified here are not stable structures but rather commonly occurring, transient conformations. For example, in the “stiffest” isotype, βIVa, the CK-2 DATA motif is occupied 69% of the time whereas the MAP2 motif EEEV is occupied 17% of the time (Table 3). The two motifs are simultaneously folded in the ensemble of CTT trajectories only 12% of the time. As the sequences contained within the two motifs are extended to include DATAEE and AEEEVA, respectively, this drops to 2% (Table 3). Rather, the shared similarity of the folded states indicates that these conformations are individually inducible, particularly when binding to their respective substrates. As we are dealing with a subset of all CTTs, we were unable to determine accurately if any additional structural or sequence motifs exist in the tails. A more accurate determination of CTT motifs would require a larger survey of CTTs across all identified α- and β-tubulin isotypes (see Table 1).

The role of the unstructured tubulin CTTs, particularly concerning MT stability and interactions with MAPs, is a question that has yet to be conclusively answered. Here we have shown, using a combination of REMD and PCA, that although the experimental structures of tubulin suggest the CTTs are unstructured, many contain a significant amount of transient conformations that are linked by a few highly rotatable bonds. We believe that the role of the transient secondary structure, coupled with increased CTT flexibility, may be twofold. First, global flexibility coupled with weak helical tendencies could provide a mechanism that allows the CTTs to search conformational space and hence easily bind to other proteins. This could provide a convenient mechanism for MAPs to dock to the MT surface. Second, increased flexibility could also enhance binding affinities to MAPs as a result of their ability to conform to a specific binding site. The presence of well-defined structural motifs within the CTTs may, consequently, play an important role in binding specificity, where the presence of smaller domains provides a common structural signal for binding, while the global flexibility of the CTT is maintained.

The role of β-tubulin CTTs in MT function is appealing, since subtle differences in a small stretch of amino acids could have profound consequences on MAP interactions. Unfortunately, it is extremely challenging at this time to elicit tight correlations between the results presented here and their biological consequences, as there is very little experimental data regarding CTTs explicitly in the context of tubulin isotypes. As interest in the properties of tubulin isotypes and their potential role in targeted chemotherapy treatments increases, we anticipate a significant body of data to become available soon. As a result, a clear understanding of CTT structure and function may provide the means by which we can begin to rationally design and develop novel drugs or peptide mimetics that target CTT binding sites, providing better selectivity for chemotherapies.