Un flux de travail computationnel pour la priorisation des gènes de l'hôte associés aux métabolites microbiens dans le syndrome de l'intestin irritable à prédominance constipation

Résumé

Ce protocole intègre la prédiction ciblée des métabolites microbiens, la transcriptomique de la muqueuse rectale, l'interaction protéine-protéine et l'enrichissement des voies métaboliques, le dockage moléculaire, la simulation de dynamique moléculaire et l'estimation de l'énergie libre de liaison par la méthode de mécanique moléculaire/surface de Poisson-Boltzmann (MM-PBSA), afin de produire une liste restreinte et classée de gènes hôtes candidats associés aux métabolites, ainsi que des complexes protéine-ligand priorisés structurellement, destinés à des validations expérimentales.

Résumé

Aucun pipeline computationnel normalisé n'existe pour hiérarchiser systématiquement les gènes de l'hôte et les complexes protéine-ligand associés aux métabolites microbiens à partir de bases de données publiques de nature chimique, génomique et structurale. Cet article décrit un flux de travail en huit étapes qui prend en entrée un ensemble défini par l'utilisateur de métabolites d'origine microbiote intestinal et produit une liste restreinte classée de gènes candidats associés aux métabolites, de voies biologiques enrichies et de complexes protéine-ligand hiérarchisés sur la base de critères structuraux, destinés à un suivi expérimental. Le pipeline intègre : (i) le profilage chimioinformatique des métabolites ; (ii) la prédiction multicritère de cibles candidates à l'aide d'outils d'interaction protéine-composé chimique, d'une méthode de prédiction cible basée sur le ligand et d'un programme de docking moléculaire ; (iii) l'analyse de l'expression différentielle des gènes à partir de données transcriptomiques publiques ; (iv) l'intersection entre les cibles prédites et les gènes différentiellement exprimés ; (v) la construction d'un réseau d'interactions protéine-protéine et l'enrichissement en voies biologiques ; (vi) un docking moléculaire réalisé avec un programme de docking moléculaire ; (vii) une simulation de dynamique moléculaire de 200 ns à l'aide d'un moteur de dynamique moléculaire utilisant un champ de forces protéique dédié aux simulations de dynamique moléculaire ; et (viii) l'estimation de l'énergie libre de liaison par la méthode MM-PBSA. À titre d'exemple appliqué, neuf métabolites d'origine ou modifiés par le microbiote intestinal, représentant des acides gras à chaîne courte, des acides biliaires, des métabolites dérivés du tryptophane et l'uro lithine A, ont été analysés à l'aide du jeu de données transcriptomiques public du mucosa rectal dans le syndrome de l'intestin irritable avec constipation (IBS-C), GSE36701. Le flux de travail a classé 17 gènes uniques prédits comme associés aux métabolites et différentiellement exprimés dans ce jeu de données. Les analyses de docking, de dynamique moléculaire et de MM-PBSA ont permis de hiérarchiser structurellement cinq complexes métabolite-protéine : l'acide lithocholique-VDR, l'acide lithocholique-NR1H4/FXR, l'acide ursodésoxycholique-NR1H4/FXR, la tryptamine-HTR2A (simulée dans une bicouche lipidique explicite de 1-Palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC)) et l'uro lithine A-CASP3. Le protocole est conçu pour être adaptable à d'autres ensembles de métabolites, à d'autres jeux de données transcriptomiques liés à des maladies et à d'autres classes de cibles ; toutes les sorties constituent des prédictions computationnelles génératrices d'hypothèses, nécessitant une réplication indépendante au niveau transcriptomique, une validation au niveau protéique et des tests fonctionnels de réponse au ligand avant que des conclusions causales ou thérapeutiques puissent être établies.

Introduction

Le syndrome du côlon irritable avec constipation (IBS-C) est un trouble fonctionnel gastro-intestinal fréquent, caractérisé par des douleurs abdominales récurrentes, des modifications des habitudes intestinales, des ballonnements et de la constipation, avec une prévalence mondiale estimée à environ 10–15 % de la population générale^1,2. Les traitements pharmacologiques actuels, notamment les sécrétagogues, les prokinétiques et les antispasmodiques, peuvent améliorer certains symptômes chez un sous-groupe de patients ; toutefois, la réponse au traitement reste hétérogène et une rémission durable est rarement obtenue, reflétant la pathobiologie complexe et multifactorielle de cette affection^1,3,4. Une compréhension mécanistique plus complète de la façon dont les signaux microbiens intestinaux sont transduits au niveau muqueux est donc nécessaire afin de formuler des hypothèses testables pour de nouvelles cibles thérapeutiques.

La microflore intestinale contribue à l'homéostasie gastro-intestinale inférieure par la production et la biotransformation de métabolites chimiquement diversifiés, notamment des acides gras à chaîne courte (AGCC), des acides biliaires secondaires, des composés dérivés du tryptophane et des métabolites dérivés des polyphénols tels que les urolithines^5,6,7,8. Ces molécules communiquent avec les cellules de l'hôte par l'intermédiaire d'un large éventail, encore incomplètement caractérisé, de cibles moléculaires qui s'étend bien au-delà des récepteurs membranaires canoniques sensibles aux métabolites, incluant des récepteurs nucléaires, des enzymes cytosoliques, des protéines modifiant les histones, des précurseurs d'hormones peptidiques et des protéines intracellulaires de signalisation⁹. Des altérations de la composition des communautés microbiennes intestinales et des profils de métabolites ont été observées chez des patients souffrant du syndrome de l'intestin irritable (SII), fournissant ainsi une justification biologique à l'étude de la perturbation transcriptionnelle, dans la muqueuse rectale du SII-C, des gènes de l'hôte associés à la réponse aux métabolites microbiens¹⁰.

Le panel de neuf métabolites a été défini a priori afin de fournir un ensemble compact, chimiquement diversifié et biologiquement interprétable de petites molécules dérivées ou modifiées par la microflore intestinale. La sélection s'est appuyée sur cinq critères : la représentation des principales classes de métabolites microbiens impliqués dans la signalisation hôte-microbiote ; une exposition plausible ou avérée à la muqueuse intestinale distale ; la disponibilité d'identifiants PubChem non ambigus et de structures canoniques ; la taille moléculaire et la faisabilité structurale pour la prédiction de cibles et le docking basés sur les ligands ; et une plausibilité antérieure concernant une signalisation épithéliale, neuro-immunitaire, entéroendocrine, via des récepteurs nucléaires ou liée à la motilité, dans le syndrome du côlon irritable avec constipation (IBS-C). Le panel sélectionné comprenait le butyrate et le propionate en tant qu'acides gras à chaîne courte (SCFA) ; l'acide chénodésoxycholique, l'acide lithocholique et l'acide ursodésoxycholique en tant qu'acides biliaires ; la tryptamine, l'acide indole-3-propionique et l'acide indole-3-lactique en tant que métabolites dérivés du tryptophane ; et l'uellithine A en tant que métabolite de polyphénol d'origine microbienne intestinale^5,6,7,8,9,10.

La plupart des études computationnelles et expérimentales antérieures ont examiné des paires isolées de métabolite–récepteur ou de métabolite–enzyme, une approche qui ne rend pas compte de la nature distribuée et convergente de la signalisation par les métabolites microbiens à travers les voies de l'hôte^9,11. L'intégration de plusieurs étapes analytiques confère un pouvoir de filtrage mutuellement renforçant que nulle étape unique ne peut fournir indépendamment. La prédiction computationnelle de cibles contre des bases de données annotées fournit un large ensemble de protéines candidates de l'hôte pour chaque métabolite. L'intersection avec des données transcriptomiques pertinentes pour la maladie filtre considérablement cet ensemble, en conservant uniquement les candidats dont les transcrits sont modifiés dans le contexte de la maladie. Les analyses d'enrichissement en voies biologiques et de réseaux d'interactions protéine–protéine permettent ensuite de cartographier la liste réduite de candidats vers des modules biologiques connus. Le dockage moléculaire fournit une évaluation computationnelle initiale de la complémentarité des poches de liaison pour chaque complexe candidat, et une simulation supplémentaire de dynamique moléculaire (DM) de 200 ns avec une décomposition de l'énergie libre de liaison par la méthode MM-PBSA ajoute une dimension thermodynamique résolue dans le temps au classement structural, qui n'est pas accessible à partir des scores de dockage seuls. Effectuer chaque étape indépendamment, sans intégration systématique ni filtrage séquentiel, produirait des listes de candidats trop vastes pour être traitables expérimentalement et ne permettrait pas de détecter l'architecture convergente des voies.

Dans le cadre de ce protocole complet, par le terme « gène associé à un métabolite » (MAG), nous entendons un gène humain dont le produit protéique a été désigné comme cible moléculaire putative d’un ou plusieurs métabolites dérivés de la microflore intestinale par au moins une base de données de prédiction computationnelle validée, et dont le transcrit est exprimé de manière différentielle dans le jeu de données transcriptomiques pertinent pour la maladie utilisé pour illustrer le flux de travail. Cette définition opérationnelle inclut délibérément les récepteurs membranaires ainsi que les récepteurs nucléaires, les enzymes cytosoliques, les protéines de signalisation, les précurseurs d’hormones peptidiques et d’autres protéines intracellulaires. L’attribution du statut MAG ne constitue pas une preuve expérimentale qu’un métabolite se lie, forme un complexe protéine-ligand, active un récepteur, modifie l’abondance protéique ou provoque une maladie, mais plutôt une désignation générée par calcul, à vocation hypothétique, qui nécessite une validation expérimentale.

Ce protocole décrit le flux de travail computationnel complet en huit étapes (Figure 1) avec suffisamment de détails opérationnels pour permettre une réplication indépendante, une adaptation à d'autres profils de métabolites ou jeux de données de maladies, et une extension à d'autres contextes d'interaction hôte-microbiote. Le flux de travail est explicitement défini comme un cadre de génération d'hypothèses et de priorisation structurale qui s'appuie exclusivement sur des ressources publiques d'omique et de structures, sans inférer de modifications des concentrations en métabolites, des états d'activation des récepteurs, des variations d'expression protéique, de l'activité de signalisation en aval ou de la signification clinique à partir des seules sorties computationnelles. Ici, nous illustrons le protocole par un exemple concret utilisant neuf métabolites dérivés ou modifiés par le microbiote intestinal et le jeu de données transcriptomiques public GSE36701 provenant de la muqueuse rectale de patients atteints de SII-C, dans le but d'identifier des MAG et de hiérarchiser les complexes métabolite-protéine pour un suivi expérimental ultérieur.

Protocole

L'analyse a utilisé uniquement des données transcriptomiques publiques et dé-identifiées provenant de GSE36701, ainsi que des bases de données chimiques, protéiques et structurales publiques. Les bases de données ont été consultées entre janvier et mai 2026. Toute date d'accès ultérieure a été documentée dans le Tableau des matériaux séparé.

1. Conception de l'étude, exigences en matière de matériel et de logiciels

