Summary: Wordom is a versatile program for manipulation of molecular dynamics trajectories and efficient analysis of simulations. Original tools in Wordom include a procedure to evaluate significance of sampling for principal component analysis as well as modules for clustering multiple conformations and evaluation of order parameters for folding and aggregation. The program was developed with special emphasis on user-friendliness, effortless addition of new modules and efficient handling of large sets of trajectories. Availability: The Wordom program is distributed with full source code (in the C language) and documentation for usage and further development as a platform-independent package under a GPL license from http://www.biochem-caflisch.unizh.ch/wordom/ Contact: caflisch@bioc.unizh.c

Caflisch, Amedeo

Cecchini, Marco

Rao, Francesco

Seeber, Michele

Settanni, Giovanni

RERO DOC Digital Library

Vol. 23 no. 19 2007, pages 2625–2627BIOINFORMATICS APPLICATIONS NOTE doi:10.1093/bioinformatics/btm378Structural bioinformaticsWordom: a program for efficient analysis of moleculardynamics simulationsMichele Seeber, Marco Cecchini, Francesco Rao, Giovanni Settanni and Amedeo Caflisch*Department of Biochemistry, University of Zurich, CH-8057 Zurich, SwitzerlandReceived on May 9, 2007; revised on July 11, 2007; accepted on July 15, 2007Advance Access publication August 23, 2007Associate Editor: Alfonso ValenciaABSTRACTSummary: Wordom is a versatile program for manipulation ofmolecular dynamics trajectories and efficient analysis of simulations.Original tools in Wordom include a procedure to evaluate signifi-cance of sampling for principal component analysis as well asmodules for clustering multiple conformations and evaluation oforder parameters for folding and aggregation. The program wasdeveloped with special emphasis on user-friendliness, effortlessaddition of new modules and efficient handling of large sets oftrajectories.Availability: The Wordom program is distributed with full sourcecode (in the C language) and documentation for usage and furtherdevelopment as a platform-independent package under a GPLlicense from http://www.biochem-caflisch.unizh.ch/wordom/Contact: caflisch@bioc.unizh.ch1 INTRODUCTIONMolecular dynamics (MD) simulations can produce severalterabytes of data (Rueda et al., 2007) whose manipulation andanalysis is problematic in the absence of proper tools. As anexample, sampling in the millisecond time-scale generates onthe order of 107 to 108 coordinate frames (Settanni et al., 2005).In this communication, a program is presented that caneffectively access MD trajectory data and make them availablefor conversion, manipulation and analysis by different mod-ules. In this way, it is possible to implement new procedures byfocusing on the relevant equations and algorithms withouthaving to deal with the data access module. The programprovides general tools for analysis as well as special purposemodules which are not available in other programs (e.g. theefficient clustering of large sets of conformations and theevaluation of order parameters for monitoring polypeptideaggregation).2 DATA MANIPULATION AND OUTPUTWordom can perform simple conversion and manipulationtasks, such as the extraction of one or more snapshots from aCHARMM (Brooks et al., 1983), GROMACS (van der Spoelet al., 2005), or NAMD (Kale´ et al., 1999) trajectory into a pdb(Berman et al., 2000) or crd file (Brooks et al., 1983), adding asingle coordinate set from a pdb or crd file to a trajectory,merging segments of trajectories, calculating the averagestructure along a trajectory, extracting xyz coordinates ofselected atoms, showing/modifying the content of a trajectoryfile header, etc. For most of these tools an atom selectionmechanism is implemented. Notably, all of these functionsare accessible with a single command line. Moreover, thecoordinate and trajectory files generated by Wordom can bevisualized by widely distributed graphics packages, e.g. VMD(Humphrey et al., 1996) and Pymol (DeLano, 2002).3 ANALYSIS OF MD SIMULATIONSWordom provides modules for the analysis of geometricalfeatures (distances, contacts, dihedrals, hydrogen bonds), aswell as evaluation of radius of gyration, DRMS (distance rootmean square) and RMSD (root mean square deviation) from areference structure. RMSD can also be computed with respectto the preceding structure, to highlight conformational changesalong the trajectory. More complex types of analysis include aprincipal component analysis (PCA) module, the computationof quasiharmonic entropy, several clustering algorithms and theevaluation of order parameters to monitor folding (i.e. numberof native and non-native contacts) and amyloid aggregation.3.1 PCA analysisPCA is a useful technique for separating large amplitudemotions from irrelevant fluctuations along an MD trajectory(Amadei et al., 1993). The approach consists of diagonalizingthe covariance matrix of atomic fluctuations with respect to theaverage structure, thus obtaining a set of eigenvectors andcorresponding eigenvalues which represent the axes of maximalvariance in the protein motion. An original feature in the PCAmodule is a procedure for iterative PCA evaluation alongtrajectory intervals of increasing length (e.g. snapshots 1–1000,snapshots 1–2000, etc.) to evaluate the significance of sampling.The eigenvectors computed at each time interval are projectedinto the ones of the preceding (i.e. shorter) time interval.A value close to 1 along the diagonal of the projection matrixindicates convergence of the corresponding eigenvector(s)suggesting that trajectory elongation is not necessary for the*To whom correspondence should be addressed. The Author 2007. Published by Oxford University Press. All rights reserved. For Permissions, please email: journals.permissions@oxfordjournals.org 2625calculation of those eigenvectors. The PCA module in Wordomcan also be used to generate pdb files where the last field of eachatom entry (which usually contains the crystallographic-factor) is set equal to the length of the vector defined bythe components of the selected eigenvector along the corres-ponding atom. This quantitative information on the proteinsegment(s) involved in the most relevant concerted motioncan be displayed directly with available graphics programs(Fig.1, left).3.2 Quasiharmonic entropyThe calculation of the quasiharmonic entropy is based on thediagonalization of the mass-weighted covariance matrix of theatomic fluctuations (Andricioaei and Karplus, 2001). The3n 6 non-zero eigenvalues i are related to the quasiharmonicfrequencies !i ¼ﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃkT=ip(where T is the simulation tempera-ture and k the Boltzmann constant), which in turn yield thevibrational entropy S ¼ kP3n6i h!i=kTe h!i=kT1 ln 1 e h!i=kT ,the quasiharmonic vibrational energy E ¼ P3n6i ð h!i2 þ h!ieh!i=kT1Þand the vibrational specific heat CV ¼ @E=@T. The entropymeasured in this way is directly related to the conformationalspace explored by the system along the MD simulation andprovides a wide range of applications from estimation of theentropic effects of mutations in the native state of proteins tothe determination of the relative entropy of ligands in a bindingpocket (Thorpe and Brooks, 2004).3.3 Clustering algorithmsClustering can be performed using two metrics (DRMS orRMSD) and three algorithms: hierarchical (Johnson, 1967),quality threshold (Heyer et al., 1999) or leader-like (Hartigan,1975). These algorithms have different CPU-time requirementsand yield different clustering performances, hierarchical andquality threshold being slower but more precise and leader-likebeing much faster but dependent on the frame order. For theformer, it is possible to first take into account a subset of thesnapshots and then cluster the remaining snapshots to reducememory and CPU-time requirements. In the case of RMSD-based clustering, superposition can be carried out beforecomputing the distance, depending on whether the relativeposition of the structures is relevant. It is also possible tosuperpose the structures according to a subset of atoms andcluster them relatively to another subset, which is useful foranalyzing multiple poses obtained by ligand docking into asingle protein structure.3.4 Order parametersWordom efficiently evaluates the number of contacts, i.e.C-pairs with a distance below a given threshold (input value).The subset of contacts which are present in the native structureis a useful variable to monitor the progress of (un)folding(Chan and Dill, 1998).The polar and nematic order parameters are used to studyamyloid formation (Cecchini et al., 2004). They are defined asP1 ¼ 1NPNi¼1 z^i  d^ and P2 ¼ 1NPNi¼132 ðz^i  d^Þ2  12 ; where d^ is aunit vector pointing in the preferred direction of alignment, z^i isa suitably defined molecular vector and N is the numberof peptides. These order parameters capture differentorientational properties of the system and yield useful andcomplementary information. The polar P1 describes thepolarity of the system, i.e. how much the molecular vectors(z^i) point in the same direction, and discriminates betweenparallel, antiparallel or mixed ordered aggregates. The nematicP2 describes the orientational order of the system anddiscriminates between ordered and disordered conformations.The evaluation of P1 and P2 facilitates the analysis ofsimulations of peptide aggregation (Cecchini et al., 2004) bymeasuring the degree of structural order during the fibrilformation process (Fig.1, right).4 CONCLUSIONSAnalyzing the ever growing amounts of MD trajectory data canbecome the rate-limiting step in MD simulation studies.The Wordom program facilitates this task. Manipulations ofmultiple trajectory files and simple geometric analyses requirea single line of input, while input files for clustering and PCAanalysis consist of about five lines. Therefore, for efficientlyextracting information from massive sets of MD data Wordomis better suited than other programs which require complicatedinput files and/or the use of a graphical interface. Moreover,thanks to low-level data access functions, adding new analysismodules requires only the coding of the core function allowingfor a steadily growing library of algorithms introduced by thecommunity of users.ACKNOWLEDGEMENTSThe authors thank F. Raimondi and F. Fanelli for Figure. 1,left and helpful discussions, A. Cavalli, P. Kolb and D. Westfor technical help and R. Friedman and F. Marchand forcomments to the manuscript. This work was supported bygrants from the SNF and NCCR-neuro to A.C.Fig. 1. (Left) Ribbon diagram of the -subunit of protein Gq coloredaccording to the length (red: large; green: intermediate; blue, small) ofthe atomic component of the first eigenvector on a 6-ns MD run by theWordom PCA tool and displayed by Pymol (DeLano, 2002). (Right)Time series of P1 (blue) and P2 (black) (Cecchini et al., 2004) along anMD trajectory of three Q15KLVFFA21 heptapeptides from theAlzheimer’s A-peptide. Insets show representative snapshots ofmixed and parallel -sheet aggregates corresponding to low and highvalues of P1, respectively.M.Seeber et al.2626Conflict of Interest: none declared.REFERENCESAmadei,A. et al. (1993) Essential dynamics of proteins. Proteins: Struct. Funct.Genet., 17, 412–425.Andricioaei,I. and Karplus,M. (2001) On the calculation of entropy fromcovariance matrices of the atomic fluctuations. J. Chem Phys., 115,6289–6292.Berman,H.M. et al. (2000) The protein data bank. Nucleic Acids Res., 28,235–242.Brooks,B.R. et al. (1983) CHARMM: a program for macromolecular energy,minimization, and dynamics calculations. J. Comput Chem., 4, 187–217.Cecchini,M. et al. (2004) Replica exchange molecular dynamics simulations ofamyloid peptide aggregation. J. Chem Phys., 121, 10748–10756.Chan,H.S. and Dill,K.A. (1998) Protein folding in the landscape perspective:Chevron plots and non-Arrhenius kinetics. Proteins: Struct. Funct.Bioinformatics, 30, 2–33.De Lano,W. (2002) The PyMOL Molecular Graphics System, San Carlos,CA, USA.Hartigan,J. (1975) Clustering Algorithms, New York, USA.Heyer,L.J. et al. (1999) Exploring expression data: identification and analysis ofcoexpressed Genes. Genome Res., 9, 1106–1115.Humphrey,W. et al. (1996) VMD–Visual molecular dynamics. J. Mol. Graph.Model., 14, 33–38.Johnson,S. (1967) Hierarchical clustering schemes. Psychometrika, 32, 241–254.Kale´,L. et al. (1999) Namd2: greater scalability for parallel molecular dynamics.J. Comput. Phys., 283–312.Rueda,M. et al. (2007) A consensus view of protein dynamics. Proc. Natl Acad.Sci. USA., 104, 796–801.Settanni,G. et al. (2005) -Value analysis by molecular dynamics simulations ofreversible folding. Proc. Natl Acad. Sci. USA, 102, 628–633.Thorpe,I. and Brooks ,C.L.III (2004) The coupling of structural fluctuations tohydride transfer in dihydrofolate reductase. Proteins: Struct. Funct.Bioinformatics, 57, 444–457.van der Spoel,D. et al. (2005) Gromacs: Fast, flexible, and free. J. Comput Chem.,26, 1701–1718.Wordom2627

Wordom: a program for efficient analysis of molecular dynamics simulations

Abstract

Similar works

Full text

Available Versions

RERO DOC Digital Library