Découverte :
Champs de force, énergie potentielle et autres joyeusetés en modélisation moléculaire

Introduction

Cet article, d'une très longue série, a pour but de donner plus de détails et peut-être - qui sait - vous titiller suffisamment pour enjamber la barrière physique (au sens littéral) et rejoindre le fabuleux monde d'Oz de la bioinformatique structurale.

Avant de commencer, je vous conseille de prendre du papier, un crayon, un thé ou un café. Ces breuvages serviront à illustrer mes exemples et faire une pause dans votre lecture.

La bioinformatique structurale, aussi connue sous le nom barbare de biologie structurale computationnelle, a pour objectif de comprendre les facteurs d'influence qui déterminent la fonction biologique. Quoi ? Vous faites ça tous les jours sur votre bureau et avec votre petit ordinateur ? Je vais vous faire déjanter tout de suite. Comprendre ces facteurs implique dès le départ d'avoir une bonne connaissance des éléments susceptibles de jouer un rôle. Ainsi et comme vous le savez très bien : fonction biologique = structure. Toutefois, on ne s'intéresse pas exclusivement à la structure finale. On va regarder comment cette structure se forme et pour cela, on peut regarder les acides aminés conservés au cours de l'évolution et responsables de certaines catégories de repliements structuraux ou encore comparer les molécules homologues et non-homologues ce qui apporte également beaucoup d'informations car elles renseignent sur la convergence structurale comme c'est le cas des protéines appartenant aux familles de repliements dont certains membres ont des scores de similarités médiocres...

Depuis longtemps, on sait que l'information du repliement est contenue dans la séquence protéique ou nucléique que l'on peut obtenir facilement aujourd'hui. La base de ces repliements réside dans les propriétés physico-chimiques liées aux interactions entre atomes. Le nombre de repliements semble donc infini aux vues des nombreuses conformations possibles. Heureusement pour nous, certains repliements agissent comme une ancre pour ne former qu'un nombre limité de conformations stables et leurs détections facilitent le travail de détermination de la fonction biologique.

Ainsi, la bioinformatique structurale utilise les modèles 1D (séquence) et 3D (structure). L'histoire pourrait s'arrêter là mais ce serait réducteur de penser que le travail se limite à faire un blast et un RMSD. En effet, la particularité de s'intéresser au Vivant est justement...qu'il est Vivant, donc modulable. Cette modularité au niveau structural a pour conséquence d'avoir des modules qui adoptent un repliement caractéristique indépendamment du reste de la séquence et de manière indépendante. De plus, certains repliements s'associent entre eux pour former un domaine structural fonctionnel mais parfois un module seul suffit pour permettre au repliement d'être fonctionnel. L'exemple le plus parlant est un ARN de transfert, dont les modules s'associent entre eux pour former un domaine spécifique d'un acide aminé : si bien que même si l'on reconnaît un module par une matrice poids-positions, cela n'indique en aucun cas son influence sur la structure 3D. La prédiction et la simulation sont donc des solutions alternatives à l'expérimentation tout azimut (qui est chère, longue et souvent hasardeuse dans sa réalisation).

