Structural Bioinformatics Library
Template C++ / Python API for developping structural bioinformatics applications.
User Manual

Molecular_covalent_structure

Authors: F. Cazals and A. Chevallier and T. Dreyfus

Introduction

Describing a molecule: topology, geometry, biochemistry

The description of a molecule encompasses several aspects:

  • topology defined as a molecular graph: the vertices of the graph are the particles namely atoms or pseudo-atoms for coarse grain models; the edges of the graph are the covalent bonds connecting the particles.

    Note that the graph may not be connected, which corresponds to several molecules.

    See also Biomolecules and the PDB format for a discussion of topology issues with respect to molecular data formats.

    The topological representation of a molecule allows one to:

    • iterate on molecules,
    • iterate on particles, bonds, bond angles and torsion angles,
    • report tuples (pairs, triples, quadruples) used to compute potential energies,

    The topology is handled by the class SBL::CSB::T_Molecular_covalent_structure

  • geometry: the conformation i.e. the embedding in three dimensions of the molecular graph. The geometric representation of a molecule allows one to convert back and forth between Cartesian and internal coordinate representations, and to compute partial derivatives of potential energies wrt coordinates when performing energy minimization. The default representation is in cartesian coordinates with the class SBL::Models::T_Geometric_conformation_traits::Conformation . These coordinates can be converted into internal coordinates as explained in the package Molecular_coordinates

  • biochemistry: the properties qualifying the particles, from which potential energies are computed, see the package Molecular_potential_energy. It is handled by the concept ParticleInfo , as detailed below.

Amongst these three aspects, this package solely handles the topological ones.

Covalent data structures: default and optimized

When the covalent structure is used to compute the potential energy and the gradient of the energy of a large number of conformations, it becomes necessary to optimize its traversal – and the calculation of internal coordinates. The base data structure uses a graph for representing it, so that it is easy to visit particles and bonds. However, the following comments are in order:

  • visiting bond angles and torsion angles require some calculations, since they are not stored within the graph data structure;
  • supplemental filters from the parameters used for computing the energy may reduce the number of particles, bonds or angles to visit : this is particularly the case when 1-3 interactions are ignored, or when a force constant is considered null for a given bond or angle type.

For these reasons, two covalent data structures are provided:

File formats

Practically, the following file formats are used:

  • The PDB format: the reference format for biomolecules whose structure has been resolved and are deposited in the Protein Data Bank ( PDB website).

  • For molecules other than peptides and proteins, we use the MOL format – see section Beyond the PDB format and MOL format.

Implementation: default covalent structure – with graph representation

The first implementation provided uses internally a graph data structure to code particles and bonds, which comes with the following pros and cons:

  • (pro) The graph comes with topological operations
  • (cons) Access to the data members (internal coordinates in particular) are not optimized.

Three important data structures are provided by this package :

ParticleInfo concept

This concept is used to represent a particle in the covalent structure and to attach information to it.

Depending on the context, the relevant information and the necessary operations on a particle vary.

For example, the particle could belong to a small molecule or to a polypeptidic chain.

In this latter case, one may want to retrieve a particular atom (e.g. the Calpha carbon) – an irrelevant request if the molecule is not a polypeptide chain.

For this reason, the class SBL::CSB::T_Particle_info_traits< ParticleInfo > provides all the necessary types and operations depending on the model of ParticleInfo . Such a mechanism uses the C++ Partial Template Specialization : any structure can be used as a model of the concept ParticleInfo provided that the class SBL::CSB::T_Particle_info_traits is redefined for this particular model. The required operations depend on the context where it is used :

Two models of ParticleInfo are provided :

  • As a simple structure, a std::string can be used for representing the information of a particle; this can be done by including the file "SBL/CSB/Particle_info_string.hpp"
Note that some atoms in the covalent structure may have no equivalent in the ESBTL data structure, e.g. because there is no representation of this atom in the loaded file representing the protein. This happens typically for missing hydrogen atoms, or portions of long side chains unresolved in a crystal structure. To be able to identify an atom from the particle info data structure even without the ESBTL data structure, the chain identifier, residue sequence number, insertion code and atom name are also directly stored within the SBL::CSB::T_Particle_info_for_proteins . structure.


: Comment on the following case: a system containing a protein and a drug is loaded. which ParticleInfo should we use? That case is not covered by remark 1.