1. Définissez le flux de travail avant de commencer l'analyse. Utilisez huit étapes : sélection des métabolites, prédiction des cibles, analyse d'expression différentielle, chevauchement des gènes différentiellement exprimés (DEG), enrichissement en interactions protéine-protéine (PPI)/voies métaboliques, docking moléculaire, simulation de dynamique moléculaire (MD) et estimation MM-PBSA.
2. Précisez que le docking, la MD et l'estimation MM-PBSA sont uniquement des outils de priorisation structurale. N'interprétez pas ces résultats comme une preuve expérimentale de liaison, d'activation du récepteur, de modification de l'abondance protéique, d'efficacité thérapeutique ou de causalité dans la maladie.
3. Vérifiez le matériel informatique avant d'exécuter les simulations de dynamique moléculaire. Utilisez un système d'exploitation Linux 64 bits, un processeur d'au moins 6 cœurs, une plateforme de calcul accélérée par GPU disposant d'au moins 8 Go de mémoire vidéo (VRAM), ou une plateforme équivalente avec au minimum 8 Go de VRAM, 32 Go de RAM, et au moins 200 Go d'espace de stockage libre par système de dynamique moléculaire.
4. Répertoriez les logiciels essentiels : un moteur de dynamique moléculaire, un outil de docking moléculaire des ligands métabolites sur des protéines cibles, un outil de conversion de formats de fichiers chimiques¹², un générateur de ligands tridimensionnels, un outil de préparation des ligands, une boîte à outils pour la préparation des entrées de docking, un outil de docking moléculaire des ligands métabolites sur des protéines cibles, un environnement de programmation généraliste, un environnement de calcul statistique doté d'un cadre logiciel bioinformatique, et un package d'analyse d'expression génique différentielle.
5. Répertoriez les outils d'analyse structurale : un outil en ligne de construction de systèmes membranaires, un service de paramétrage de ligands compatible CHARMM, un outil de calcul de l'énergie de liaison en mécanique moléculaire/solvant continu, une bibliothèque de conversion de topologie et de paramètres moléculaires, un programme de visualisation moléculaire tridimensionnel, et un outil de visualisation moléculaire et de diagramme d'interactions bidimensionnel (2021) (voir la Table of Materials pour les liens de téléchargement et les informations de version).
6. Indiquez dans la Table of Materials/fichier d'environnement les identifiants exacts des champs de force protéiques utilisés pour les simulations de dynamique moléculaire, CGenFF, CHARMM-GUI, le cadre logiciel R/bioinformatique, et les versions de l'outil de calcul de l'énergie de liaison en mécanique moléculaire/solvant continu. Indiquez comme « non récupérable » les identifiants absents ; ne les déduisez pas.

2. Sélection des métabolites et caractérisation chimioinformatique

1. Définir le panel de métabolites avant la prédiction des cibles. Inclure le butyrate (PubChem CID : 264), le propionate (CID : 1032), l'acide chénodésoxycholique (CID : 10133), l'acide lithocholique (CID : 9903), l'acide ursodésoxycholique (CID : 31401), la tryptamine (CID : 1150), l'acide indole-3-propionique (CID : 3744), l'acide indole-3-lactique (CID : 92904) et l'uro lithine A (CID : 5488186).
2. Récupérer le système simplifié d'entrée de ligne de description moléculaire canonique (SMILES) et le CID PubChem pour chaque métabolite. Vérifier les synonymes et les structures en double avant la prédiction des cibles. Stocker les identifiants finaux dans la feuille maîtresse des métabolites.
3. Soumettre les chaînes SMILES canoniques à l'outil web de prédiction des propriétés physico-chimiques et du profil ADME¹³ (voir Tableau des matériaux). Noter la masse moléculaire, la surface polaire topologique (TPSA), le logP consensuel, le nombre de donneurs d'hydrogène, d'accepteurs d'hydrogène, de liaisons rotatives, l'absorption gastro-intestinale prédite, la prédiction de la glycoprotéine P, ainsi que les alertes Lipinski, Veber, Ghose, Egan, Muegge et PAINS.
4. Conserver les métabolites présentant une reconnaissance de structure réussie, une masse moléculaire ≤500 Da et aucune alerte PAINS. Noter tout critère non satisfait ainsi que la décision de conserver ou d'exclure le métabolite.
5. Attribuer les états d'ionisation avant la prédiction des cibles en aval et l'ajustement moléculaire (docking). Utiliser les formes carboxylates déprotonées pour le butyrate et le propionate, les formes acides carboxyliques neutres pour les acides biliaires, la forme ammonium protonée pour la tryptamine, et les formes neutres pour les autres métabolites.

3. Prédiction des cibles humaines candidates

1. Ouvrir une prédiction d'interaction chimique-protéine cible¹⁴ (voir Tableau des matériaux). Saisir le nom de chaque métabolite ou son CID PubChem, sélectionner Homo sapiens (identifiant taxonomique : 9606), et régler le score d'interaction combiné minimal à ≥0,700.
2. Hiérarchiser les canaux de preuves expérimentales et de bases de données curatées dans une prédiction d'interaction chimique-protéine cible. Télécharger le tableau complet des associations protéiques pour chaque métabolite.
3. Ouvrir un programme de docking moléculaire¹⁵ (voir Tableau des matériaux). Soumettre chaque chaîne SMILES canonique en sélectionnant Homo sapiens, et conserver les cibles ayant une probabilité ≥0,70.
4. Fusionner les résultats de la prédiction d'interaction chimique-protéine cible et du programme de docking moléculaire en un ensemble union pour chaque métabolite. Conserver toute cible répondant à l'un des seuils de base de données et supprimer les entrées de symboles génétiques en double exact.
5. Normaliser les entrées protéiques selon les symboles génétiques approuvés par le HUGO Gene Nomenclature Committee (HGNC) en utilisant le mappage des identifiants protéiques vers les symboles génétiques HGNC normalisés ou la base de données intégrée d'information génétique humaine (voir Tableau des matériaux). Résoudre les alias, les anciens symboles et les annotations d'isoformes en un seul symbole génétique par protéine.
6. Classer chaque cible comme récepteur membranaire, récepteur nucléaire, enzyme, protéine de signalisation intracellulaire, hormone peptidique, protéine liée aux hormones ou autre protéine intracellulaire. Enregistrer la classe dans le tableau des cibles.

4. Analyse du jeu de données transcriptomiques et expression différentielle des gènes

1. Accéder à GSE36701 via l'outil d'analyse différentielle d'expression génique en ligne du NCBI^16,17 (voir Tableau des Matériaux). Noter que la série de données comprend des données d'expression provenant de biopsies de la muqueuse rectale issues de groupes de patients atteints de syndrome de l'intestin irritable avec constipation (IBS-C), de syndrome de l'intestin irritable avec diarrhée prédominante (IBS-D), de syndrome de l'intestin irritable post-infectieux, ainsi que de volontaires sains¹⁸.
2. Rechercher dans GEO et ArrayExpress une cohorte indépendante pour la validation. Utiliser des combinaisons de termes tels que IBS-C, syndrome de l'intestin irritable avec prédominance de constipation, muqueuse rectale, muqueuse colique, biopsie, transcriptome, microarray et RNA-séquençage. Noter les référentiels, les termes de recherche, la date de recherche, et indiquer si une série de données de validation comparable a été identifiée.
3. Lancer l'outil d'analyse différentielle d'expression génique en ligne à partir de l'enregistrement GSE36701 (voir Tableau des Matériaux). Affecter les 18 échantillons IBS-C au groupe IBS-C, affecter les 40 volontaires sains au groupe témoin, et ne pas affecter les échantillons IBS-D et post-infectieux.
4. Effectuer l'analyse d'expression différentielle à l'aide du cadre du package d'analyse différentielle d'expression génique, avec correction du taux de fausses découvertes (FDR) de Benjamini-Hochberg¹⁹. Télécharger le tableau complet des résultats contenant l'ID de la sonde, le symbole du gène, le titre du gène, le logFC, l'AveExpr, la statistique t modérée, la valeur P brute et la valeur P ajustée.
5. Réduire les sondes à des entrées au niveau du gène. Supprimer les sondes dépourvues de symbole génique ; conserver la sonde présentant le FDR le plus faible en cas de symboles en double ; et utiliser la valeur absolue du logFC la plus élevée comme critère de départage.

5. Analyse du chevauchement cible-dég et évaluation statistique

1. Intersectionner chaque liste de cibles prédites spécifique à un métabolite avec la liste de gènes différentiellement exprimés (DEG) au seuil de FDR < 0,05. Enregistrer les gènes en commun, le métabolite d'origine, le logFC, la valeur P ajustée et le sens d'expression.
2. Fusionner les listes d'intersections spécifiques à chaque métabolite en une liste non redondante de gènes associés aux métabolites (MAG). Compter le nombre total de cibles prédites, les intersections spécifiques à chaque métabolite et le nombre total de MAG uniques.
3. Évaluer la cohérence directionnelle au niveau des sondes pour les gènes possédant plusieurs sondes. Indiquer tout gène pour lequel les sondes présentent un désaccord concernant le sens d'expression.
4. Construire le tableau de contingence du test exact de Fisher en utilisant le nombre total d'entrées regroupées par gène, le nombre total de DEG, le nombre total de cibles prédites uniques et le nombre de MAG observés. Calculer la valeur P unilatérale, le rapport de cotes (odds ratio) et l'intervalle de confiance à 95 % à l'aide de l'implémentation du test exact de Fisher.
5. Si le taux de DEG dans la population de référence dépasse 50 %, rapporter l'intersection comme descriptive plutôt que comme un enrichissement indépendamment validé. Considérer une régulation négative uniforme comme un motif directionnel descriptif, sauf si un test de directionnalité distinct est réalisé.

6. Analyse du réseau d'interactions protéine-protéine et enrichissement des voies métaboliques

1. Soumettez la liste complète et unique de MAG à la construction du réseau d'interactions protéine-protéine et à l'enrichissement des voies métaboliques²⁰ (voir Tableau des matériaux). Sélectionnez Homo sapiens et réglez le score d'interaction minimal à 0,700.
2. Exportez le réseau combiné de construction d'interactions protéine-protéine et d'enrichissement des voies métaboliques ainsi que le tableau complet des interactions. Si l'extraction textuelle produit une topologie artificiellement dense, désactivez l'extraction textuelle et conservez les canaux expérimentaux, de co-expression et issus de bases de données.
3. Générez des sous-réseaux par classe de métabolite pour les MAG associés aux AGCC, aux acides biliaires, et au tryptamine/sérotoninergiques. Utilisez les mêmes paramètres d'organisme et de confiance pour la construction du réseau d'interactions protéine-protéine et l'enrichissement des voies métaboliques.
4. Effectuez l'analyse de construction du réseau d'interactions protéine-protéine et d'enrichissement des voies métaboliques contre la Kyoto Encyclopedia of Genes and Genomes (KEGG)²¹, Reactome²² et Gene Ontology (GO) Biological Process^23,24. Appliquez un seuil de FDR de Benjamini–Hochberg (BH) <0,05 et exportez tous les tableaux d'enrichissement.

7. Docking moléculaire

1. Récupérer les structures expérimentales de récepteurs à partir de la base de données Protein Data Bank (RCSB PDB) du Research Collaboratory for Structural Bioinformatics²⁵ (voir la Table des matériels). Utiliser VDR/1DB1, NR1H4/FXR/3DCT, CASP3/2DKO et HTR2A/6A93 pour les cinq complexes protéine-ligand prioritaires.
2. Préparer chaque récepteur en conservant la chaîne A et en supprimant les molécules d'eau, les ligands co-cristallisés, les cofacteurs, les ions et les enregistrements HETATM non protéiques. Pour 6A93, supprimer le segment de fusion avec la lysozyme T4 avant la préparation du récepteur.
3. Ajouter les hydrogènes polaires, attribuer les charges de Gasteiger et sauvegarder chaque récepteur au format PDBQT à l’aide de la boîte à outils de structure moléculaire. Examiner les états de protonation des histidines du site de liaison avant la conversion en PDBQT et documenter les états choisis.
4. Générer la structure tridimensionnelle de chaque ligand à l’aide de la boîte à outils de conversion de fichiers de structure chimique. Minimiser l’énergie avec le champ de force universel (Universal Force Field, UFF) pendant 500 étapes, attribuer l’état d’ionisation à pH 7,4, attribuer les charges de Gasteiger, puis sauvegarder au format PDBQT.
5. Définir une boîte d’auto-assemblage (docking) de 25 Å x 25 Å x 25 Å centrée sur le centroïde du ligand co-cristallisé. Utiliser les coordonnées centrales (10, 19, 33) pour VDR, (137, 31, 78) pour FXR, (37, 34, 32) pour CASP3 et (12, −1, 61) pour HTR2A.
6. Effectuer un auto-assemblage moléculaire des ligands métabolites sur les protéines cibles^26,27 avec un paramètre d’exhaustivité (exhaustiveness) = 8, une graine (seed) = 42, un nombre de modes (num_modes) = 9 et une plage énergétique (energy_range) = 3 kcal/mol. Enregistrer le score Vina le plus élevé et les valeurs d’écart quadratique moyen (RMSD) pour toutes les poses.
7. Sélectionner le mode 1 pour chaque complexe prioritaire. Générer des diagrammes bidimensionnels ligand-résidus à l’aide d’un outil de visualisation moléculaire et d’un outil de diagramme d’interaction bidimensionnel, ainsi que des vues tridimensionnelles récepteur-ligand dans un programme de visualisation moléculaire tridimensionnel.
8. Effectuer des contrôles de redocking pour VDR/1DB1 et FXR/3DCT. Valider la configuration d’auto-assemblage du récepteur lorsque le RMSD des atomes lourds est <2,0 Å par rapport à la pose cristallographique.
9. Effectuer des contrôles de cross-docking en assemblant l’acide lithocholique (LCA) dans CASP3 et la tryptamine dans VDR. Comparer les scores cognates et non cognates, et noter les cas où la différence de score est <1,0 kcal/mol.