La modélisation moléculaire englobe un grand nombre de techniques pour prédire la conformation d'un système moléculaire. Elles ont en commun de partir de fichiers de structures connues et disponibles dans les bases de données internationales comme RCSB PDB ou PDBe provenant d'expériences de Cristallographie RX, RMN ou CryoEM et de rechercher une nouvelle conformation. Très utilisé en Biologie Structurale, et notamment dans l'industrie pharmaceutique, cette approche permet également de mieux comprendre les nombreuses interactions les plus probables au sein d'une molécule (avec elle-même, c'est le cas lors de repliement secondaire de protéines) ou avec son environnement ou bien encore avec un ligand. Tout cela permet l'étude dynamique des molécules.

Cet article s'intéresse particulièrement à la mécanique moléculaire, il s'agit d'une des trois méthodes que l'on peut utiliser pour faire de la modélisation de structure protéique. Les deux autres impliquent soit d’utiliser la mécanique quantique soit d’utiliser un méthode mixte entre la mécanique moléculaire et la mécanique quantique. Dans les faits, la mécanique quantique s'adresse à des petits systèmes car les calculs sont complexes. La mécanique moléculaire s'adresse aux grands systèmes car les calculs sont moins complexes. La méthode mixte combine à la fois la mécanique moléculaire et quantique : la partie de la molécule à étudier (un site actif par exemple) est simulée via mécanique quantique tandis que le reste de la molécule est simulé par mécanique moléculaire.

La première étape est de bien comprendre ce que l’on modélise. Reprenons les bases et répétez avec moi : on modélise des molécules ! Bien...mais ce n'est pas fini ! Une molécule est faite d'atomes et chaque atome interagit avec son voisin. Au lieu de s'échanger de la farine ou du sel, des atomes proches mettent en commun des électrons et bim : ça fait une interaction. Et voilà, le premier nœud à problèmes : modéliser une interaction entre atomes. La solution : utiliser un champ de force.

 

La fête des atomes voisins !Champs de force

Pour bien comprendre ce qu'est un champ de force, je vais prendre un exemple très parlant. Faites un dessin hyper-simplifié de vous même debout : un point pour la tête, deux points pour les mains, un point pour une vertèbre thoracique et un autre pour une vertèbre lombaire, deux points pour vos fesses et enfin, deux points pour vos pieds. Reliez les points ensemble : chaque pied est relié à une fesse, vos fesses sont reliées à la vertèbre lombaire. Celle-ci est reliée à la vertèbre thoracique. Et la tête et chacune des mains sont reliées à cette vertèbre. Plutôt que de décrire toutes les caractéristiques internes (vos petits os, vos organes, vos cellules), on va partir du postulat que ces points suffisent amplement. Si vous avez une bonne notion du corps humain, vous avez sans doute remarqué que votre vertèbre thoracique est à mi-distance de chacune de vos mains et que cette vertèbre est plus proche de la tête que de vos fesses. Sinon, je vous invite à consulter.

Évidement, certains atomes présentent des caractéristiques propres comme leur charge (due à leurs électrons de valence) ou leur taille (due aux nombre de nucléons). Si on va plus loin, il y a également des interactions entre atomes très caractéristiques et qui jouent un rôle dans les forces qui s'exercent sur chacun des atomes liés. De même que des forces s'exercent entre atomes distants, on parle alors d'interactions entre atomes non liés.

Si on applique cela à la modélisation du corps humain, on peut dire que le point représentant la tête dispose de caractéristiques différentes de celui qui représente une main ou un pied. De même, on peut dire que la tête est liée à la vertèbre thoracique mais il existe des forces qui la lie de manière distante aux pieds. Ouf pour nous !

Comme vous vous en doutez, les échanges qui ont lieu au sein d'une interaction sont d'ordres énergétiques. Ainsi, les forces qui s'exercent sur les atomes ont une action papable sur des valeurs physiques réelles comme des distances ou des angles. Autrement dit, une force va apporter une énergie particulière pour modifier l'interaction elle-même : c'est l'énergie potentielle. Celle-ci, a, par définition, le potentiel d'être transformée en énergie cinétique donc en mouvement. Et la manière de représenter la variation est d'avoir recours à une fonction analytique. Pour les adeptes des mathématiques, il s'agit de la fonction d'une variable qui est développable en série au voisinage d'un intervalle donné. C'est-à-dire que sur un intervalle donné, la fonction générale est égale à la somme de fonctions simples.

Reprenons mon exemple. Votre schéma de base peut s'appliquer à tous les types d'humanoïdes . Toutefois, de petites différences sont notables : votre voisin sera plus grand ou plus petit que vous par exemple. On a donc des variables qu'il faut prendre en compte pour décrire les interactions entre chacun des points de votre schéma.

Pour résumer : un champ de force est la somme des forces (répulsions et attractions) qui s'exercent sur une interaction, autrement dit sur les atomes qui composent l'interaction. Chacune de ces forces se décrit par son énergie potentielle qui dépend de constantes décrites de manières empiriques et de variables qui dépendent du modèle structural (du schéma que vous avez fait).

Il existe de nombreux champs de force et chacun diffère par la valeur des constantes décrivant le type d'interaction qui lui est associé. Elles sont fournies par la spectroscopie infrarouge, la RMN ou par des calculs quantiques.

Pour repartir de mon exemple, chacun des traits que vous avez fait entre les points de votre schéma est un champ de force. Réjouissant !

Coordonnées cartésiennes et coordonnées internes

Maintenant vient la douloureuse question de la représentation. Doit-on représenter l'interaction ou les atomes qui la composent ? Vous allez me dire : c'est la même chose. En fait, non ! On peut soit décrire la molécule par la position selon les 3 coordonnées x, y et z en angström (10-12m) de chacun des atomes, on utilise alors des coordonnées cartésiennes, soit décrire la géométrie de la molécule en exprimant les longueurs de liaisons, les angles de valances et de torsion, on parle alors de coordonnées internes ou Z-matrice ou encore appelé modes de vibrations internes et que pour éviter la redondance des expressions matricielles, le premier atome est l'origine du référentiel.

On a vu qu'une énergie potentielle peut engendrer une modification d'un des paramètres de l'interaction. Ces modifications sont les coordonnées internes et sont réparties en deux grandes catégories : les vibrations d'allongements (appelé aussi vibrations de valence ou modes d'extension ou stretching) et les vibrations de déformations (appelé aussi mode de déformation ou bending). Parmi les stretching, il y a les vibrations symétriques (Vs) et asymétriques (Vas). Parmi les bending, il y a les déformations dans le plan de molécule : la rotation plane (notée β et appelée aussi rocking) et le cisaillement (noté δ et appelé aussi scissoring) ainsi que les déformations hors du plan de la molécule : le balancement (noté ω et appelé aussi wagging) et la torsion (notée τ et appelée aussi twisting).

Dans un système 3D, un objet peut être représenté selon sa position (x, y et z) et sa quantité de mouvement. C'est ce que l'on appelle les degrés de liberté. Ainsi, on a 3 degrés pour la position et 3 degrés pour la rotation chez les molécules linéaires (2 chez les linéaires, car un axe est bloqué) et un certain nombre de coordonnées internes. On peut associer une coordonnée interne à un degré de liberté et puisqu'un certain nombre a déjà été comptabilisé, on peut dénombrer les coordonnées internes pour chaque type de molécule. On a donc 3n-5 coordonnées internes (n étant le nombre d'atomes) pour une molécule linéaire et 3n-6 coordonnées internes pour une molécule non-linéaire. Par exemple, le dioxyde de carbone (CO2) est une molécule linéaire, une de ces rotations est bloquée par la nature des liaisons qui engagent le Carbone : elle n'a donc que 4 coordonnées internes. En revanche, le méthane (CH4) est une molécule non-linéaire, elle a 9 coordonnées internes.

Aujourd'hui, le format de fichier PDB - le plus utilisé - contient exclusivement des données sous forme de coordonnées cartésiennes. Vous pourrez retrouver la notice de ce format ici.

Un petit exemple tiré de la structure d'une protéase pancréatite bovine. J'ai volontairement extrait uniquement les coordonnées cartésiennes d'une proline de cette molécule.

A droite : le fichier pdb. A gauche : la représentation du fichier dans VMD

A droite : le fichier pdb. A gauche : la représentation du fichier dans VMD

Pour les coordonnées internes, on peut utiliser des logiciels (VDM) ou des packages R (Rpdb) ou encore Python (Bio.PDB). Finalement, j'ai utilisé VDM. J'ai obtenu toutes les distances entre atomes (7), chacun des angles (8) et des angles dièdres (6). Cela correspond à la géométrie de notre petite proline. Du coup, comme je vous ai parlé de 3n-6 coordonnées internes utiles, on a 21 coordonnées internes calculées par VMD soit au final, 6 valeurs redondantes. À l'échelle d'une petite molécule, ça fait peut mais quand on considère des macromolécules de plusieurs milliers d'atomes, ça change la donne... et notamment au niveau format de fichier et de la compression de données. Dans les faits, cela ne sert strictement à rien d'avoir la géométrie pour un acide aminé, il est plus intéressant d'obtenir les angles dièdres du squelette carboné d'une protéine (backbone carbone) et de les représenter via un diagramme de Ramachandran. Je reviendrai dessus pour un prochain article.

Calcul de l'énergie potentielle en mécanique moléculaire

En mécanique moléculaire, les noyaux et les électrons forment une seule entité : un point ponctuel. Ils sont tous sphériques, de charge nette (c'est-à-dire qu'elle est la même pour tous les atomes) et leurs tailles varient selon l'atome.

On a vu plus haut que l'énergie potentielle est en fait l'énergie nécessaire à apporter pour qu'une interaction soit modifiée dans un espace tri-dimensionnel. Un atome dans une molécule présente de nombreuses interactions : des interactions entre atomes liés et non-liés. Chacun de ces termes peut être représenté par une énergie potentielle. Au sein d’interactions entre atomes liés, on a une énergie potentielle de déformation de liaison (stretch), de déformation d'angles (bend) et de torsion d'angles dièdres (dihedral angle). Dans le cas d’interactions entre atomes non-liés, on a une énergie potentielle électrostatique (potentiel de Coulomb), des liaisons hydrogènes et des interactions de Van Der Waals (potentiel de Lennard-Jones). Notez que dans certains cas, une énergie potentielle pour les angles de torsion impropres peut être utilisées en plus !
À des fins de rapidité de calcul, il existe pour chaque champ de force des termes de couplage qui associent deux potentiels pour certaines interactions typiques comme par exemple : une élongation des liaisons quand l'angle diminue (stretch-bend) ou lorsque une conformation est éclipsée (stretch-torsion).

Au final, cela donne ces équations, avec en vert, les paramètres variables et en rouge, les paramètres empiriques.

Avec les indicateurs :

Ka : constance de force de liaison

l et l0 : distance interatomique mesurée et à l'équilibre

Kb : constante de force d'angle

θ et θ0 : angle mesuré et à l'équilibre

V : constante de force de torsion associé à n

n : ordre de séries de Fourier

Φ : angle mesuré dièdre

γ : angle de phase du dièdre

qi et qj : charges partielles de l'atome i et j

ε : constante diélectrique

rij : distance interatomique entre les atomes i et j

Aij : profondeur de puits du potentiel approché

r*ij : distance interatomique au minimum d'énergie

Et enfin la fameuse équation représentant l'énergie potentielle totale :

Le solvant est pris en compte par le choix de la valeur de la constante diélectrique.

Petit détail qui a son importance : il est exclu de vouloir comparer l'énergie potentielle de deux molécules même si elles présentent une forte similarité de conformation car de nombreux paramètres d'un champ de force à l'autre ne sont pas superposables. Enfin, sachez que selon le logiciel de calcul de champs de force utilisé, les unités utilisées ne sont pas les mêmes : en kcal/mol/A pour Charmm ou en kJ/mol/nm pour GROMOS.

Si vous êtes encore vivants et curieux d'en savoir un peu plus, je vous invite au prochain article sur l'analyse des conformations et de la minimisation. Par ailleurs, pour les aficionados du code et autres adeptes de la gâchette du clavier, sachez que BioPython dispose d'un ensemble de packages que l'on utilise pour l'étude des structures PDB : Bio.PDB. De toute manière, les principaux logiciels de visualisation de structures sont en Python...

Voilà, c'est tout pour le moment. Vous venez de voir comment fonctionnent globalement les programmes et les logiciels de champs de force en modélisation moléculaire. La prochaine fois, je vous montrerai comment on peut s'en servir pour rechercher une conformation particulière.

  • À propos de
  • Issu du master de Biologie Structurale Intégrative et Bioinformatique de Strasbourg, je suis un passionné des sciences qui s'intéresse à tous les sujets pouvant permettre de mieux comprendre les mécanismes du Vivant. Cette curiosité m'a poussé à devenir technicien de laboratoire puis bioinformaticien dans les domaines de la biologie structurale et de l'analyse NGS.

6 commentaires sur “Champs de force, énergie potentielle et autres joyeusetés en modélisation moléculaire

  1. Très bon article ! enfin de la modélisation moléculaire sur bioinfo-fr :).

    Concernant les champs de forces, je pense qu\'il est intéressant de préciser les plus classiques :
    - GROMOS (J. Comp. Chem., 2004)
    - AMBER (J. Am. Chem. Soc., 1995)
    - OPLS (J. Phys. Chem. B., 2001)
    - CHARMM (J. Phys. Chem. B., 1998)

    A ne pas confondre avec les programme de modélisation moléculaire (Gromacs, Amber, Charmm) qui ont souvent le même nom que les champs de force.

    Pour ma part, je vous conseille GROMACS comme programme pour faire de la MD, très simple d\'utilisation et il fonctionne avec tous les champs de forces classiques.

    • De la... MD ? C\'est un raccourci pour modélisation 😮 ?

      • Ah oui pardon, j\'ai trop le nez dedans. MD c\'est pour Molecular Dynamics (dynamique moléculaire).

  2. Article très clair et intéressant !

    Par ailleurs, suis-je la seule à ne plus pouvoir voir les équations dans la partie \"Calcul de l\'énergie potentielle en mécanique moléculaire\" ?

    • Hmmm, en effet, il y a un soucis d\'affichage.
      On va voir ce qu\'on peut faire pour corriger ça.

    • Ça semble venir d\'un problème d\'accès au site latex.codecogs.com, on ne va pas pouvoir corriger a, désolé :/
      Si le site continue de dysfonctionner, on demandera à l\'auteur d\'en trouver un autre ou bien d\'insérer un screenshot.

      Merci du retours!

Laisser un commentaire