Gromologist: a Gromacs-Oriented Utility Library for Structure and Topology Manipulation
Despite the increasing automation of workflows for the preparation of systems for molecular dynamics simulations, the custom editing of molecular topologies to accommodate non-standard modifications remains a daunting task even for experienced users. To alleviate this issue, we created Gromologist, a utility library that provides the simulation community with a toolbox of primitive operations, as well as useful repetitive procedures identified during years of research. The library has been developed in response to users' feedback, and will continue to grow to include more use cases, thorough automatic testing and support for a broader spectrum of rare features. The program is available at gitlab.com/KomBioMol/gromologist and via Python's pip.
Introduction<!>Features Topology manipulation<!>Structure manipulation<!>2/5<!>Simultaneous processing of topology and structure<!>Checking and printing<!>Dependencies<!>Outlook
<p>Molecular dynamics (MD) simulations are a cornerstone of modern-day computational chemistry and biology, and in the last 40 years a range of general-purpose codes have been developed to make simulations more reproducible and accessible to the community [1][2][3][4][5][6][7]. Although the preparation of typical biomacromolecular systems, e.g. proteins, lipid membranes or nucleic acids containing standard residues, has become extremely streamlined and automated [8][9][10], the corresponding topology files -containing a full description of the system's connectivity and energetics -can be unreadable and not amenable to simple manual manipulation. As a result, any non-automated departure from these standards, such as removing atoms or adding bonds, can become a painstakingly complex process, in particular for novice users not aware of nuances and caveats of these fixed file formats.</p><p>To deal with these routine but non-standard tasks, high-level utility libraries are needed that can accommodate the increasing complexity of inputs and dynamically respond to popular requests from the community. One notable example of such a library, Parmed, has already been incorporated into Amber and contains many OpenMM-specific functions [11], but its primary goal is to provide easy access to general properties and allow for smooth interconversion between standard file formats. Another one, pmx, serves to prepare hybrid topologies for "alchemical" simulations and the analysis of thus obtained results [12]. Many routine tasks, though, are not covered by any popular library, leaving plenty of room for new developments.</p><p>In an attempt to fill that gap, we present Gromologist, a Pythonic utility library oriented at custom manipulation of Gromacs topologies. For many years now, Gromacs has been not only the most broadly used but also one of the fastest-growing MD engines in the field of molecular simulations [13,1] and has a well-established user base, whichcombined with the ease of format interconversion allowed by other tools -makes Gromologist a broadly applicable tool for all types of non-standard operations including, but not limited to, adding and removing atoms and bonds, introducing amino acid mutations in both the structure and topology, checking and listing force field parameter compatibilities, looking for mismatches between structures and topologies, merging files, manipulating alchemical states, modifying force field types or combination rules, parameter optimization or automated editing of structures. In addition, we maintain a strict zero-dependency policy (except for a few very specific features), so that the library is conveniently lightweight and can be easily installed with any modern Python distributions without causing dependency conflicts.</p><!><p>To handle Gromacs topologies, Gromologist implements a hierarchical representation of sections, subsections and entries. Sections represent high-level abstractions, such as individual molecule definitions, force field parameter sets or headers/footers. Subsections correspond to individual fragments of the topology file, delimited with the [ subsection name ] syntax, and directly handle most operations. Entries correspond to individual lines and provide a low-level access to all key features of the contents of a topology.</p><p>Splitting and merging. By default, topologies are saved as a single standalone file, including all .itp files referenced therein, but splitting to individual .itp is also supported. If needed, files can be made lightweight by keeping (a) only molecules listed in the [ system ] subsection and/or (b) only parameters used by the system.</p><p>Bonds and atoms. Primitive operations on topologies include adding, removing or swapping atoms in molecules while preserving a correct numbering scheme across all sections. Adding or removing bonds is also supported, including new bonds between two molecules, in which case the two become merged into a single molecule entry, and all new bonded terms (pairs, angles, dihedrals) are automatically identified and added, while existing ones are preserved and renumbered correctly. This allows for easy creation of chemical adducts or cyclic molecules, as well as editing of protonation states or simple chemical modifications without the need of extending the basic residue library (.rtp/.hdb files). Another heavily requested feature was the automated introduction of ,,special bonds" -disulfides and transition metal coordination bonds -that are not always correctly guessed by default Gromacs tools. Finally, Gromologist also provides a low-level interface for the mutation library, and can be easily incorporated into custom modification workflows.</p><p>Mutations. In this vein, a mutations module allows to introduce standard amino acid mutations to the topology without the need to go through the full process of topology generation. This has proven particularly useful in the study of heavily glycosylated systems, such as the SARS-CoV-2 Spike [14], where ready-made topologies and systems were often shared publicly but introducing even minor changes (such as a point mutation) would entail going through the same complicated server-based pipeline, risking introducing additional errors at this stage. This idea is also extended to protonation states, so that individual titratable residues can be protonated or deprotonated easily and in a consistent manner. As long as standard atom name conventions are followed, the mutation module is force-field agnostic, and the user can point to a specific .rtp file containing the respective residue definitions.</p><p>Alchemistry. Another set of utilities is related to the management of alchemical states. Gromologist can add or remove selected alchemical states, or swap them -a useful feature e.g. for the equilibration of B-states that helps avoid the computational overhead associated with the use of the free energy code.</p><p>Force fields. Finally, some features are meant to facilitate force field development and the implementation of selected methods such as solute tempering (REST2) [15]. These include the possibility of automatically cloning types along with all their bonded interaction parameters and renaming all selected atoms to the new type, or by automatically adding NBFIX-type modifications just by specifying the deviation from the standard Lorentz-Berthelot rules [16,17]. An additional module implements a generic strategy for multiple dihedral fitting to QM data, akin to that available in the FFTK module of VMD for NAMD topologies [18].</p><!><p>The structure module can deal with both PDB and GRO files, interconverting between them when necessary. Common operations include constructing new atoms from existing ones, interpolating between structures, atom/residue renumbering according to custom rules or generation of virtual sites, e.g. upgrading 3-point water molecules to 4-point ones.</p><!><p>In addition, many often-required features of the PDB files can be filled in automatically. Chain assignment, an information typically lost in PDB to GRO conversion can be inferred based on a simple distance criterion, similarly to the CONECT entries required by some analysis programs. Elements can be inferred from atom names, and beta-factors can be set to reflect external data sets, including an option to smooth the values out spatially.</p><p>To assist with these operations, Gromologist implements a robust selection language akin to that of VMD [19], allowing for logical operations, selections "within X of Y", and a number of predefined molecule classes such as "protein", "nucleic" or "solvent". Periodic boundary condition (PBC) treatment is also implemented for any generalized triclinic box.</p><!><p>In Gromologist, topologies can be matched with a compatible structure, so that where relevant, topology operations can be simultaneously performed on the structure attribute (and vice versa), facilitating the work with simulation-ready systems, and fully leveraging the information available in both files. Examples of such procedures include atom addition/removal, or introduction of mutations. When a topology is available, chain assignment can be made using molecule definitions therein. Moreover, when Gromacs often only indicates an unspecified mismatch between e.g. the number of atoms in structure and topology, Gromologist provides a specific list of all unmatched atom names, facilitating the identification of missing atoms or residues in either one.</p><!><p>The last area of applications of Gromologist covers extracting information about the topologies and structures, or making it better visible to the user. On the most basic level, one can list molecules in the system, atoms in the molecule, as well as selected bonded terms (bonds, 1-4 pairs, angles, dihedrals) using atom names or atom types. Charges and masses of both the whole system and individual molecules are easy-to-access attributes. For a protein or a nucleic acid structure, the sequence can be printed chain by chain, and per-residue missing atoms identified, helping identify cases in which residues have to be rebuilt. Gradually diagnostic functions are being incorporated to quickly identify structural issues prior to simulation, like one currently allowing for validation of all chiral centres in a protein.</p><p>When required, atom names can be explicitly stated in the comments of bonded interaction entries to facilitate debugging, and numerical values of the parameters can be explicitly included in the topology. In the same vein, fields set using the define syntax can be set to their explicit values. Parameters whose numerical definitions are missing can also be identified with a single function call.</p><!><p>To be applicable as broadly as possible, Gromologist uses no dependencies other than base Python libraries. Although a few optional features, such as smoothing of beta-factors or dihedral optimization, do require numpy and scipy, these make up extremely rare use cases and therefore do not form part of the official dependency list.</p><!><p>Given the existing set of capabilities, Gromologist is becoming an easily extensible platform for the introduction of further convenience functions, editing tools and workflows. Further development of the platform could benefit the Gromacs community by rapidly implementing utilities tailored to new or popular protocols, and a closer collaboration with the Gromacs development team is expected to directly address these issues. Moreover, extensive testing and some code reorganization will be required to ensure that the library covers all features supported by Gromacs, including more cryptic ones. In the long run, however, we expect Gromologist to become a stable element of the Gromacs environment and a valuable contribution to the broad simulation community.</p><p>3/5</p>