8. Simulation de dynamique moléculaire

1. Générez les paramètres du ligand à l'aide d'un service de paramétrisation de ligands compatible CHARMM²⁸ (voir Tableau des matériaux). Examinez tous les scores de pénalité et signalez tout paramètre présentant une pénalité >50.
2. Convertissez les fichiers de flux du ligand en fichiers .itp et .prm compatibles avec le moteur de dynamique moléculaire à l'aide d'un script de conversion de topologie de champ de force. Combinez les fichiers de topologie du ligand et de la protéine pour chaque complexe.
3. Appliquez une redistribution des masses des hydrogènes à l'aide d'une bibliothèque de conversion de topologie et de paramètres moléculaires. Générez des topologies aqueuses en utilisant le champ de force protéique destiné aux simulations de dynamique moléculaire²⁹ et un modèle explicite d'eau à trois sites³⁰.
4. Solvatez les complexes aqueux dans une boîte dodécaédrique avec un espacement minimal entre le soluté et le bord de 1,2 nm au moins. Neutralisez les systèmes et ajoutez du NaCl jusqu'à une concentration de 0,15 M.
5. Construisez le système membranaire tryptamine-HTR2A à l'aide d'un outil en ligne de construction de systèmes membranaires^31,32,33 (voir Tableau des matériaux). Utilisez des coordonnées du récepteur alignées sur une base de données d'orientations membrane-protéine³⁴ (voir Tableau des matériaux), une bicouche pure de POPC, des couches d'eau de 22,5 Å et une concentration de NaCl de 0,15 M.
6. Minimisez l'énergie de tous les systèmes par la méthode de la plus forte pente, jusqu'à 50 000 pas. Vérifiez la convergence lorsque Fmax <1000 kJ·mol^-1·nm^-1 avant l'équilibration.
7. Équilibrez les systèmes aqueux selon les ensembles à nombre de particules, volume et température constants (NVT) et à nombre de particules, pression et température constants (NPT). Équilibrez le système membranaire à l'aide d'un protocole en ligne en six étapes de préparation et d'équilibration multistade du système moléculaire, avec un relâchement progressif des contraintes.
8. Effectuez une simulation de dynamique moléculaire productive de 200 ns pour les cinq complexes. Utilisez un pas de temps de 4 fs avec la redistribution des masses des hydrogènes (HMR), un thermostat V-rescale à 310 K, un barostat de Parrinello-Rahman à 1 bar, l'électrostatique Particle Mesh Ewald (PME)³⁵ et des contraintes LINCS³⁶.
9. Analysez les trajectoires finales à l'aide d'utilitaires d'analyse de trajectoires de dynamique moléculaire. Calculez le RMSD du squelette, la fluctuation quadratique moyenne des atomes Cα (RMSF), le rayon de giration, la surface accessible au solvant (SASA) et les liaisons hydrogène entre la protéine et le ligand, en utilisant les 150 ns finales comme fenêtre d'analyse principale.

9. Estimation de l'énergie libre de liaison par MM-PBSA

1. Extraire des instantanés de trajectoire pour l'analyse MM-PBSA. Utiliser 2 001 cadres pour chaque complexe en solution aqueuse et 201 cadres traités pour le sous-système HTR2A intégré dans la membrane.
2. Exécuter l'outil de calcul de l'énergie de liaison en mécanique moléculaire/solvant continu³⁷ avec solvatation de Poisson-Boltzmann, une constante diélectrique interne = 1, une constante diélectrique externe = 80, une solvatation non polaire basée sur la SASA, et sans correction d'entropie. Indiquer l'énergie libre moyenne de liaison et l'écart type.
3. Effectuer une décomposition par résidu pour les cinq complexes. Indiquer les résidus stabilisants et déstabilisants dont les contributions absolues sont ≥0,5 kcal·mol⁻¹.

    

Résultats

Cibles associées aux métabolites candidates

Les neuf métabolites ont produit des ensembles cibles prédits hétérogènes selon un programme de prédiction d'interactions chimique-protéine et un programme de docking moléculaire. Le propionate, la tryptamine, les acides biliaires et l'uro lithine A ont révélé plusieurs cibles ayant une pertinence connue dans la signalisation gastro-intestinale. Le paysage cible prédit comprenait des récepteurs membranaires canoniques, des récepteurs nucléaires, des enzymes intracellulaires, des protéines de signalisation et des protéines associées aux hormones peptidiques. Les résultats ultérieurs sont donc décrits comme des gènes associés aux métabolites (MAG) plutôt que comme des découvertes limitées aux récepteurs (Tableau 1).

Évaluation comparative par rapport aux interactions métabolite-protéine rapportées

Afin de comparer les résultats de prédiction de cibles à la connaissance expérimentale existante, les relations prédites entre métabolites et cibles ont été classées en trois niveaux de preuves : (i) interactions directes ou au niveau de la classe proche entre métabolite et protéine étayées expérimentalement, dans lesquelles le métabolite ou un métabolite endogène étroitement apparenté a été signalé comme se liant à la protéine codée, l'activant, l'inhibant ou la régulant fonctionnellement ; (ii) interactions étayées par voie métabolique ou par classe de cible, dans lesquelles la cible prédite appartient à une voie métabolique sensible aux métabolites ou à une famille de récepteurs établie, mais pour lesquelles les preuves directes concernant la paire exacte métabolite-protéine sont limitées ; et (iii) associations exclusivement computationnelles pour lesquelles aucune interaction expérimentale directe n'a été identifiée dans la littérature examinée. Ce processus de comparaison a été utilisé pour contextualiser, et non valider, les MAG prédits.

Plusieurs prédictions ont confirmé des données biologiques précédemment rapportées. L'interaction propionate-FFAR2 a été considérée comme soutenue expérimentalement, car FFAR2/GPR43 est un récepteur canonique des acides gras à chaîne courte. L'interaction butyrate-HDAC3 a été classée comme soutenue expérimentalement ou par voie indirecte, car le butyrate est un inhibiteur reconnu des histones désacétylases, et la prédiction recouvrait un membre de la famille des HDAC. Les prédictions associées aux acides biliaires impliquant NR1H4/FXR et VDR ont été jugées étayées par la biologie établie des récepteurs nucléaires des acides biliaires, en particulier pour les acides biliaires hydrophobes comme l'acide lithocholique (LCA) ; les prédictions concernant FXR associées à l'acide ursodésoxycholique (UDCA) ont été interprétées avec prudence, car l'UDCA est généralement un ligand de FXR plus faible ou dépendant du contexte. Les prédictions associées à la tryptamine et concernant HTR1B, HTR2A, HTR2B et HTR6 ont été classées comme soutenues par la voie sérotoninergique, mais non confirmées comme des interactions directes spécifiques aux récepteurs, car la tryptamine est une monoamine dérivée du tryptophane microbien, et les récepteurs de la sérotonine sont des régulateurs établis de la motilité et de la sécrétion gastro-intestinales. L'interaction urolithine A-CASP3 a été considérée comme soutenue par voie indirecte en raison de liens publiés entre l'urolithine A et les réponses apoptotiques ou liées aux caspases, mais non par une preuve directe de liaison à CASP3. Les interactions acide indole-3-lactique-KYAT1 et acide indole-3-propionique-KYAT1 ont été conservées comme hypothèses uniquement computationnelles, car la littérature générale soutient la signalisation de l'hôte par les dérivés microbiens de l'indole, mais pas une liaison directe de ces métabolites précis à KYAT1^7,8,38,39,40.