The properties required for the class used are twofold:
  • comparable: this ensures that each instance of ParticleInfo uniquely identifies a particle;

  • streamable: this ensures that each particle in the covalent structure has an associated name, which is e.g. used when printing the covalent structure as a graph using Graphviz .

The concept can be represented by a simple string representing the name of the particle. When used for proteins, the name is built from the chain identifier, the residue sequence number, the insertion code, and the atom name. In doing so, it is possible to find particles in the covalent structure from biochemical information.
A more complete data structure is provided as SBL::CSB::Particle_info_for_proteins and defines three features :
  • it stores separatly all the information required for identifying the particle; this allows to use the covalent structure in the twin package Molecular_potential_energy where access to the biochemical information is required for associating the force field parameters;


Representing covalent structures: main data structure

Basically, the class SBL::CSB::T_Molecular_covalent_structure is a wrapper for a boost graph with accessors, modifiers and iterators over the graph. More precisely, it uses the boost subgraph data structure to offer the possibility of clustering the covalent structure, e.g by polypeptidic chains or by residues. There is a unique template parameter ParticleInfo representing the information attached to a particle, that is a simple string by default. An object of type ParticleInfo should be ordered, streamable and default constructible. Four main functionnality are available.

Building the covalent structure.

Adding a particle to the covalent structure is done using the method SBL::CSB::T_Molecular_covalent_structure::add_particle : this method takes a particle info representing that particle.

Adding a bond is done using the method SBL::CSB::T_Molecular_covalent_structure::add_bond : this method takes two vertices of the graph, and adds, if possible, an edge between the two vertices. It also takes optionnaly a bond type (simple by default, or double) to represent the covalence of the bond.

Establishing a correspondence with an embedded geometry specified by cartesian coordinates.

It is possible to map a geometric model over the covalent structure by simply mapping the vertices of the graph to the geometric representation of the particles.

A vertex is termed embedded if it is endowed with Cartesian coordinates. A covalent structure is termed embedded if all its vertices are embedded, and partially embedded otherwise.


The notion of partial embedding is especially important when loading polypeptide chains from PDB files, which in general do not contain coordinates for all atoms (hydrogen atoms, heavy atoms located in flexible regions).

In the SBL, a geometric model is represented by a conformation that is a D-dimensional point consisting of the cartesian coordinates of the particles of a molecule. By default, the ith inserted vertex is mapped to the ith particle in the conformation.

Given an embedded vertex, it is then possible to access x, y and z coordinates from the covalent structure from the methods SBL::CSB::T_Molecular_covalent_structure::get_x, SBL::CSB::T_Molecular_covalent_structure::get_y and SBL::CSB::T_Molecular_covalent_structure::get_z .

Iterating over the covalent structure.

It is possible to iterate over the particles, bonds, bond angles and torsion angles of the covalent structure. There are two modes of iteration : over all the entities, or over all the entities that are mapped with a geometric model. For example, if ones want to compute internal coordinates, it is only possible to do so from particles that are mapped. To switch between the modes, each iteration method has as argument an optional boolean tag (by default, it iterates over mapped particles).

Consider a set of particles associated with one or several connected components of the covalent structure. The subgraph induced by these particles is returned by SBL::CSB::T_Molecular_covalent_structure::get_sub_structure . The type of this graph is SBL::CSB::T_Molecular_covalent_structure . Note that stored information are not duplicated : modifying the sub-structure modifies also the parent covalent structure, and modifying the parent covalent structure modifies the sub-structure.

As an application, one can iterate on the individual molecules of a molecular system, as these correspond to the connected components of the covalent structure. To do so, one proceeds as follows :

Displaying the covalent structure.

The covalent structure can be dumped in a file using Graphviz through the method SBL::CSB::T_Molecular_covalent_structure::print . For large systems (e.g. large proteins), the high number of particles, the Graphviz software may take a time to create the output. Also note that the neato command from Graphviz has a better rendering than the more classical dot command in this case.

Building covalent structures

Building a covalent structure requires three pieces of information:

  • Attributes for vertices (type of vertices).

  • Attributes for edges (type of covalent bonds).

  • The topology itself. For example, when building the subgraph associated with an amino-acid, the topology of this amino-acid is required. this information can be implicit as in the case of amnino-acids and polypeptide chains, or can be explicit i.e. provided in the input file. This accounts for the two types of builders provided thereafter.