Par conséquent, le Tableau 1 distingue la nomination computationnelle de cibles du niveau de soutien expérimental ou de voie préalable. Il indique également, pour chaque cible, la source de prédiction (prédiction d'interaction chimique-protéine, programme de docking moléculaire, ou les deux), le score d'interaction combiné pour la prédiction d'interaction chimique-protéine, et la probabilité du programme de docking moléculaire lorsque la cible a été identifiée par un tel programme. Les cibles prédites sans preuve expérimentale directe préalable sont décrites comme des gènes candidats associés à des métabolites nécessitant une validation indépendante au niveau protéique et par réponse au ligand.

Recouvrement entre les gènes cibles prédits et les gènes différentiellement exprimés dans le SII-C

L'intersection des listes de cibles prédites combinées par union et des résultats d'expression différentielle au niveau des gènes a identifié 17 gènes uniques prédits comme associés à des métabolites, qui étaient exprimés de manière significativement différentielle dans la comparaison IBS-C par rapport aux volontaires sains. Les 17 gènes étaient tous régulés à la baisse. Ce groupe comprenait des récepteurs membranaires et nucléaires (CASR, FFAR2, GPR68, HTR1B, HTR2A, HTR2B, HTR6, NR1H4, TBXA2R, VDR) ainsi que des protéines non réceptrices (CASP3, GCG, GNAQ, GPHN, HDAC3, KYAT1, MLN) (Tableau 1, Figure 2A,B).

Tous les 17 génomes microbiens assemblés (MAG) présentaient un taux de faux positifs (FDR) inférieur à 0,05 ; 16 des 17 ont atteint le seuil de FDR plus strict < 0,001, le gène restant (HTR1B) étant significatif au seuil FDR < 0,05. Sept des 17 cibles (CASP3, GCG, GNAQ, GPHN, GPR68, HDAC3, TBXA2R) ont satisfait à la fois au FDR < 0,001 et un logarithme absolu₂ changement d'expression supérieur à 1,0 (intervalle de logFC de −1,34 à −1,10), indiquant une régulation négative forte et constante pour ce sous-ensemble. Les autres cibles ont montré une régulation négative modérée mais statistiquement significative (|logFC| compris entre 0,45 et 0,97). Ce profil descriptif uniforme a été interprété avec prudence, compte tenu des caractéristiques d'expression à l'échelle du génome de la base de données (voir évaluation statistique ci-dessous).

Évaluation statistique du recouvrement entre les gènes différentiellement exprimés et les cibles

Afin d'évaluer formellement la signification statistique du recoupement de 17 gènes, un test exact de Fisher unilatéral a été appliqué en utilisant les 17 gènes cibles prédits comme ensemble interrogé et l'ensemble des 18 296 entrées uniques regroupées par gène détectées dans GSE36701 comme fonds génomique. Parmi ce fonds, 17 296 gènes (94,5 %) étaient exprimés de manière différentielle à un FDR < 0,05, reflétant une suppression transcriptionnelle quasi universelle dans la comparaison de la muqueuse rectale de patients atteints de SCI-C. Les 17 gènes cibles prédits appartenaient tous aux gènes exprimés de manière différentielle (recoupement observé : 17/17, soit 100 %). Étant donné un taux de 94,5 % de gènes différentiellement exprimés dans le fonds, le recoupement attendu pour tout ensemble aléatoire de 17 gènes est de 16,1 gène. Le test exact de Fisher a donné une valeur p = 0,384 avec un odds ratio corrigé pour la continuité de 2,03 (intervalle de confiance à 95 % : 0,12–33,73), ce qui n'était pas statistiquement significatif au seuil α = 0,05 (Figure 3A–C).

Ce résultat indique que le chevauchement observé de 17 sur 17 ne dépasse pas le chevauchement attendu par hasard selon le profil d'expression génomique de ce jeu de données. En conséquence, ces résultats sont interprétés comme un modèle descriptif directionnel, dans lequel les 17 cibles prédites étaient toutes de manière constante et significative sous-régulées dans le tissu muqueux rectal de patients atteints de SII-C, plutôt que comme une preuve d'enrichissement statistique ou de validation indépendante par rapport à un fond génomique. Un test formel d'enrichissement nécessiterait une réplication dans des jeux de données transcriptomiques présentant des profils de différentiel d'expression plus sélectifs, dans lesquels nettement moins de la moitié des gènes atteignent la significativité. Il convient de souligner que la sous-régulation uniforme des 17 gènes chevauchants constitue une observation descriptive plutôt qu'un résultat statistique validé indépendamment, car le fond de gènes différentiellement exprimés dans ce jeu de données est lui-même principalement sous-régulé ; une direction descendante commune parmi les gènes chevauchants était donc attendue et n'a pas fait l'objet d'un test formel de directionnalité. Cette direction uniforme ne doit donc pas être interprétée comme une preuve statistique indépendante de régulation coordonnée et spécifique au métabolite.

Modèles spécifiques aux métabolites

Le propionate présentait le plus grand nombre de gènes chevauchants, notamment CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN et TBXA2R, suggérant une implication possible des voies de signalisation sensibles aux acides gras à chaîne courte et associées à Gq. Le butyrate partageait des gènes communs avec HDAC3, ce qui est cohérent avec la biologie des histones désacétylases associée au butyrate, bien que la simple diminution de l’ARNm n’établisse pas formellement une modification de la sensibilité au butyrate. Les chevauchements liés aux acides biliaires incluaient les récepteurs nucléaires VDR et NR1H4, tous deux reconnus comme des effecteurs du signal des acides biliaires dans l’intestin^38,39. La tryptamine partageait des gènes communs avec HTR1B, HTR2A, HTR2B et HTR6, impliquant une implication potentielle du système sérotoninergique, un système jouant un rôle bien établi dans la motilité et la sécrétion gastro-intestinales⁴⁰. L’acide indole-3-lactique et l’acide indole-3-propionique présentaient un chevauchement avec KYAT1, et l’urolithine A avec CASP3.

Enrichissement des voies

L'analyse d'enrichissement fonctionnel des 17 gènes communs a révélé des voies liées à la signalisation en aval des récepteurs couplés aux protéines G (RCPG), à la signalisation Gαq, à la liaison des ligands aux RCPG, aux synapses sérotoninergiques, à l'interaction entre ligands neuroactifs et récepteurs, à la transduction du signal calcique, à la signalisation de l'AMPc et à la sécrétion d'hormones peptidiques. Ces résultats sont cohérents avec la composition du jeu de gènes et soutiennent sa cohérence biologique, mais reflètent l'annotation fonctionnelle des gènes soumis plutôt qu'une preuve indépendante d'activité au niveau des voies.

Structure du réseau d'interactions protéine-protéine

La construction du réseau d'interactions protéine-protéine et l'analyse d'enrichissement des voies ont été interprétées à partir de trois réseaux complémentaires. Dans le méta-réseau combiné de 17 gènes (Réseau 1), la structure la plus évidente appuyée par annotation était un composant de signalisation GPCR/Gαq centré sur GNAQ, reliant GNAQ à des gènes associés aux récepteurs, notamment TBXA2R, CASR, HTR2A et HTR2B. Une connectivité limitée des récepteurs de la sérotonine était également conservée, principalement entre HTR2A et HTR2B, tandis que plusieurs autres gènes restaient isolés ou faiblement connectés au seuil de confiance sélectionné. Le réseau spécifique au propionate (Réseau 2) présentait une topologie plus restreinte, GNAQ conservant des liaisons appuyées par annotation avec CASR et TBXA2R, tandis que FFAR2, GPR68, GCG, GPHN et MLN étaient isolés ou faiblement connectés. Le réseau tryptamine/sérotonine (Réseau 3) incluait HTR1B, HTR2A, HTR2B et HTR6 ; au sein de ce sous-ensemble, HTR2A et HTR2B présentaient la connexion principale appuyée par annotation, tandis que HTR1B et HTR6 n'étaient pas directement connectés au seuil choisi (Figure 4A–C).

Le docking moléculaire

Un dockage moléculaire a été réalisé sur cinq complexes métabolite-protéine sélectionnés. Les paires acide biliaire-récepteur nucléaire ont montré des scores Vina plus favorables que l'urolythine A-CASP3 et la tryptamine-HTR2A. Le complexe LCA-VDR a obtenu le meilleur score avec −10,0 kcal/mol, suivi par LCA-NR1H4/FXR (−9,9 kcal/mol) et UDCA-NR1H4/FXR (−9,4 kcal/mol). L'urolythine A-CASP3 et la tryptamine-HTR2A ont obtenu des scores inférieurs mais néanmoins raisonnables de −7,1 kcal/mol (Tableau 2).

Pour le complexe LCA-VDR (ID PDB : 1DB1), la conformation prédite était stabilisée par une liaison hydrogène conventionnelle entre l'oxygène carboxylate de la LCA et Ser278 (4,29 Å), ainsi que par d'importants contacts hydrophobes impliquant Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233 et His397, et des contacts additionnels de type van der Waals avec Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275 et Phe150. La conformation classée en tête présentait un score Vina de −10,0 kcal/mol, une taille de cavité de 2055 ų et un centre de grille situé à (10, 19, 33) (Tableau 3, Figure 5A,B).

Pour le complexe LCA-NR1H4/FXR (ID PDB : 3DCT), le score de docking de −9,9 kcal/mol était associé à des liaisons hydrogène prédites impliquant His294 et Ile335, une interaction π-Sigma avec His294, ainsi qu'à des contacts hydrophobes Alkyle ou π-Alkyle impliquant Met290, Met328, Ala291, Leu287, Ile352 et His447, des contacts supplémentaires de type van der Waals contribuant à l'accommodation du squelette stéroïdien dans la poche du FXR (Tableau 4, Figure 6A,B).

La pose prédite du complexe UDCA-NR1H4/FXR (PDB ID : 3DCT) présentait une liaison hydrogène classique avec His447 (3,66 Å), une autre liaison hydrogène avec Gly322 (3,46 Å), une interaction π-anion avec Val325 (4,96 Å) et une liaison carbone-hydrogène avec Trp469 (4,51 Å). La carte des interactions a également révélé des contacts donneur-donneur défavorables avec Arg395 (3,89 Å) et Gln396 (3,40 Å), suggérant que le score Vina inférieur de l'UDCA par rapport à celui de la LCA dans le même site de liaison pourrait être dû à une géométrie locale ou à des interactions électrostatiques moins favorables (Tableau 5, Figure 7A,B).

Dans le complexe urolithine A-CASP3 (ID PDB : 2DKO), le mode de liaison prédit présentait des liaisons hydrogène conventionnelles avec Gln161 (3,78 et 4,19 Å), Ser120 (3,95 Å) et Arg207 (3,05 et 3,77 Å), et était davantage stabilisé par des interactions π-cation avec Arg207, une liaison hydrogène π-donneur avec Cys163, ainsi que des contacts π-alkyle et de van der Waals supplémentaires impliquant Arg64, Ala162, His121, Ser205 et Trp206 (Tableau 6, Figure 8A,B).

Pour le complexe tryptamine-HTR2A (ID PDB : 6A93), la conformation prédite était stabilisée par un pont salin électrostatique entre l'amine protonée de la tryptamine et Asp155, l'aspartate conservé de l'hélice transmembranaire 3 (D3.32 selon la numérotation de Ballesteros-Weinstein) qui ancre l'amine protonée des ligands aminergiques dans les récepteurs de la sérotonine et apparentés^41,42,43, ainsi que par des liaisons hydrogène avec Thr160 et Ser159, des contacts aromatiques avec Phe340 et Trp336, et des interactions π-alcoyle avec Val156 et Ile163. Des contacts supplémentaires de type van der Waals avec Tyr370, Phe339, Ser242, Phe243, Phe332 et Leu123 soutenaient un mode de liaison au niveau du site orthostérique (Tableau 7, Figure 9A,B).

Validation du protocole d'arrimage

Afin d'évaluer la fiabilité du protocole d'assemblage moléculaire, deux expériences de contrôle complémentaires ont été réalisées. Pour les contrôles de réassemblage (positifs), les ligands co-cristallisés ont été extraits de leurs structures cristallographiques de référence et réassemblés dans leurs sites de fixation natifs. La pose prédite la mieux classée pour l'analogue de la vitamine D, VDX, dans VDR/1DB1 présentait un écart de 0,87 Å par rapport à la position cristallographique, et le ligand co-cristallisé WAY-362450 dans FXR/3DCT présentait un écart de 1,79 Å ; ces deux valeurs étant inférieures au seuil d'acceptation conventionnel de 2,0 Å, elles confirment la validité géométrique du protocole d'assemblage pour ces systèmes récepteurs (Figure 10A,B). Pour les contrôles d'assemblage croisé (négatifs), l'acide lithocholique a été assemblé dans la caspase-3 (2DKO), une protéase à cystéine dont il n'est pas un ligand connu, produisant un score prédit (−8,3 kcal/mol) 1,7 kcal/mol plus faible que celui obtenu avec sa cible naturelle VDR (−10,0 kcal/mol), ce qui est compatible avec une sélectivité prévue du site de fixation. L'assemblage de la tryptamine dans VDR a donné un score prédit de −6,4 kcal/mol comparé à −7,1 kcal/mol sur sa cible naturelle HTR2A, une différence de 0,7 kcal/mol qui se situe dans la marge d'incertitude rapportée pour les scores d'assemblage moléculaire des ligands métaboliques sur les protéines cibles, indiquant ainsi une sélectivité prédite modeste pour ce ligand plus petit (Figure 10C). Dans leur ensemble, ces contrôles indiquent que le protocole d'assemblage reproduit les géométries de liaison connues et permet de distinguer les paires homologues des paires non homologues dans les conditions testées, tout en restant des prédictions computationnelles qui ne remplacent pas les mesures expérimentales d'affinité (Table 8).

Simulation de dynamique moléculaire

Des simulations de dynamique moléculaire ont été réalisées pour les cinq complexes prioritaires sur des trajectoires de production de 200 ns. Les quatre complexes solubles et les complexes récepteurs nucléaires ont été simulés dans un solvant aqueux explicite, tandis que le complexe tryptamine-HTR2A a été simulé dans un doublet lipidique explicite de POPC afin de fournir un environnement membranaire physiologiquement approprié pour ce récepteur couplé aux protéines G. Les analyses ont évalué la stabilité dynamique des poses d'ancrage dans des conditions dépendantes du temps et ont permis de comparer le comportement structural relatif entre les complexes (Tableau 9).

Le profil RMSD du complexe LCA-VDR/1DB1 a montré une courte période d'équilibration durant les 10 premières ns, suivie d'un plateau stable, avec des fluctuations principalement comprises entre 0,20 et 0,28 nm (Figure 11A). Les valeurs de RMSF étaient faibles, et les fluctuations de la charpente étaient < 0,15 nm pour la plupart des résidus (Figure 11B). L'analyse des liaisons hydrogène a révélé un réseau persistant de 2 à 5 liaisons hydrogène, avec des augmentations occasionnelles jusqu'à 7 (Figure 11C). Le rayon de giration (Rg) a été maintenu dans la plage de 1,25 à 1,75 nm, et la surface accessible au solvant (SASA) a été maintenue aux environs de 130 nm² (Figure 11D,).

Le complexe urolithine A-CASP3/2DKO a présenté une activité dynamique plus importante. Le RMSD a initialement augmenté, puis oscillé entre 0,4 et 0,7 nm, avec un bref épisode de forte déviation vers 165 ns (Figure 12A). L'analyse RMSF a révélé une grande mobilité au niveau des résidus, avec les plus grandes fluctuations dans la région de la boucle flexible autour du résidu 175 (Figure 12B). L'analyse des liaisons hydrogène a mis en évidence un réseau initial étendu d'environ 2 à 5 liaisons durant les 30 à 40 premières ns, suivi principalement par 0 à 2 liaisons intermittentes (Figure 12C). Les profils correspondants du rayon de giration et de la surface accessible au solvant (SASA) sont présentés dans la Figure 12D,E.

Pour les systèmes biliaires NR1H4/FXR (3DCT), le profil du RMSD de la structure principale est resté dans une plage relativement étroite pendant la majeure partie de la trajectoire (Figure 13A), tandis que le profil du RMSF a montré une mobilité réduite dans les régions centrales et des fluctuations plus importantes dans les régions flexibles (Figure 13B). Le complexe LCA-3DCT a maintenu environ trois à quatre liaisons hydrogène persistantes tout au long de la trajectoire, alors que le complexe UDCA-3DCT a présenté des fluctuations plus marquées des liaisons hydrogène et une réduction du nombre de ces liaisons après environ 125 ns. Les profils du rayon de giration pour les systèmes liés au LCA et à l'UDCA sont présentés respectivement dans les Figure 13C,D, et les profils correspondants de SASA sont illustrés dans les Figure 13E,F.

Dynamique moléculaire de la membrane du complexe tryptamine-HTR2A

Le complexe tryptamine-HTR2A/6A93 a été simulé pendant 200 ns dans un bicouche lipidique explicite de POPC composée de 258 molécules lipidiques, un modèle d'eau à trois sites explicite et 0,15 M de NaCl, pour une taille totale du système d'environ 100 925 atomes^33,44,45. Le récepteur est resté stablement intégré dans la bicouche tout au long de la trajectoire (Figure 14). Le RMSD de la chaîne principale est passé d'environ 0,10 nm à un plateau stable aux alentours de 0,15–0,20 nm durant les 100 premières ns et est resté stable par la suite, avec toutes les valeurs inférieures à 0,25 nm, indiquant que le récepteur a conservé une conformation stable dans l'environnement membranaire sans dépliement global (Figure 15A). L'analyse RMSF par résidu a montré de faibles fluctuations dans le cœur hélicoïdal transmembranaire, avec une mobilité plus élevée attendue dans les boucles et les régions terminales, en accord avec la flexibilité typique des RCPG (Figure 15B). Le rayon de giration a été étroitement confiné entre environ 2,06 et 2,12 nm, et la SASA a fluctué dans une bande étroite sans dérive progressive, confirmant ainsi le maintien du faisceau transmembranaire compact (Figure 15C,D).

Les liaisons hydrogène entre la protéine et le ligand ont été maintenues tout au long de la trajectoire (Figure 15E), avec des fluctuations importantes du nombre de liaisons hydrogène, allant de 1 à 3. Afin d'évaluer spécifiquement la persistance de l'interaction ionique clé, la distance minimale entre l'azote ammonium protoné de la tryptamine et les atomes d'oxygène carboxylate de l'Asp155 (D3.32) a été surveillée tout au long de l'ensemble de la trajectoire. Cette distance est restée étroitement distribuée autour d'une moyenne de 0,270 nm (minimum 0,247 nm, maximum 0,424 nm), et le contact du pont salin (< 0,4 nm) a été maintenu pendant 99,9 % de la simulation, avec seulement deux brèves excursions transitoires, sans aucun événement de dissociation prolongé (Figure 16). Ces résultats suggèrent que l'interaction ionique conservée impliquant l'Asp155 était suffisante pour stabiliser la tryptamine au sein du site orthostérique de l'HTR2A tout au long de la simulation en milieu membranaire.

Énergie libre de liaison MM-PBSA et décomposition par résidu

Une analyse MM-PBSA a été réalisée afin d'ajouter une couche supplémentaire de priorisation énergétique aux cinq complexes (Tableau 10). Pour les quatre complexes aqueux, la décomposition par résidu a identifié les principaux contributeurs énergétiques de chaque mode de liaison prédit. Dans le complexe LCA-VDR/1DB1, le ligand et Gln317 ont apporté une contribution favorable, tandis que Trp286 a présenté une contribution défavorable. Dans le complexe urolithine A-CASP3/2DKO, Arg64 et Arg207 ont montré des contributions par résidu fortement négatives, indiquant une stabilisation polaire ou électrostatique importante ; néanmoins, la trajectoire correspondante est restée très dynamique, ce qui démontre que des énergétiques favorables au niveau des résidus ne garantissent pas à elles seules une stabilité durable du complexe. Pour les systèmes 3DCT, la liaison de LCA était principalement due à Arg331, tandis que la liaison de UDCA impliquait un réseau énergétique plus étendu comprenant Glu326, Asp394, Arg395, Arg441 et Asp470. Dans l'ensemble des quatre systèmes aqueux, la décomposition MM-PBSA a confirmé la priorisation relative des complexes basés sur la LCA.

Pour le complexe tryptamine-HTR2A/6A93 intégré dans la membrane, une analyse MM-PBSA a été réalisée sur le sous-système protéine-ligand extrait de la trajectoire dans la bicouche^46,47. Des contributions favorables ont été observées pour le ligand et Asp155 (D3.32), qui s'est avéré de loin le résidu le plus important en termes de stabilisation, en accord avec l'interaction par pont salin identifiée dans les analyses de docking et de distance au cours de la trajectoire. Trp137 a présenté la contribution défavorable la plus importante au niveau individuel parmi les résidus environnant la poche orthostérique (Ser86, Phe87, Phe133, Phe140, Phe141, Val156, Ser159, Thr160, Ile163, Val167, Tyr171), qui forment ensemble le réseau de contacts aromatiques et polaires tapissant la poche de liaison. Ces valeurs représentent des estimations computationnelles relatives utilisées pour la hiérarchisation structurale et ne correspondent pas à des affinités de liaison expérimentales.

Figure 1 : Flux de travail computationnel pour la priorisation de gènes de l'hôte associés aux métabolites dans le SII-C. Représentation schématique du flux de travail en huit étapes intégrant la sélection des métabolites, la prédiction de cibles, l'expression différentielle transcriptomique, l'analyse de chevauchement, l'enrichissement en réseau et en voies, le docking moléculaire, la simulation de dynamique moléculaire et l'analyse de l'énergie libre de liaison MM-PBSA. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 2 : Analyse différentielle de l'expression et recouvrement entre métabolites cibles dans la muqueuse d'IBS-C. (A) Graphique volcan représentant l'expression différentielle au niveau des gènes dans GSE36701. Points bleus : gènes significativement sous-régulés ; points rouges : gènes significativement sur-régulés ; points gris : gènes non significatifs. Les gènes associés aux métabolites présentant un recouvrement sont indiqués. (B) Diagramme de Venn montrant le recouvrement entre 330 cibles uniques de métabolites prédites et les gènes sous-régulés dans GSE36701 ; 17 gènes sont communs. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 3 : Évaluation statistique des 17 gènes cibles de métabolites prédits par rapport à GSE36701. (A) Variation logarithmique (log2) de l'expression pour chacun des 17 gènes, colorés selon le niveau de significativité. (B) Taux d'expression différentielle des gènes de fond par rapport aux cibles prédites, analysé par le test exact de Fisher. (C) Un tableau de contingence 2×2 utilisé pour le test exact de Fisher. Les 17 cibles étaient toutes significativement régulées à la baisse ; l'intersection est interprétée comme un motif descriptif directionnel plutôt qu'un enrichissement statistique. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 4 : Construction du réseau d'interaction protéine-protéine composite et réseaux d'enrichissement des voies métaboliques pour les gènes associés aux métabolites communs. (A) Réseau 1 : méta-réseau combiné des 17 gènes. (B) Réseau 2 : réseau spécifique au propionate comprenant huit gènes (CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN, TBXA2R). (C) Réseau 3 : réseau de la tryptamine/sérotonine comprenant quatre gènes (HTR1B, HTR2A, HTR2B, HTR6). Les réseaux ont été générés pour Homo sapiens au minimum, à l'aide de la construction de réseaux d'interaction protéine-protéine et d'un seuil de confiance d'enrichissement des voies ≥ 0,700. Les arêtes représentent une association fonctionnelle étayée par des annotations Veuillez cliquer ici pour consulter une version agrandie de cette figure.

Figure 5 : Représentation structurale tridimensionnelle et bidimensionnelle de l'acide lithocholique en complexe avec le VDR (ID PDB : 1DB1). (A) Représentation tridimensionnelle en surface et en forme de bande avec l'acide lithocholique représenté sous forme de sphères. (B) Carte d'interactions bidimensionnelle montrant la liaison hydrogène avec la sérine 278 et les contacts hydrophobes environnants ainsi que les interactions de van der Waals. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 6. Représentation tridimensionnelle et bidimensionnelle de la structure de l'acide lithocholique en complexe avec NR1H4/FXR (ID PDB : 3DCT). (A) Représentation tridimensionnelle en surface et en forme schématique. (B) Carte d'interaction bidimensionnelle montrant les liaisons hydrogène avec His294 et Ile335, une interaction π-Sigma, ainsi que les contacts environnants. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 7 : Représentation structurale tridimensionnelle et bidimensionnelle de l'acide ursodésoxycholique en complexe avec NR1H4/FXR (ID PDB : 3DCT). (A) Représentation 3D en surface et sous forme de schéma fil de fer. (B) Carte d'interactions bidimensionnelle montrant des liaisons hydrogène avec His447 et Gly322, une interaction π-anion avec Val325, une liaison carbone-hydrogène avec Trp469, ainsi que des contacts donneur-donneur défavorables avec Arg395 et Gln396. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 8 : Représentation structurale tridimensionnelle et bidimensionnelle de l'acide urolithine A en complexe avec CASP3 (ID PDB : 2DKO). (A) Représentation tridimensionnelle en surface et en schéma fil de fer. (B) Carte d'interactions bidimensionnelle montrant les liaisons hydrogène avec Gln161, Ser120 et Arg207, les interactions π-cation avec Arg207, une liaison hydrogène de type donneur π avec Cys163, ainsi que les contacts environnants. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 9 : Représentation structurale tridimensionnelle et bidimensionnelle de la tryptamine en complexe avec HTR2A (ID PDB : 6A93). (A) Représentation 3D en surface et sous forme de schéma générée à l’aide d’un logiciel de visualisation moléculaire tridimensionnel. (B) Carte d’interaction bidimensionnelle générée à l’aide d’un outil de visualisation moléculaire et d’un diagramme d’interaction 2D, illustrant le pont salin Asp155 et d’autres interactions au niveau du site de liaison. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 10 : Validation du protocole de docking. (A,B) Redocking des ligands co-cristallisés dans VDR/1DB1 (RMSD 0,87 Å) et FXR/3DCT (RMSD 1,79 Å) ; les poses cristallographiques et redockées sont superposées, toutes deux en dessous du seuil d'acceptation de 2,0 Å. (C) Sélectivité en cross-docking : scores Vina homologues par rapport aux scores hétérologues pour l'acide lithocholique et la tryptamine. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 11. Analyse de la trajectoire de dynamique moléculaire du complexe LCA-VDR/1DB1 sur 200 ns. (A) Profil de RMSD. (B) Profil de RMSF. (C) Nombre de liaisons hydrogène. (D) Profil du rayon de giration. (E) Profil de SASA. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 12 : Analyse de la trajectoire de dynamique moléculaire du complexe urolithine A-CASP3/2DKO sur 200 ns. (A) Profil du RMSD montrant des fluctuations conformationnelles importantes et un événement transitoire de forte déviation vers 165 ns. (B) Profil du RMSF montrant une flexibilité importante au niveau des résidus vers le résidu 175. (C) Nombre de liaisons hydrogène. (D) Profil du rayon de giration. (E) Profil de la SASA. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 13 : Analyse de la trajectoire de dynamique moléculaire des systèmes d'acides biliaires NR1H4/FXR (3DCT) sur 200 ns. (A) Profil RMSD du squelette du complexe 3DCT. (B) Profil RMSF du squelette. (C) Profil du rayon de giration pour 3DCT-LCA. (D) Profil du rayon de giration pour 3DCT-UDCA. (E) Profil SASA pour 3DCT-LCA. (F) Profil SASA pour 3DCT-UDCA. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 14 : Le complexe tryptamine-HTR2A intégré dans une bicouche lipidique POPC explicite. Le récepteur est représenté sous forme de schéma traversant la bicouche, les lipides POPC sous forme de lignes avec les groupements phosphate des têtes polaires mis en évidence, et la tryptamine située dans le site orthostérique. L'eau est représentée au-dessus et au-dessous de la membrane. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 15 : Analyse de la trajectoire de dynamique moléculaire du complexe tryptamine-HTR2A/6A93 sur 200 ns dans une bicouche lipidique POPC explicite. (A) Profil RMSD du squelette. (B) Profil RMSF par résidu. (C) Profil du rayon de giration. (D) Profil de la surface accessible au solvant (SASA). (E) Nombre de liaisons hydrogène entre la protéine et le ligand. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Figure 16 : Persistance de l'interaction ionique entre la tryptamine et Asp155 (D3.32) au cours de la trajectoire de 200 ns dans la membrane. La distance minimale entre l'azote ammonium de la tryptamine et les atomes d'oxygène carboxylate d'Asp155 est tracée en fonction du temps ; la ligne en pointillés indique le seuil de contact de pont salin à 0,4 nm. Le contact a été maintenu pendant 99,9 % de la simulation. Veuillez cliquer ici pour visualiser une version agrandie de cette figure.

Symbole du gène	Métabolite(s) d'origine	Classification fonctionnelle	log2FC	FDR (valeur-P ajustée)	Niveau de signification
GCG	Propionate	Protéine liée à l'hormone peptidique	−1,342	1,97e−7	FDR <0,001 & |logFC| > 1
HDAC3	Butyrate	Enzyme	−1,234	2,44e−6	FDR <0,001 & |logFC| > 1
CASP3	Urolithine A	Enzyme	−1,198	6,66e−7	FDR <0,001 & |logFC| > 1
GPR68	Propionate	Récepteur membranaire	−1,137	4,35e−6	FDR <0,001 & |logFC| > 1
GNAQ	Propionate	Protéine de signalisation intracellulaire	−1,122	1,05e−6	FDR <0,001 & |logFC| > 1
GPHN	Propionate	Autre protéine intracellulaire	−1,109	1,13e−6	FDR <0,001 & |logFC| > 1
TBXA2R	Propionate	Récepteur membranaire	−1,104	4,04e−7	FDR <0,001 & |logFC| > 1
HTR6	Tryptamine	Récepteur membranaire	−0,967	2,17e−5	FDR <0,001
VDR	Acide lithocholique	Récepteur nucléaire	−0,942	5,73e−7	FDR <0,001
HTR2A	Tryptamine	Récepteur membranaire	−0,937	4,99e−6	FDR <0,001
FFAR2	Propionate	Récepteur membranaire	−0,889	1,44e−4	FDR <0,001
NR1H4	Acide lithocholique / Acide ursodésoxycholique	Récepteur nucléaire	−0,861	3,68e−6	FDR <0,001
HTR2B	Tryptamine	Récepteur membranaire	−0,702	1,29e−4	FDR <0,001
MLN	Propionate	Protéine liée à l'hormone peptidique	−0,605	7,39e−5	FDR <0,001
KYAT1	Acide indole-3-lactique / Acide indole-3-propionique	Enzyme	−0,530	3,61e−4	FDR <0,001
CASR	Propionate	Récepteur membranaire	−0,483	4,05e−4	FDR <0,001
HTR1B	Tryptamine	Récepteur membranaire	−0,455	3,18e−2	FDR <0,05

Tableau 1 : Gènes cibles associés aux métabolites prédits chevauchant les gènes différentiellement exprimés dans le jeu de données de la muqueuse rectale de l'IBS-C. Tous les gènes chevauchants répertoriés étaient régulés à la baisse. Le tableau 1 est fourni séparément sous forme de tableau de feuille de calcul et indique, pour chaque cible, le ou les métabolites d'origine, la catégorie fonctionnelle, la source de prédiction de la cible (une prédiction d'interaction chimique-protéine, un programme d'ancrage moléculaire, ou les deux), un score d'interaction combiné provenant d'une prédiction d'interaction chimique-protéine, et, lorsque disponible, la probabilité du programme d'ancrage moléculaire, le niveau de prédiction, le log2 de changement d'expression, le FDR et le niveau de significativité de l'expression. Source : les valeurs d'expression génique ont été obtenues à partir du tableau de différentielle d'expression GSE36701 avec les gènes regroupés (sonde par gène ayant le FDR le plus bas). Les sources de prédiction des cibles et les valeurs de confiance ont été compilées à partir des sorties d'une prédiction d'interaction chimique-protéine et d'un programme d'ancrage moléculaire, en utilisant comme seuils un score d'interaction combiné ≥ 0,700 pour la prédiction d'interaction chimique-protéine et une probabilité ≥ 0,70 pour le programme d'ancrage moléculaire. Les scores de la prédiction d'interaction chimique-protéine sont des scores combinés sur une échelle de 0 à 1 ; STP désigne la probabilité du programme d'ancrage moléculaire. Niveau 1 = soutien strict de la prédiction d'interaction chimique-protéine ; niveau 1+ = soutien strict de la prédiction d'interaction chimique-protéine corroboré par un programme d'ancrage moléculaire.

Complexe	Protéine (ID PDB)	Ligand	Score Vina (kcal/mol)	Taille de la cavité (A^3)	Centre de la grille X,Y,Z (A)	Boîte de recherche (A)
LCA-VDR	VDR (1DB1)	Acide lithocholique	−10.0	2055	10, 19, 33	25 x 25 x 25
LCA-NR1H4/FXR	NR1H4/FXR (3DCT)	Acide lithocholique	−9.9	3395	137, 31, 78	25 x 25 x 25
UDCA-NR1H4/FXR	NR1H4/FXR (3DCT)	Acide ursodésoxycholique	−9.4	3395	137, 31, 78	25 x 25 x 25
Urolithine A-CASP3	CASP3 (2DKO)	Urolithine A	−7.1	233	37, 34, 32	25 x 25 x 25
Tryptamine-HTR2A	HTR2A (6A93)	Tryptamine	−7.1	3238	12, −1, 61	25 x 25 x 25

Tableau 2 : Résultats du dockage moléculaire : meilleurs scores de dockage moléculaire des ligands métabolites aux protéines cibles et paramètres des cavités pour les cinq complexes protéine-ligand prioritaires. La taille de la cavité est indiquée en ų. Source : Docking_Validation/Results/Docking_Validation_Results.xlsx, feuille « Original_Docking_Scores ». Un dockage moléculaire des ligands métabolites aux protéines cibles ; exhaustivité = 8, graine = 42 (fixe), num_modes = 9 pour tous les complexes ; la conformation classée au premier rang (mode 1) est indiquée.

Type d'interaction	Résidu(s)	Distance (Å)	Remarques
Liaison hydrogène conventionnelle	Ser278	4,29	Oxygène carboxylate de LCA
Contact hydrophobe / Pi-Alkyle	Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233, His397	-
Contact de Van der Waals	Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275, Phe150	-

Tableau 3 : Modes de liaison générés pour le docking de l'acide lithocholique avec le VDR (ID PDB : 1DB1). Source : outil de visualisation moléculaire et de diagramme d'interactions ligand-résidu en deux dimensions, diagrammes d'interaction ligand-résidu 2D, tel que rapporté dans les Résultats du manuscrit (Docking moléculaire). « - » indique qu'aucune valeur de distance n'a été rapportée individuellement pour ce contact.

Type d'interaction	Résidu(s)	Distance (Å)	Remarques
Liaison hydrogène	His294	-
Liaison hydrogène	Ile335	-
Interaction Pi-Sigma	His294	-
Alkyle / Pi-Alkyle (hydrophobe)	Met290, Met328, Ala291, Leu287, Ile352, His447	-
Contact de Van der Waals	Résidus supplémentaires du site (non spécifiés individuellement dans la source)	-	Permet l'accommodation du squelette stéroïdique

Tableau 4 : Modes de liaison générés pour le dockage de l'acide lithocholique avec NR1H4/FXR (ID PDB : 3DCT).

Source : visualisation moléculaire et outil de diagramme d'interaction bidimensionnel, diagrammes d'interaction ligand-résidu 2D, tels que rapportés dans les résultats du manuscrit (Dockage moléculaire). « - » indique qu'aucune valeur de distance n'a été rapportée individuellement pour ce contact.

Type d'interaction	Résidu(s)	Distance (Å)	Notes
Liaison hydrogène classique	His447	3,66
Liaison hydrogène	Gly322	3,46
Interaction pi-anion	Val325	4,96
Liaison carbone-hydrogène	Trp469	4,51
Contact donneur-donneur défavorable	Arg395	3,89
Contact donneur-donneur défavorable	Gln396	3,40

Tableau 5 : Modes de liaison générés pour le docking de l'acide ursodésoxycholique avec NR1H4/FXR (ID PDB : 3DCT). Source : outil de visualisation moléculaire et de diagramme d'interaction ligand-résidu en deux dimensions (diagrammes d'interaction ligand-résidu 2D), tel que décrit dans les Résultats du manuscrit (Docking moléculaire). Le signe « - » indique qu'aucune valeur de distance n'a été rapportée individuellement pour ce contact.

Type d'interaction	Résidu(s)	Distance (Å)	Remarques
Liaison hydrogène conventionnelle	Gln161	3,78
Liaison hydrogène conventionnelle	Gln161	4,19	deuxième contact
Liaison hydrogène conventionnelle	Ser120	3,95
Liaison hydrogène conventionnelle	Arg207	3,05
Liaison hydrogène conventionnelle	Arg207	3,77	deuxième contact
Interaction Pi-cation	Arg207	-
Liaison hydrogène donneur Pi	Cys163	-
Contact Pi-alkyle / van der Waals	Arg64, Ala162, His121, Ser205, Trp206	-

Tableau 6 : Modes de liaison générés pour le docking de l'urothine A avec CASP3 (ID PDB : 2DKO). Source : outil de visualisation moléculaire et diagrammes d'interaction ligand-résidu en deux dimensions, tels qu'indiqués dans les Résultats du manuscrit (Docking moléculaire). Le symbole « - » indique qu'aucune valeur de distance n'a été rapportée individuellement pour ce contact.

Type d'interaction	Résidu(s)	Distance (A)	Remarques
Pont salin électrostatique	Asp155 (D3.32)	-	amine protonée de la tryptamine
Liaison hydrogène	Thr160	-
Liaison hydrogène	Ser159	-
Contact aromatique	Phe340, Trp336	-
Interaction Pi-Alkyle	Val156, Ile163	-
Contact de Van der Waals	Tyr370, Phe339, Ser242, Phe243, Phe332, Leu123	-

Tableau 7 : Modes de liaison générés pour l'ancrage de la tryptamine avec HTR2A (PDB ID : 6A93). Source : outil de visualisation moléculaire et de diagramme d'interaction ligand-résidu en deux dimensions, diagrammes d'interaction ligand-résidu 2D, tel que rapporté dans les Résultats du manuscrit (Ancrage moléculaire). « - » indique qu'une valeur de distance n'a pas été rapportée individuellement pour ce contact.

(A) Validation par redockage (contrôles positifs)
ID PDB	Protéine	Ligand co-cristallisé	Score Vina (kcal/mol)	RMSD (Å)	Seuil (Å)	Résultat
1DB1	VDR	VDX (analogue de la vitamine D)	−13.0	0.87	2.0	PASS
3DCT	FXR	WAY-362450 (064)	−11.9	1.79	2.0	PASS
(B) Validation par cross-dockage (contrôles négatifs)
Ligand	Cible homologue (PDB)	Score homologue (kcal/mol)	Cible non homologue (PDB)	Score non homologue (kcal/mol)	Delta (kcal/mol)	Sélectivité
Acide lithocholique	VDR (1DB1)	−10.0	CASP3 (2DKO)	−8.3	1.7	Confirmée
Tryptamine	HTR2A (6A93)	−7.1	VDR (1DB1)	−6.4	0.7	Modeste (dans l'incertitude de Vina +/−0,5–1,0)

Tableau 8 : Résultats de la validation du protocole d'ancrage : valeurs RMSD du réancrage (témoins positifs) et scores d'ancrage croisé (témoins négatifs). Source : Docking_Validation/Results/Docking_Validation_Results.xlsx et Docking_Validation/Logs/*.log (un ancrage moléculaire de ligands métabolites sur des protéines cibles, exhaustivité = 8, graine = 42, boîte de 25 Å × 25 Å × 25 Å). RMSD calculé par correspondance des atomes lourds et des noms d'atomes (aucune superposition).

Complexe	RMSD (nm), moyenne + / – SD (intervalle)	Rg (nm), moyenne + / – SD (intervalle)	SASA (nm^2), moyenne + / – SD (intervalle)	liaisons hydrogène, moyenne + / – SD (intervalle)	RMSF (nm), moyenne (max)
LCA-VDR/1DB1	0.230 + / – 0.025 (0.167–0.296)	1.889 + / − 0.009 (1.863–1.919)	130.4 + / − 2.3 (122.3–137.4)	1.9 + / − 0.9 (0–7)	0.093 (max 0.600 au résidu 120)
LCA-NR1H4/FXR/3DCT	0.190 + / – 0.020 (0.135–0.281)	1.824 + / − 0.008 (1.804–1.849)	129.7 + / − 2.3 (121.9–138.1)	3.8 + / − 0.7 (1–6)	0.113 (max 0.298)
UDCA-NR1H4/FXR/3DCT	0.190 + / – 0.020 (0.135–0.281)	1.834 + / − 0.013 (1.809–1.921)	131.0 + / −3.4 (121.6–143.5)	1.1 + / − 1.1 (0–5)	0.113 (max 0.298)
Urolithin A-CASP3/2DKO	0.521 + / – 0.058 (0.244–0.755)	1.892 + / − 0.024 (1.839–1.984)	134.9 + / − 3.0 (126.4–146.4)	0.6 + / − 0.7 (0–3)	1.172 (max 2.532 au résidu 175)
Tryptamine-HTR2A/6A93 (membrane)	0.177 + / –0.017 (0.131–0.227)	2.089 + / − 0.007 (2.070–2.116)	165.1 + / − 2.7 (156.–172.7)	1.7 + / − 0.7 (0–4)	0.090 (max 0.319)

Tableau 9 : Résumé du comportement en simulation de dynamique moléculaire de 200 ns pour les cinq complexes protéine-ligand prioritaires, incluant le système tryptamine-HTR2A intégré dans la membrane. Source : fichiers d'utilitaires d'analyse de trajectoire de dynamique moléculaire (.xvg) — gmx rms, gmx gyrate, gmx sasa, gmx hbond, gmx rmsf — calculés sur les 150 ns finales (50–200 ns) de chaque phase de production de 200 ns, conformément à l'étape 8.8 du protocole. RMSD/Rg ajustés sur le squelette peptidique ; rayon de la sonde SASA de 0,14 nm ; seuil donneur-accepteur de liaison hydrogène à 0,35 nm / 30 °. LCA-3DCT et UDCA-3DCT partagent une trajectoire commune du squelette peptidique (RMSD, RMSF), avec des valeurs spécifiques au ligand pour Rg/SASA/les liaisons hydrogène.

Tryptamine-HTR2A/6A93 (membrane) — décomposition quantitative par résidu
Résidu	Contribution totale à ddG (kcal/mol), moyenne + / − DS	Sens
Asp155 (D3.32)	−89.94 + / − 6.81	Stabilisant (dominant)
Tryptamine (ligand)	−13.01 + / − 6.22	Stabilisant
Tyr171	13.62 + / − 4.54	Déstabilisant
Val167	23.32 + / − 3.96	Déstabilisant
Val156	20.03 + / − 3.81	Déstabilisant
Thr160	4.86 + / − 3.64	Déstabilisant
Ser159	24.16 + / − 3.48	Déstabilisant
Ser86	24.48 + / − 3.65	Déstabilisant
Phe87	35.18 + / − 4.04	Déstabilisant
Phe133	32.80 + / −3.70	Déstabilisant
Phe140	30.63 + / − 3.84	Déstabilisant
Phe141	35.25 + / − 3.55	Déstabilisant
Ile163	27.64 + / − 3.71	Déstabilisant
Trp137	53.77 + / − 4.32	Déstabilisant (le plus défavorable)
Autres complexes — résidus identifiés dans la décomposition par résidu (qualitatif)
Complexe	Résidu	Sens
LCA-VDR/1DB1	Ligand (LCA)	Favorable
LCA-VDR/1DB1	Gln317	Favorable
LCA-VDR/1DB1	Trp286	Défavorable
LCA-NR1H4/FXR/3DCT	Arg331	Favorable (dominant)
UDCA-NR1H4/FXR/3DCT	Glu326	Réseau mixte/distribué
UDCA-NR1H4/FXR/3DCT	Asp394	Réseau mixte/distribué
UDCA-NR1H4/FXR/3DCT	Arg395	Réseau mixte/distribué
UDCA-NR1H4/FXR/3DCT	Arg441	Réseau mixte/distribué
UDCA-NR1H4/FXR/3DCT	Asp470	Réseau mixte/distribué
Urolithin A-CASP3/2DKO	Arg64	Fortement favorable (polaire/électrostatique)
Urolithin A-CASP3/2DKO	Arg207	Fortement favorable (polaire/électrostatique)

Tableau 10 : Décomposition MM-PBSA par résidu RÉSUMÉ COURT : résidus stabilisants et déstabilisants (contribution absolue ≥ 0,5 kcal mol⁻¹) pour chacun des cinq complexes protéine-ligand prioritaires, incluant le système tryptamine-HTR2A intégré à la membrane. Source : Membrane Simulation/03_MMPBSA/results/FINAL_DECOMP_MMPBSA.dat (outil de calcul de l'énergie de liaison par mécanique moléculaire/solvant continu Generalized Born (GB), décomposition par résidu, « Complex : Total Energy Decomposition »). Les numéros de résidus ont été convertis à partir de la numérotation interne du système construit par CHARMM-GUI (décalage +68) vers la numérotation PDB originale 6A93 utilisée ailleurs dans cet article.

Source : Données de simulation de dynamique moléculaire antérieures/1DB1,2KD0, LCA & UDCA_3DCT}/mmpbsa_*/Decomposition_NORMAL_GB_Complex_TDC*.svg et résultats de l'article (énergie libre de liaison MM-PBSA et décomposition par résidu). Ces quatre complexes ne possèdent pas de sortie numérique par résidu au format .dat/.csv dans le répertoire du projet (seules des illustrations SVG rendues, contenant du texte vectoriel non extractible par machine, sont disponibles) ; seules l'identité des résidus et la direction favorable ou défavorable sont indiquées dans le texte de l'article. Les contributions exactes en kcal/mol pour ces quatre complexes ne sont pas disponibles dans le dépôt source.

Discussion

Cette étude computationnelle exploratoire démontre un flux de travail intégré et reproductible permettant de hiérarchiser les gènes de l'hôte et les complexes protéine-ligand associés aux métabolites microbiens, appliqué ici à un jeu de données transcriptomiques public de la muqueuse rectale de patients atteints de SII-C. À l'aide de ce flux de travail, un sous-ensemble de gènes prédits comme associés aux métabolites microbiens présentait un chevauchement avec des gènes régulièrement sous-régulés dans le jeu de données, s'agrégeant au sein de voies associées aux RCPG, à la signalisation sérotoninergique, à la signalisation calcique, aux récepteurs-ligands neuroactifs et aux récepteurs nucléaires, des systèmes de plus en plus impliqués dans la communication entre le microbiote et l'hôte^48,49. Ces résultats doivent être interprétés strictement comme générant des hypothèses : l'analyse ne mesure ni les concentrations des métabolites microbiens, ni l'abondance des protéines réceptrices, ni la liaison du ligand, ni l'activation du récepteur, ni la signalisation en aval, ni la motilité, ni la sécrétion, ni les réponses à la douleur, ni les résultats cliniques. La conclusion la plus solide qui puisse être tirée est que les gènes et voies identifiés constituent des candidats pour une validation expérimentale, plutôt que des mécanismes établis de la maladie.

L'importance cruciale de ce protocole, par rapport aux travaux antérieurs examinant des paires métabolite-récepteur isolées, réside dans son intégration de la prédiction ciblée, des données transcriptomiques publiques, de l'analyse de réseaux, du docking avec des contrôles de validation, de la dynamique moléculaire et de la méthode MM-PBSA au sein d'un pipeline séquentiel de priorisation. Chaque étape affine et contextualise l'ensemble de candidats produit par l'étape précédente, et ce filtrage séquentiel rend la liste finale de candidats expérimentalement abordable. Le module GPCR centré sur GNAQ et le module associé aux récepteurs de la sérotonine identifiés ici sont biologiquement plausibles, compte tenu du rôle de la voie de signalisation Gq dans l'activation de la phospholipase C, la production d'inositol 1,4,5-trisphosphate, la mobilisation du calcium, la sécrétion et la fonction entéroendocrine, ainsi que des rôles établis des signaux dérivés des acides gras à chaîne courte et du tryptophane dans l'homéostasie muqueuse, et du signal sérotoninergique dans la motilité gastro-intestinale, la sécrétion, la sensibilité viscérale et la communication intestin-cerveau^10,18,50,51.

Une caractéristique méthodologique essentielle de cette étude est le traitement du récepteur HTR2A intégré dans la membrane. Étant donné qu'une simulation en phase soluble ne peut pas reproduire l'environnement lipidique qui gouverne le comportement conformationnel d'un récepteur couplé aux protéines G, le complexe tryptamine-HTR2A a été simulé dans un bicouche POPC explicite. Dans cet environnement membranaire, le récepteur est resté structurellement stable tout au long de la trajectoire complète de 200 ns, et le pont salin entre l'ammonium de la tryptamine et Asp155 (D3.32) a été maintenu pratiquement durant toute la simulation. Le fait que trois lignes indépendantes de preuves — la pose d'ancrage, la distance de contact persistante au cours de la trajectoire et la contribution dominante par résidu selon l'analyse MM-PBSA — convergent vers la même interaction conservée au niveau de D3.32 confère une cohérence interne au mode de liaison de la tryptamine prédit, qui reproduit la géométrie de liaison canonique des ligands aminergiques aux récepteurs de la sérotonine.

Il existe certaines questions méthodologiques à prendre en compte lors de la reproduction de ce flux de travail. Les erreurs dans la structure canonique ou les composés interférant dans les tests pan-essais (composés signalés comme PAINS) se propagent par la prédiction ciblée et le docking, ce qui exige une sélection précise des métabolites et une curatation chimioinformatique rigoureuse. Réduisez au minimum les ensembles cibles induits par le bruit en appliquant systématiquement les critères de confiance (prédiction de cibles par interaction chimique-protéine : ≥0,700 ; programme de docking moléculaire : ≥0,70 ; construction du réseau d'interactions protéine-protéine et enrichissement des voies : ≥0,700). Les cibles prédites doivent être regroupées par catégorie fonctionnelle afin d'éviter de caractériser à tort tous les gènes associés aux métabolites comme des récepteurs. Un prétraitement précis des structures PDB, la minimisation de l'énergie du ligand et le positionnement de la grille autour des résidus de liaison connus sont des aspects essentiels du docking, tandis que les contrôles de redocking et de cross-docking introduits ici fournissent une mesure objective de la justesse de la méthodologie de docking. La marge de reproductibilité en dynamique moléculaire est définie par une combinaison de paramétrisation du champ de forces, de solvatation adéquate ou de construction de membrane, d'équilibrage étagé et d'échantillonnage suffisant en phase de production.

Les adaptations et étapes de dépannage typiques incluent l'assouplissement des seuils si la prédiction cible ne retourne aucun résultat, la vérification de la cohérence directionnelle au niveau des sondes pour les gènes comportant plusieurs sondes, ainsi que l'interprétation de nœuds isolés dans la construction de réseaux d'interactions protéine-protéine et l'enrichissement de voies comme dépendant du seuil plutôt que biologiquement non pertinents. Pour les récepteurs membranaires, une simulation explicite du doublet lipidique doit être utilisée plutôt qu'une simulation aqueuse, comme illustré par la méthode HTR2A décrite ici. Lorsqu'une décomposition énergétique par résidu est nécessaire, le calcul doit être effectué à l'aide d'un moteur capable de décomposition, et la numérotation des résidus indiquée doit être mise en cohérence avec la numérotation native du récepteur afin d'éviter toute ambiguïté. Nous recommandons de considérer les résultats d'enrichissement de voies principalement comme un contexte organisationnel pour la liste de candidats, plutôt que comme une validation au niveau des voies. Mécaniquement, un enrichissement des termes GPCR, sérotoninergiques ou de signalisation calcique se produira chaque fois que la liste de gènes contient plusieurs gènes de récepteurs de la sérotonine, indépendamment de la co-régulation au niveau protéique. Les valeurs RMSD, Rg et RMSF du système HTR2A intégré à la membrane doivent être interprétées en tenant compte du doublet lipidique : une diminution de Rg dans la trajectoire ultérieure peut refléter une adaptation conformationnelle induite par le doublet du faisceau transmembranaire plutôt qu'un déroulement global, et les liaisons hydrogène persistantes entre le ligand et la protéine doivent être interprétées conjointement avec la stabilité globale du RMSD.

Les limites de cette étude sont importantes et restreignent l'interprétation. La recherche s'appuyait sur un seul jeu de données public relativement restreint, et une recherche dans les principaux dépôts transcriptomiques publics (outil web d'analyse différentielle d'expression génique et ArrayExpress) n'a pas permis d'identifier un jeu de données indépendant de transcriptome de la muqueuse rectale chez les patients atteints de SCI-C, de conception et de plateforme comparables, pouvant servir de cohorte de réplication au moment de l'analyse. L'absence de réplication transcriptomique indépendante constitue une limitation majeure, et aucune affirmation formulée dans cet article ne doit être interprétée comme une validation externe des résultats obtenus à partir d'un seul jeu de données. Le jeu de données révèle une expression différentielle quasi universelle (environ 94,5 % des gènes sont significatifs, dont la grande majorité sont sous-régulés), une caractéristique qui rend les statistiques classiques d'enrichissement non informatives et empêche de tirer des conclusions sur la spécificité de la sous-régulation des gènes cibles par rapport au fond génomique ; l'intersection est donc rapportée comme un motif directionnel descriptif plutôt que comme un enrichissement statistique. Le transcriptome muqueux en bloc ne permet pas de distinguer une régulation génique réelle des modifications de la composition cellulaire. L'expression des ARNm ne détermine pas l'abondance des protéines ni la réponse fonctionnelle. Les bases de données de prédiction de cibles souffrent de biais d'annotation, et les résultats des simulations d'ancrage, de dynamique moléculaire (MD) et de MM-PBSA dépendent du choix du champ de force, de la paramétrisation du ligand, de la position initiale, de la durée de la simulation et de la qualité de l'échantillonnage. Les identifiants exacts de version ou de construction pour certains composants de serveurs web et de paquets, notamment le service de paramétrisation de ligands compatible CHARMM, CHARMM-GUI et l'environnement de calcul statistique

Les versions des packages et des outils de calcul d'énergie de liaison en mécanique moléculaire/solvant continu n'ont pas pu être entièrement récupérées à partir de l'archive du projet et doivent être indiquées comme disponibles dans le Tableau des matériaux séparé. Les valeurs MM-PBSA sont des estimations relatives, ne contiennent pas de terme explicite d'entropie configurationnelle et ne doivent pas être interprétées comme des affinités expérimentales. L'étude ne comporte pas de données métabolomiques et ne permet pas de déterminer si la disponibilité des ligands est modifiée dans le SII-C, ni si les modifications d'expression observées sont des causes, des conséquences, des réponses compensatoires ou des corrélations sans lien direct.

Les applications futures de cette méthode devraient inclure une réplication indépendante au niveau du transcriptome, une validation par réaction de polymérisation en chaîne quantitative (qPCR) et au niveau protéique, une localisation cellulaire par transcriptomique unicellulaire ou spatiale, un profilage métabolomique des classes de métabolites pertinentes, ainsi que des tests fonctionnels de réponse aux ligands dans des colonoïdes dérivés de patients, des explants de muqueuse ou des modèles comparables. Des comparaisons avec des cohortes de syndrome de l'intestin irritable (SII) à prédominance de diarrhée, de SII mixte, de maladie inflammatoire de l'intestin et de constipation non liée au SII^1,2 permettraient d'établir la spécificité de la maladie. Pour la composante structurale, la réplication des trajectoires de dynamique moléculaire (MD), la réalisation d'analyses de sensibilité avec des positions initiales alternatives, ainsi qu'une documentation complète du dépôt des fichiers de topologie, des trajectoires et des entrées et sorties de MM-PBSA renforceraient davantage la reproductibilité. Des tests expérimentaux de réponse aux ligands restent nécessaires pour déterminer si les complexes prioritaires sont fonctionnellement pertinents ; les résultats actuels ne soutiennent aucune affirmation clinique ou thérapeutique.

Matériaux

Liste des matériaux utilisés dans cet article

Nom	Entreprise	Numéro de catalogue	Commentaires
AutoDock Vina	Scripps Research / open source	v1.2.7; https://vina.scripps.edu/	Accolement moléculaire de ligands métabolites à des protéines cibles.
CGenFF/ParamChem	SilcsBio / University of Maryland	v4.6; https://cgenff.com/	Paramétrisation du champ de forces des ligands pour la dynamique moléculaire.
Champ de forces CHARMM36m	Développeurs CHARMM / open source	CHARMM36m; https://www.charmm.org/charmm/resources/charmm-force-fields/	Champ de forces protéique utilisé pour les simulations de dynamique moléculaire.
CHARMM-GUI Membrane Builder	CHARMM-GUI / Lehigh University	Serveur web ; version exacte non récupérable ; https://www.charmm-gui.org/?doc=input/membrane	Construction et configuration de l'équilibration du système membranaire POPC explicite.
Discovery Studio Visualizer	BIOVIA (Dassault Systèmes)	2021; https://discover.3ds.com/discovery-studio-visualizer-download	Analyse des interactions ligand-résidu en deux dimensions.
GEO2R	NCBI Gene Expression Omnibus	Outil web ; accédé de janvier à mai 2026 ; https://www.ncbi.nlm.nih.gov/geo/geo2r/	Analyse différentielle d'expression de GSE36701.
GeneCards	Institut Weizmann des sciences	Base de données web ; accédée de janvier à mai 2026 ; https://www.genecards.org/	Vérification des symboles génétiques et des informations génétiques lors de la standardisation des cibles.
gmx_MMPBSA	Open source (Valdés-Tresanco et al.)	1.5.x; https://valdes-tresanco-ms.github.io/gmx_MMPBSA/	Estimation de l'énergie libre de liaison par MM-PBSA et décomposition par résidu.
GROMACS	Équipe de développement de GROMACS / open source	2024.2; https://www.gromacs.org/	Moteur de simulation de dynamique moléculaire.
Jeu de données transcriptomiques GSE36701	NCBI Gene Expression Omnibus	GSE36701; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE36701	Jeu de données public d'expression de la muqueuse rectale dans le SCI-C.
Open Babel	Open source	3.2.0; https://openbabel.org/	Conversion de formats de fichiers chimiques, génération de ligands tridimensionnels et préparation des ligands.
Base de données OPM	Université du Michigan	Base de données web ; accédée de janvier à mai 2026 ; https://opm.phar.umich.edu/	Coordonnées d'orientation des protéines dans les membranes utilisées pour aligner HTR2A.
ParmEd	Développeurs ParmEd / open source	4.x; https://parmed.github.io/ParmEd/html/index.html	Repartitionnement des masses des hydrogènes et traitement des topologies de simulation moléculaire.
PyMOL	Schrödinger / open source	2.x; https://www.pymol.org/	Visualisation structurale tridimensionnelle et figures récepteur-ligand.
Banque de données RCSB Protein Data Bank	RCSB PDB	Base de données web ; accédée de janvier à mai 2026 ; https://www.rcsb.org/	Source des structures protéiques expérimentales et des coordonnées PDB.
STITCH	Consortium STITCH (EMBL)	v5.0; https://stitch.embl.de/	Prédiction des cibles d'interaction chimique-protéine.
STRING	Consortium STRING / ELIXIR	v12.0; https://version-12-0.string-db.org/	Construction de réseaux d'interaction protéine-protéine et enrichissement des voies.
SwissADME	Institut suisse de bioinformatique SIB / Université de Lausanne	Outil web ; accédé de janvier à mai 2026 ; https://www.swissadme.ch/	Descripteurs chémo-informatiques, prédictions pharmacocinétiques et évaluation des PAINS.
SwissTargetPrediction	Institut suisse de bioinformatique SIB / Université de Lausanne	Outil web ; accédé de janvier à mai 2026 ; https://www.swisstargetprediction.ch/	Prédiction basée sur les ligands des cibles protéiques humaines.
Cartographie d'identifiants UniProt	Consortium UniProt	Service web ; accédé de janvier à mai 2026 ; https://www.uniprot.org/id-mapping	Association des identifiants protéiques à des symboles génétiques normalisés approuvés par HGNC.
NVIDIA RTX 3080	NVIDIA Corporation	RTX 3080; ≥8 Go de mémoire vidéo ; version CUDA/pilote non précisée dans le manuscrit	Unité de traitement graphique compatible CUDA utilisée pour les simulations de dynamique moléculaire.
GPU compatible CUDA	NVIDIA Corporation	Version du toolkit CUDA non précisée dans le manuscrit ; ≥8 Go de mémoire vidéo	GPU compatible CUDA avec ≥8 Go de mémoire vidéo ; la station de travail nécessitait également ≥32 Go de RAM et un processeur 6 cœurs.
Ubuntu Linux	Canonical Ltd. / open source	22.04 LTS	Système d'exploitation Linux 64 bits.
Python 3.9	Fondation Python	3.9	Environnement de programmation généraliste utilisé pour les scripts de flux de travail et l'analyse.
Gene Expression Omnibus (GEO)	NCBI / Bibliothèque nationale de médecine des États-Unis	Référentiel web public ; aucune version logicielle précisée dans le manuscrit	Référentiel public de données de génomique fonctionnelle.
AutoDockTools/MGLTools	Laboratoire de graphisme moléculaire, Scripps Research	1.5.7	Boîte à outils de préparation des structures moléculaires et des entrées d'accolement.
Outils d'analyse GROMACS	Équipe de développement de GROMACS / open source	2024.2	Utilitaires d'analyse des trajectoires de dynamique moléculaire.
Protocole en six étapes CHARMM-GUI	CHARMM-GUI / Lehigh University	Protocole web ; version exacte non récupérable	Flux de travail web multistade de préparation et d'équilibration du système moléculaire.
cgenff_charmm2gmx_py3.py	Script de conversion open source ; source non précisée dans le manuscrit	Version non précisée dans le manuscrit	Script de conversion de topologie de champ de forces.
Python	Fondation Python	3.9	Environnement de programmation généraliste.
SciPy	Communauté SciPy / open source	Version non précisée dans le manuscrit	Bibliothèque de calcul scientifique.
scipy.stats.fisher_exact	Communauté SciPy / open source	Version SciPy non précisée dans le manuscrit	Implémentation du test exact de Fisher’s.
R	Fondation R pour le calcul statistique	4.3.x	Environnement de calcul statistique.
Bioconductor	Projet Bioconductor / open source	3.18	Infrastructure logicielle en bioinformatique.
limma	Projet Bioconductor / open source	Version non précisée dans le manuscrit	Pack d'analyse d'expression génique différentielle.
Procédure de Benjamini–Hochberg	Méthode statistique	Non applicable (procédure statistique)	Méthode de correction du taux de fausses découvertes.
NVIDIA RTX 3080	NVIDIA Corporation	RTX 3080; ≥8 Go de mémoire vidéo ; version CUDA/pilote non précisée dans le manuscrit	Unité de traitement graphique disposant d'au moins 8 Go de mémoire vidéo.
GPU compatible CUDA	NVIDIA Corporation	Version du toolkit CUDA non précisée dans le manuscrit ; ≥8 Go de mémoire vidéo	Unité de traitement graphique prenant en charge le calcul parallèle à usage général.
Ubuntu Linux 22.04 LTS	Canonical Ltd. / open source	22.04 LTS	Système d'exploitation Linux 64 bits.
TIP3P	Développeurs du champ de forces CHARMM / open source	TIP3P ; aucune version logicielle applicable	Modèle d'eau explicite à trois sites.
MM/PBSA	Développeurs de gmx_MMPBSA / open source	gmx_MMPBSA 1.5.x	Méthode d'énergie de liaison mécanique moléculaire/Poisson–Boltzmann surface-area.