In the context of potential energy calculations – see Molecular_potential_energy, molecular dynamics software (e.g Gromacs) usually use parameter files including the topology of those molecules (e.g itp file format). We deliberately choose to separate the topology from those parameters, allowing to perform independently analysis and calculations over the covalent structure. In particular, the covalent structure is represented by a graph that can be used as input of combinatorial packages in CADS such as Betti_numbers or Earth_mover_distance


Implicit builders for polypeptide chains. The class SBL::CSB::T_Molecular_covalent_structure_builder_for_proteins is a functor building the covalent structure from an arbitrary data structure SBL::CSB::T_Molecular_covalent_structure_builder_for_proteins::Structure defining the polypeptidic chains and their residues. More precisely, a protein is represented by a mapping from a chain identifier to the chain itself, each chain being represented by a list of pairs (residue name, residue id). The graph of the covalent structure can be constructed at different scales : (0) with all atoms from all residues, (1) with only the heavy atoms, (2) with only the carbon alpha of each residue, or (3) with pseudo-atoms from the Martini coarse grain model. The scale factor is determined at the construction of the builder and is 0 by default. The only template parameter Covalentstructure is the covalent structure type.

The nomemclature for atoms in PDB files can be found here: IMGT / Amino acids: 3D model and atoms nomenclature. Should the PDB parser find an atom whose naming is unknown, a warning is issued.


For sequence alignement purposes, residues in PDB files are identified by their sequence number and their insertion code. In particular, sequence numbers are not necessarily successive, and several residues with the same sequence number and different insertion codes may coexist. For this reason, in the builder, residues are sorted by the atomic serial number of their first atom. In addition, a geometric predicate is used to check if the atom C of a residue makes a bond with the atom N of the next residue in the sequence. Practically, it is checked that the distance between these atoms is less than 2 $\AA$. If it is not the case, no bond is created in the covalent structure and the covalent structure has multiple connected components.


Explicit builder for other molecules.

Since there is no particular difficulty nor operation for loading molecules from the MOL format, building the covalent structure from a MOL file is directly done in the loader SBL::CSB::T_Molecular_covalent_structure_loader_from_MOL .

Loading atomic resolution covalent structures

Loader from PDB files.

The class SBL::IO::T_Molecular_covalent_structure_loader_from_PDB is a loader as defined in the Module_base package : it loads an input PDB , uses a builder to build a covalent structure, and maps the vertices of the covalent structure to the loaded atoms by order of their atomic sequence identifier.

The following comments are in order:

Specific atomic naming conventions.

  • Some atom names may have different conventions. However, it is important for further calculations (e.g potential energy calculation) to have a unique convention for each atom name. When required, we fixed the convention for the naming to be the one used within force field parameters. For example, hydrogen atoms bonded to the first carbon in the N-terminal cap called ACE, or the C-terminal cap called NME can be names XHH3 or HH3X, where X is replaced by 1, 2 or 3. In Amber force fields, it is named HH3X. To handle such cases, we fixed our convention on HH3X. This label is set by the operator Aliases defined in a specialization of the class SBL::CSB::T_Particle_info_traits .

Loader from MOL files.

The class SBL::IO::T_Molecular_covalent_structure_loader_from_MOL is a loader as defined in the Module_base package : it simply loads an input MOL file, and directly builds the covalent structure from the loaded data and maps the vertices to the laoded particles. It has a unique template parameter MolecularCovalentStructure , that is the representation of the covalent structure, by default SBL::CSB::T_Molecular_covalent_structure <> (see package Molecular_covalent_structure)

Loading coarse grain covalent structures

It is possible to set the coarse level used by the builder to create the covalent structure. This class has four template parameters, all having a default type :

Implementation: advanced features

Canonical representations

Manipulating molecular coordinates (package Molecular_coordinates) or computing potential energies (package Molecular_potential_energy) benefits from a canonical representation of internal coordinates.

The method SBL::CSB::T_Molecular_covalent_structure::get_canonical_rep returns the cannonized representation of a bond, a bond angle or a dihedral angle.


Canonical representation of a bond. Such a representation is obtained by ordering its two vertices – that is the particle_rep used to represent the vertices.


Canonical representation of a valence angle. Such a representation is obtained by ordering the first and third vertices – since the middle one is imposed.


Canonical representation of a dihedral angle. This case covers two sub-cases (see also the package Molecular_coordinates):

  • Proper angle: second and third vertices sorted;
  • Improper angle: central atom first, remaining three sorted.

Implementation: advanced features

The class SBL::CSB::T_Molecular_covalent_structure_optimal< ParticleInfo > is an optimized version of the covalent structure data structure. In particular, all entities (particles, bonds and angles) that are used in the computation of the energy are stored in a vector, ensuring a constant time access each time.

In order to use this class, once has to first build a covalent structure with the class SBL::CSB::T_Molecular_covalent_structure, then use the class SBL::CSB::T_Molecular_potential_energy_structure_parameters_traits_optimal to build the optimized version of the covalent structure data structure. See the package Molecular_potential_energy for examples and more details.

The functionalities depending on the graph structure which are not used in the potential energy calculation are unavailable with the optimized data structure. This includes using the graph data structure itself, dynamically changing the graph (adding particles or bonds), and printing the covalent structure in the .dot format.


Examples

The following example shows how to load a PDB file and build the corresponding covalent structure at a given scale. It then output a dot file in the Graphviz format to check the covalent structure.

#include <iostream>
#include <sstream>
#include <iterator>
#include <algorithm>
#include <SBL/Models/Atom_with_hierarchical_info_traits.hpp>
#include <SBL/CSB/Particle_info_for_proteins.hpp>
#include <SBL/CSB/Molecular_covalent_structure.hpp>
#include <SBL/CSB/Molecular_covalent_structure_builder_for_proteins.hpp>
#include <SBL/IO/Molecular_covalent_structure_loader_from_PDB.hpp>
int main(int argc, char *argv[])
{
if(argc < 2)
return -1;
//Loads the covalent structure from a PDB file.
loader.add_input_file_name(argv[1]);
if(argc > 2)
loader.set_coarse_level(atoi(argv[2]));
loader.set_loaded_water(false);//(true);
bool b = loader.load(true, std::cout);
if(!b)
return -1;
std::ofstream out("structure.dot");
S.print(out);
out.close();
std::cout << "Number of bonds (theoretical vs traversal) : " << boost::num_edges(S.get_graph()) << " " << std::distance(S.bonds_begin(true), S.bonds_end(true)) << "/" << std::distance(S.bonds_begin(false), S.bonds_end(false)) << std::endl;
//Compute the theoretical number of bond angles.
unsigned nb = 0;
for(Covalent_structure::Particles_iterator it = S.particles_begin(false); it != S.particles_end(false); it++)
{
unsigned n = boost::degree(*it, S.get_graph());
if(n > 1)
nb += n*(n-1)/2;
}
std::cout << "Number of bond angles (theoretical vs traversal) : " << nb << " " << std::distance(S.bond_angles_begin(true), S.bond_angles_end(true)) << "/" << std::distance(S.bond_angles_begin(false), S.bond_angles_end(false)) << std::endl;
//Compute the theoretical number of proper dihedrals.
nb = 0;
for(Covalent_structure::Bonds_iterator it = S.bonds_begin(false); it != S.bonds_end(false); it++)
{
unsigned n1 = boost::degree(boost::source(*it, S.get_graph()), S.get_graph()) - 1;
unsigned n2 = boost::degree(boost::target(*it, S.get_graph()), S.get_graph()) - 1;
if(n1 > 0 && n2 > 0)
nb += n1*n2;
}
std::cout << "Number of proper dihedrals (theoretical vs traversal) : " << nb << " " << std::distance(S.dihedral_angles_begin(true), S.dihedral_angles_end(true)) << "/" << std::distance(S.dihedral_angles_begin(false), S.dihedral_angles_end(false)) << std::endl;
//Compute the theoretical number of improper angles.
nb = 0;
for(Covalent_structure::Particles_iterator it = S.particles_begin(false); it != S.particles_end(false); it++)
{
unsigned n = boost::degree(*it, S.get_graph());
if(n > 2)
nb += n*(n-1)*(n-2)/6;
}
std::vector<Covalent_structure::Particle_rep> ccs;
S.get_molecules(std::back_inserter(ccs));
std::cout << "Number of molecules : " << ccs.size() << std::endl;
for(unsigned i = 0; i < ccs.size(); i++)
{
S.get_molecule(ccs[i], mol);
std::cout << "Molecule " << (i+1) << " (particles / bonds) : " << std::distance(mol.particles_begin(true), mol.particles_end(true)) << "/" << std::distance(mol.particles_begin(false), mol.particles_end(false)) << " / " << std::distance(mol.bonds_begin(true), mol.bonds_end(true)) << "/" << std::distance(mol.bonds_begin(false), mol.bonds_end(false)) << std::endl;
std::cout << "First atom: " << mol.get_graph()[p].atom_name << " " << mol.get_graph()[p].res_name << " " << mol.get_graph()[p].res_id << mol.get_graph()[p].ins_code << " " << mol.get_graph()[p].chain_id << std::endl;
next++;
while(next != mol.particles_end(false))
{
jt++;
next++;
}
p = *jt;
std::cout << "Last atom: " << mol.get_graph()[p].atom_name << " " << mol.get_graph()[p].res_name << " " << mol.get_graph()[p].res_id << mol.get_graph()[p].ins_code << " " << mol.get_graph()[p].chain_id << std::endl;
}
//Special feature : look for a particle from info :
std::pair<bool, Covalent_structure::Particle_rep> pp = Particle_info_traits::Finder()(S, 'A', 1);
return 0;
}
Builds a covalent structure from a given structure.
Definition: Molecular_covalent_structure_builder_for_proteins.hpp:85
Self & get_molecule(Particle_rep p, Self &S)
Builds a new covalent structure containing only the particles of the molecule containing p;.
Definition: Molecular_covalent_structure.hpp:1520
Covalent_structure_graph_traits::vertex_descriptor Particle_rep
Definition: Molecular_covalent_structure.hpp:146
Bond_angles_iterator bond_angles_end(bool is_embedded=true) const
Iterate over the bond angles (u, v, w) so that u and v are the largest indices lower than w such that...
Definition: Molecular_covalent_structure.hpp:1428
Particles_iterator particles_begin(bool is_embedded=true) const
Definition: Molecular_covalent_structure.hpp:1379
Dihedral_angles_iterator dihedral_angles_end(bool is_embedded=true) const
Iterate over the torsion angles (u, v, w, t) so that u, v and w are the largest indices lower than t ...
Definition: Molecular_covalent_structure.hpp:1442
boost::transform_iterator< Make_valid, Bonds_iterator_filtered > Bonds_iterator
Definition: Molecular_covalent_structure.hpp:165
const Covalent_structure_graph & get_graph(void) const
Definition: Molecular_covalent_structure.hpp:985
Particles_iterator particles_end(bool is_embedded=true) const
Definition: Molecular_covalent_structure.hpp:1386
OutputIterator get_molecules(OutputIterator out) const
Computes the connected components of this graph and returns one vertex for each connected component.
Definition: Molecular_covalent_structure.hpp:1464
Bonds_iterator bonds_end(bool is_embedded=true) const
Definition: Molecular_covalent_structure.hpp:1400
boost::filter_iterator< Is_embedded, Particles_iterator_base > Particles_iterator
Definition: Molecular_covalent_structure.hpp:160
Bonds_iterator bonds_begin(bool is_embedded=true) const
Definition: Molecular_covalent_structure.hpp:1393
Dihedral_angles_iterator dihedral_angles_begin(bool is_embedded=true) const
Iterate over the torsion angles (u, v, w, t) so that u, v and w are the largest indices lower than t ...
Definition: Molecular_covalent_structure.hpp:1435
std::ostream & print(std::ostream &out) const
Definition: Molecular_covalent_structure.hpp:1529
Bond_angles_iterator bond_angles_begin(bool is_embedded=true) const
Iterate over the bond angles (u, v, w) so that u and v are the largest indices lower than w such that...
Definition: Molecular_covalent_structure.hpp:1421
Loader for covalent structures from PDB files using ESBTL.
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:119
static void set_allow_incomplete_chains(bool with_incomplete_chains)
Set the load incomplete residues tag.
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:496
const MolecularCovalentStructure & get_covalent_structure(unsigned i) const
ith covalent structure (const).
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:553
static void set_coarse_level(unsigned coarse_level)
Set the coarse level.
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:463
static void add_input_file_name(const std::string &file_name)
Add manually the name of an input PDB file.
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:454
bool load(unsigned verbose=false, std::ostream &out=std::cout)
Loads the data.
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:799
static void set_loaded_water(bool with_water)
Set the loaded water tag.
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:480
static void set_without_disulfide_bonds(bool without_ss)
Set the without disulfide bonds tag.
Definition: Molecular_covalent_structure_loader_from_PDB.hpp:504
Concept for manipulating particle info attached to a covalent structure.
Definition: Particle_info_traits.hpp:70
Traits class defining atoms traits (biophysical and geometric properties). Traits class defining atom...
Definition: Atom_with_hierarchical_info_traits.hpp:198
p
Definition: extract-wales-graph.py:74
n
Definition: extract-wales-graph.py:70