Un flujo de trabajo computacional para priorizar los genes del huésped asociados a metabolitos microbianos en el síndrome del intestino irritable con predominio de estreñimiento

Resumen

Este protocolo integra la predicción de dianas de metabolitos microbianos, la transcriptómica de la mucosa rectal, la interacción proteína-proteína y la enriquecimiento de vías, el acoplamiento molecular, la simulación de dinámica molecular y la estimación de la energía libre de unión mediante mecánica molecular/superficie de Poisson-Boltzmann (MM-PBSA) para generar una lista corta clasificada, con fines generadores de hipótesis, de genes huésped candidatos asociados a metabolitos y complejos proteína-ligando con prioridad estructural para su validación experimental.

Resumen

No existe una canalización computacional estandarizada para priorizar sistemáticamente genes del huésped asociados a metabolitos microbianos y complejos proteína-ligando a partir de bases de datos públicas disponibles sobre aspectos químicos, genómicos y estructurales. En este artículo se describe un flujo de trabajo de ocho etapas que acepta un conjunto definido por el usuario de metabolitos derivados de la microbiota intestinal y genera una lista corta clasificada de genes candidatos del huésped asociados a metabolitos, vías biológicas enriquecidas y complejos proteína-ligando priorizados estructuralmente, listos para su validación experimental. La canalización integra: (i) perfilado quimioinformático de metabolitos; (ii) predicción de dianas candidatas mediante múltiples bases de datos, utilizando herramientas de predicción de dianas basadas en interacciones proteína-químico y en ligandos, así como un programa de acoplamiento molecular; (iii) análisis de expresión génica diferencial de datos transcriptómicos disponibles públicamente; (iv) superposición entre dianas y genes diferencialmente expresados; (v) construcción de redes de interacción proteína-proteína y enriquecimiento de vías; (vi) acoplamiento molecular con un programa de acoplamiento molecular; (vii) simulación de dinámica molecular de 200 ns utilizando un motor de dinámica molecular con un campo de fuerza proteico empleado en simulaciones de dinámica molecular; y (viii) estimación de la energía libre de unión mediante MM-PBSA. Como ejemplo desarrollado, se procesaron nueve metabolitos derivados o modificados por la microbiota intestinal, representativos de ácidos grasos de cadena corta, ácidos biliares, metabolitos derivados del triptófano y urolitina A, utilizando el conjunto de datos transcriptómicos públicos de la mucosa rectal en IBS-C GSE36701. El flujo de trabajo clasificó 17 genes únicos asociados a metabolitos predichos que mostraron expresión diferencial en este conjunto de datos. Los análisis de acoplamiento molecular, simulación de dinámica molecular y MM-PBSA priorizaron estructuralmente cinco complejos metabolito-proteína: ácido litocólico-VDR, ácido litocólico-NR1H4/FXR, ácido ursodesoxicólico-NR1H4/FXR, triptamina-HTR2A (simulada en una bicapa lipídica explícita de 1-palmitoil-2-oleoil-sn-glicerol-3-fosfocolina (POPC)) y urolitina A-CASP3. El protocolo está diseñado para adaptarse a otros conjuntos de metabolitos, conjuntos de datos transcriptómicos de enfermedades y clases de dianas; todas las salidas son predicciones computacionales generadoras de hipótesis que requieren replicación transcriptómica independiente, validación a nivel proteico y ensayos funcionales de respuesta al ligando antes de que puedan extraerse conclusiones causales o terapéuticas.

Introducción

El síndrome del intestino irritable con estreñimiento (SII-E) es un trastorno gastrointestinal funcional prevalente caracterizado por dolor abdominal recurrente, alteraciones en los hábitos intestinales, distensión abdominal y estreñimiento, con una prevalencia mundial estimada de aproximadamente el 10–15 % de la población general^1,2. Las terapias farmacológicas actuales, incluidos secretagogos, procinéticos y antiespasmódicos, pueden mejorar los síntomas individuales en un subgrupo de pacientes; sin embargo, la respuesta al tratamiento sigue siendo heterogénea y rara vez se logra una remisión duradera, lo que refleja la compleja patobiología multifactorial de esta afección^1,3,4. Por lo tanto, se requiere una comprensión mecanicista más completa de cómo se transducen a nivel mucoso las señales microbianas intestinales para generar hipótesis comprobables sobre nuevos objetivos terapéuticos.

La microbiota intestinal contribuye a la homeostasis gastrointestinal inferior mediante la producción y biotransformación de metabolitos químicamente diversos, incluyendo ácidos grasos de cadena corta (AGCC), ácidos biliares secundarios, compuestos derivados del triptófano y metabolitos derivados de polifenoles como las urolitinas^5,6,7,8. Estas moléculas se comunican con las células del huésped a través de un amplio y aún incompletamente caracterizado repertorio de dianas moleculares que va más allá de los receptores de membrana canónicos sensibles a metabolitos, e incluye receptores nucleares, enzimas citosólicas, proteínas que modifican histonas, precursores de hormonas peptídicas y proteínas intracelulares de señalización⁹. Se ha documentado una composición alterada de la comunidad microbiana intestinal y perfiles modificados de metabolitos en pacientes con síndrome del intestino irritable (SII), lo que proporciona una justificación biológica para investigar si los genes del huésped asociados con la respuesta a metabolitos microbianos están transcripcionalmente alterados en la mucosa rectal de pacientes con SII-C¹⁰.

El panel de nueve metabolitos se definió a priori para proporcionar un conjunto compacto, químicamente diverso y biológicamente interpretable de pequeñas moléculas derivadas de la microbiota intestinal o modificadas por esta. La selección se basó en cinco criterios: representación de las principales clases de metabolitos microbianos implicados en la señalización huésped-microbiota; exposición conocida o plausible en la mucosa intestinal distal; disponibilidad de identificadores inequívocos en PubChem y estructuras canónicas; tamaño molecular y trazabilidad estructural adecuados para la predicción y acoplamiento de dianas basados en ligandos; y plausibilidad previa de participación en señalización epitelial, neuroinmune, enteroendocrina, mediada por receptores nucleares o relacionada con la motilidad en el SII-C. El panel seleccionado incluyó butirato y propionato como ácidos grasos de cadena corta (AGCC); ácido quenodesoxicólico, ácido litocólico y ácido ursodesoxicólico como ácidos biliares; triptamina, ácido indol-3-propiónico y ácido indol-3-láctico como metabolitos derivados del triptófano; y urolitina A como metabolito de polifenoles derivado de la microbiota intestinal^5,6,7,8,9,10.

La mayoría de las investigaciones computacionales y experimentales previas han examinado pares individuales de metabolito–receptor o metabolito–enzima de forma aislada, un enfoque que no captura la naturaleza distribuida y convergente de la señalización por metabolitos microbianos a través de las vías del huésped^9,11. La integración de múltiples etapas analíticas confiere un poder de filtrado mutuamente reforzante que ninguna etapa individual puede proporcionar por separado. La predicción computacional de dianas frente a bases de datos curadas produce un conjunto amplio de proteínas candidatas del huésped para cada metabolito. La intersección con datos transcriptómicos relevantes para la enfermedad filtra sustancialmente este conjunto, conservando únicamente los candidatos cuyos transcritos se alteran en el contexto de la enfermedad. Los análisis de enriquecimiento de vías y de redes de interacción proteína–proteína luego asignan la lista reducida de candidatos a módulos biológicos conocidos. El acoplamiento molecular proporciona una evaluación computacional inicial de la complementariedad del bolsillo de unión para cada complejo candidato, y una simulación adicional de dinámica molecular (MD) de 200 ns con descomposición de la energía libre de unión MM-PBSA aporta una dimensión termodinámica resuelta en el tiempo a la priorización estructural, que no está disponible únicamente a partir de las puntuaciones de acoplamiento. Realizar cada paso de forma independiente, sin integración sistemática y filtrado secuencial, daría como resultado listas de candidatos demasiado amplias para ser manejables experimentalmente y no permitiría detectar la arquitectura convergente de las vías.

Dentro del marco de todo este protocolo, por el término «gen asociado a metabolitos» (GAM), entendemos un gen humano cuyo producto proteico ha sido propuesto como un blanco molecular putativo de uno o más metabolitos derivados de la microbiota intestinal por al menos una base de datos computacionales curada de predicciones, y cuyo transcrito se expresa de forma diferencial en el conjunto de datos transcriptómicos relevante para la enfermedad utilizado para demostrar el flujo de trabajo. Esta definición operativa incluye deliberadamente receptores de membrana y receptores nucleares, enzimas citosólicas, proteínas de señalización, precursores de hormonas peptídicas y otras proteínas intracelulares. La designación GAM no constituye evidencia experimental de que un metabolito se una, forme un complejo proteína-ligando, active un receptor, cambie la abundancia proteica o cause enfermedad, sino una propuesta generada computacionalmente con fines hipotéticos que requiere validación experimental.

Este protocolo describe el flujo de trabajo computacional completo de ocho etapas (Figura 1) con suficiente detalle operativo para permitir su réplica independiente, su adaptación a otros paneles de metabolitos o conjuntos de datos de enfermedades, y su extensión a otros contextos de interacción huésped-microbiota. El flujo de trabajo se define explícitamente como un marco para generación de hipótesis y priorización estructural que opera exclusivamente con recursos públicos de ómica y estructurales, y que no infiere concentraciones alteradas de metabolitos, estados de activación de receptores, cambios en la expresión proteica, actividad de señalización downstream ni significancia clínica a partir únicamente de los resultados computacionales. Aquí, demostramos el protocolo como un ejemplo práctico utilizando nueve metabolitos derivados o modificados por la microbiota intestinal y el conjunto de datos transcriptómicos públicos de la mucosa rectal en IBS-C GSE36701, con el objetivo de identificar MAGs y priorizar complejos metabolito-proteína para seguimiento experimental posterior.

Protocolo

El análisis utilizó únicamente datos transcriptómicos públicos y desidentificados de GSE36701, así como bases de datos químicas, de proteínas y estructurales de acceso público. Las bases de datos se consultaron entre enero y mayo de 2026. Cualquier fecha de acceso posterior se documentó en la Tabla separada de materiales.

1. Diseño del estudio, requisitos de hardware y software

1. Defina el flujo de trabajo antes de comenzar el análisis. Utilice ocho etapas: selección de metabolitos, predicción de dianas, análisis de expresión diferencial, solapamiento de genes diferencialmente expresados (DEG), enriquecimiento de interacciones proteína-proteína (PPI)/vías metabólicas, acoplamiento molecular, simulación de dinámica molecular (MD) y estimación MM-PBSA.
2. Registre que el acoplamiento molecular, la dinámica molecular (MD) y la estimación MM-PBSA son únicamente herramientas de priorización estructural. No interprete estos resultados como evidencia experimental de unión, activación del receptor, cambio en la abundancia proteica, eficacia terapéutica o causalidad de enfermedad.
3. Confirme el hardware de cómputo antes de ejecutar simulaciones de dinámica molecular (MD). Utilice un sistema operativo Linux de 64 bits, una CPU de 6 núcleos o superior, una plataforma con aceleración GPU con al menos 8 GB de VRAM, o una plataforma equivalente con aceleración GPU con al menos 8 GB de VRAM, un mínimo de 32 GB de RAM y al menos 200 GB de almacenamiento libre por cada sistema de dinámica molecular.
4. Registre el software principal: un motor de dinámica molecular, un programa de acoplamiento molecular de ligandos metabolitos a proteínas diana, conversión de formatos de archivos químicos¹², generación tridimensional de ligandos, preparación de ligandos, un kit de herramientas para preparar entradas de acoplamiento molecular, un programa de acoplamiento molecular de ligandos metabolitos a proteínas diana, un entorno de programación de propósito general, un entorno de computación estadística con marco de software de bioinformática y un paquete de análisis de expresión génica diferencial.
5. Registre las herramientas de análisis estructural: herramienta web para la construcción de sistemas de membrana, servicio de parametrización de ligandos compatible con CHARMM, herramienta de cálculo de energía de unión mediante mecánica molecular/disolvente continuo, biblioteca de conversión de topología y parámetros moleculares, programa de visualización molecular tridimensional y herramienta de visualización molecular y diagramas de interacción bidimensionales 2021 (véase la Tabla de Materiales para enlaces de descarga e información de versiones).
6. Registre en la Tabla de Materiales/archivo de entorno los identificadores exactos del campo de fuerza proteico utilizado en las simulaciones de dinámica molecular, CGenFF, CHARMM-GUI, marco de software R/bioinformática y la herramienta de cálculo de energía de unión mediante mecánica molecular/disolvente continuo. Marque los identificadores ausentes como “no recuperable”; no los infiera.

2. Selección de metabolitos y caracterización quimioinformática

1. Defina el panel de metabolitos antes de la predicción de objetivos. Incluya butirato (PubChem CID: 264), propionato (CID: 1032), ácido quenodesoxicólico (CID: 10133), ácido litocólico (CID: 9903), ácido ursodesoxicólico (CID: 31401), triptamina (CID: 1150), ácido indol-3-propiónico (CID: 3744), ácido indol-3-láctico (CID: 92904) y urolitina A (CID: 5488186).
2. Recupere el sistema de entrada simplificado de líneas moleculares canónicas (SMILES) y el CID de PubChem para cada metabolito. Verifique los sinónimos y estructuras duplicadas antes de la predicción de objetivos. Almacene los identificadores finales en la hoja maestra de metabolitos.
3. Envíe las cadenas SMILES canónicas a la herramienta web de predicción de propiedades fisicoquímicas y ADME¹³ (consulte la Tabla de materiales). Registre el peso molecular, el área superficial polar topológica (TPSA), el logP consenso, donadores de enlaces de hidrógeno, aceptores de enlaces de hidrógeno, enlaces rotatorios, absorción gastrointestinal predicha, predicción de glicoproteína P y alertas de Lipinski, Veber, Ghose, Egan, Muegge y PAINS.
4. Mantenga los metabolitos con reconocimiento exitoso de estructura, peso molecular ≤500 Da y sin alertas de PAINS. Registre cualquier criterio fallido y la decisión de mantener o excluir el metabolito.
5. Asigne estados de ionización antes de la predicción de objetivos y el acoplamiento molecular. Utilice carboxilatos desprotonados para el butirato y el propionato, formas neutras de ácidos carboxílicos para los ácidos biliares, amonio protonado para la triptamina y formas neutras para los metabolitos restantes.

3. Predicción de dianas humanas candidatas

1. Abra un programa de predicción de interacciones químico-proteína¹⁴ (consulte la Tabla de materiales). Ingrese el nombre de cada metabolito o el CID de PubChem, seleccione Homo sapiens (ID taxonómico: 9606) y establezca el puntaje mínimo combinado de interacción en ≥0,700.
2. Priorice los canales de evidencia experimental y de bases de datos curadas en la predicción de interacciones químico-proteína. Descargue la tabla completa de asociaciones proteicas para cada metabolito.
3. Abra un programa de acoplamiento molecular¹⁵ (consulte la Tabla de materiales). Envíe cada cadena SMILES canónica con Homo sapiens seleccionado y conserve los blancos con una probabilidad ≥0,70.
4. Combine las salidas de la predicción de interacciones químico-proteína y del programa de acoplamiento molecular como un conjunto unión para cada metabolito. Mantenga cualquier blanco que cumpla con al menos uno de los umbrales de la base de datos y elimine las entradas duplicadas exactas de símbolos génicos.
5. Normalice las entradas proteicas a los símbolos génicos aprobados por el Comité de Nomenclatura Genética HUGO (HGNC) mediante la asignación de identificadores proteicos a símbolos génicos estandarizados HGNC o a la base de datos integrada de información genética humana (consulte la Tabla de materiales). Resuelva los alias, símbolos obsoletos y anotaciones de isoformas a un único símbolo génico por proteína.
6. Clasifique cada blanco como receptor de membrana, receptor nuclear, enzima, proteína de señalización intracelular, hormona peptídica, proteína relacionada con hormonas u otra proteína intracelular. Registre la clase en la tabla de blancos.

4. Conjunto de datos transcriptómicos y análisis de expresión diferencial de genes

1. Acceda a GSE36701 a través de la herramienta de análisis de expresión génica diferencial basada en web del NCBI^16,17 (consulte la Tabla de Materiales). Registre que el conjunto de datos incluye datos de expresión de biopsias de mucosa rectal de grupos con síndrome del intestino irritable con estreñimiento (IBS-C), síndrome del intestino irritable con predominio de diarrea (IBS-D), síndrome del intestino irritable postinfeccioso y voluntarios sanos¹⁸.
2. Busque en GEO y ArrayExpress una cohorte independiente para validación. Utilice combinaciones de los términos síndrome del intestino irritable con estreñimiento (IBS-C), síndrome del intestino irritable con predominio de estreñimiento, mucosa rectal, mucosa cólica, biopsia, transcriptoma, microarreglo y RNA-seq. Registre los repositorios, los términos de búsqueda, la fecha de búsqueda y si se identificó un conjunto de datos de validación comparable.
3. Inicie la herramienta web de análisis de expresión génica diferencial desde el registro de GSE36701 (consulte la Tabla de Materiales). Asigne las 18 muestras de IBS-C al grupo de IBS-C, asigne los 40 voluntarios sanos al grupo control y deje sin asignar las muestras de IBS-D y de síndrome del intestino irritable postinfeccioso.
4. Realice el análisis de expresión diferencial utilizando el entorno del paquete de análisis de expresión génica diferencial con corrección de la tasa de falsos descubrimientos (FDR) de Benjamini-Hochberg¹⁹. Descargue la tabla completa de resultados con el ID del sonda, símbolo del gen, título del gen, logFC, AveExpr, estadístico t moderado, valor P crudo y valor P ajustado.
5. Consolide las sondas en entradas a nivel de gen. Elimine las sondas que carezcan de símbolos de gen; conserve la sonda con el FDR más bajo para los símbolos duplicados; y utilice el logFC absoluto mayor como criterio de desempate.

5. Análisis de superposición entre diana y degradación y evaluación estadística

1. Intersecte cada lista predicha específica de metabolitos con la lista de DEG a nivel génico con FDR < 0,05. Registre los genes superpuestos, el metabolito de origen, logFC, valor P ajustado y la dirección de expresión.
2. Combine las listas de superposición específicas de metabolitos en una lista no redundante de MAG. Cuente el número total de objetivos predichos, superposiciones específicas de metabolitos y el número total de MAG únicos.
3. Evalue la consistencia direccional a nivel de sonda para genes con múltiples sondas. Marque cualquier gen cuyas sondas muestren desacuerdo en la dirección de expresión.
4. Construya la tabla de contingencia del test exacto de Fisher utilizando el número total de entradas colapsadas por gen, el número total de DEG, el número total de objetivos predichos únicos y los MAG observados. Calcule el valor P unilateral, la razón de momios y el intervalo de confianza del 95 % mediante la implementación del test exacto de Fisher.
5. Si la tasa de DEG en el fondo supera el 50 %, informe la superposición como descriptiva y no como un enriquecimiento validado independientemente. Considere la regulación uniforme hacia abajo como un patrón direccional descriptivo, a menos que se realice una prueba de direccionalidad independiente.

6. Análisis de la red de interacciones proteína-proteína y enriquecimiento de vías

1. Envíe la lista completa y única de MAG al análisis de red de interacción proteína-proteína y enriquecimiento de vías²⁰ (consulte la Tabla de materiales). Seleccione Homo sapiens y establezca el puntaje mínimo de interacción en 0,700.
2. Exporte la red combinada de interacción proteína-proteína y enriquecimiento de vías, junto con la tabla completa de interacciones. Si la minería de texto produce una topología artificialmente densa, desactive la minería de texto y conserve los canales experimentales, de coexpresión y de bases de datos.
3. Genere subredes por clase de metabolito para los MAG asociados a ácidos grasos de cadena corta, asociados a ácidos biliares y asociados a triptamina/serotoninérgicos. Utilice el mismo organismo y configuración de confianza en el análisis de red de interacción proteína-proteína y enriquecimiento de vías.
4. Realice el análisis de red de interacción proteína-proteína y enriquecimiento de vías frente a Kyoto Encyclopedia of Genes and Genomes (KEGG)²¹, Reactome²² y Gene Ontology (GO) Biological Process^23,24. Aplique una tasa de falsos descubrimientos (FDR) de Benjamini–Hochberg (BH) <0,05 y exporte todas las tablas de enriquecimiento.

7. Acoplamiento molecular

1. Recuperar estructuras de receptores determinadas experimentalmente desde la base de datos RCSB PDB (Research Collaboratory for Structural Bioinformatics Protein Data Bank)²⁵ (véase la Tabla de materiales). Utilizar VDR/1DB1, NR1H4/FXR/3DCT, CASP3/2DKO y HTR2A/6A93 para los cinco complejos proteína-ligando prioritarios.
2. Preparar cada receptor conservando la cadena A y eliminando moléculas de agua, ligandos co-cristalizados, cofactores, iones y registros HETATM no proteicos. Para el 6A93, eliminar el segmento de fusión de la lisozima T4 antes de la preparación del receptor.
3. Añadir hidrógenos polares, asignar cargas de Gasteiger y guardar cada receptor en formato PDBQT mediante el kit de herramientas de estructura molecular. Inspeccionar los estados de protonación de las histidinas en el sitio de unión antes de la conversión a PDBQT y documentar los estados seleccionados.
4. Generar la estructura tridimensional de cada ligando utilizando la herramienta de conversión de archivos de estructura química. Minimizar la energía con el campo de fuerza universal (UFF) durante 500 pasos, asignar el estado de ionización a pH 7,4, asignar cargas de Gasteiger y guardar en formato PDBQT.
5. Definir una caja de acoplamiento de 25 Å x 25 Å x 25 Å centrada en el centroide del ligando co-cristalizado. Utilizar los centros (10, 19, 33) para VDR, (137, 31, 78) para FXR, (37, 34, 32) para CASP3 y (12, −1, 61) para HTR2A.
6. Realizar un acoplamiento molecular de ligandos metabolitos a las proteínas diana^26,27 con exhaustividad = 8, semilla = 42, num_modes = 9 y energy_range = 3 kcal/mol. Registrar la puntuación Vina mejor clasificada y los valores de desviación cuadrática media (RMSD) para todas las conformaciones.
7. Seleccionar el modo 1 para cada complejo prioritario. Generar diagramas bidimensionales de ligando-residuo mediante una herramienta de visualización molecular y una herramienta de diagrama de interacciones bidimensional, y vistas tridimensionales del complejo receptor-ligando en un programa de visualización molecular tridimensional.
8. Realizar controles de reacoplamiento para VDR/1DB1 y FXR/3DCT. Aceptar la configuración de acoplamiento del receptor cuando el RMSD de los átomos pesados sea <2,0 Å en relación con la conformación cristalográfica.
9. Realizar controles de acoplamiento cruzado acoplando el ácido litocólico (LCA) a CASP3 y la triptamina a VDR. Comparar las puntuaciones cognadas y no cognadas y registrar los casos en los que la diferencia de puntuación sea <1,0 kcal/mol.

8. Simulación de dinámica molecular

1. Genere parámetros del ligando con el servicio de parametrización de ligandos compatible con CHARMM²⁸ (consulte la Tabla de materiales). Revise todas las puntuaciones de penalización y marque cualquier parámetro con una penalización >50.
2. Convierta los archivos de flujo del ligando a archivos .itp y .prm compatibles con el motor de dinámica molecular mediante el script de conversión de topología del campo de fuerza. Combine los archivos de topología del ligando y de la proteína para cada complejo.
3. Aplique la redistribución de masa de hidrógeno utilizando la biblioteca de conversión de parámetros y topología molecular. Genere topologías acuosas utilizando el campo de fuerza de la proteína empleado en las simulaciones de dinámica molecular²⁹ y un modelo explícito de agua de tres sitios³⁰.
4. Solvate los complejos acuosos en una caja dodecaédrica con un espacio libre entre el soluto y el borde de al menos 1,2 nm. Neutralice los sistemas y agregue NaCl hasta alcanzar una concentración de 0,15 M.
5. Construya el sistema de membrana triptamina-HTR2A utilizando una herramienta web para la construcción de sistemas de membrana^31,32,33 (consulte la Tabla de materiales). Utilice coordenadas del receptor alineadas con la base de datos de orientación de proteínas de membrana³⁴ (consulte la Tabla de materiales), una bicapa de POPC pura, capas de agua de 22,5 Å y NaCl 0,15 M.
6. Minimice la energía de todos los sistemas mediante el método de descenso más empinado durante hasta 50.000 pasos. Confirme la convergencia cuando Fmax <1000 kJmol^-1nm^-1 antes del equilibrio.
7. Equilibre los sistemas acuosos mediante los ensambles de número constante de partículas, volumen y temperatura (NVT) y número constante de partículas, presión y temperatura (NPT). Equilibre el sistema de membrana utilizando el flujo de trabajo web multietapa de seis pasos para la preparación y equilibrado de sistemas moleculares, con restricciones liberadas gradualmente.
8. Realice simulaciones de dinámica molecular productivas de 200 ns para los cinco complejos. Utilice un paso de tiempo de 4 fs con redistribución de masa de hidrógeno (HMR), termostato V-rescale a 310 K, barostato Parrinello-Rahman a 1 bar, electrostática por Ewald en malla de partículas (PME)³⁵ y restricciones LINCS³⁶.
9. Analice las trayectorias finales con utilidades de análisis de trayectorias de dinámica molecular. Calcule la desviación cuadrática media (RMSD) del esqueleto, la fluctuación cuadrática media de Cα (RMSF), el radio de giro, el área de superficie accesible al solvente (SASA) y los enlaces de hidrógeno entre proteína y ligando, utilizando los últimos 150 ns como ventana principal de análisis.

9. Estimación de la energía libre de enlace mediante MM-PBSA

1. Extraiga instantáneas de la trayectoria para el análisis MM-PBSA. Utilice 2 001 fotogramas para cada complejo acuoso y 201 fotogramas procesados para el subsistema de HTR2A embebido en la membrana.
2. Realice el cálculo de energía de unión mediante la herramienta de mecánica molecular/disolvente continuo³⁷ con solvatación de Poisson-Boltzmann, constante dieléctrica interna = 1, constante dieléctrica externa = 80, solvatación no polar basada en SASA y sin corrección de entropía. Informe la energía libre media de unión y la desviación estándar.
3. Realice la descomposición por residuo para los cinco complejos. Informe los residuos estabilizantes e inestabilizantes con contribuciones absolutas ≥0,5 kcalmol⁻¹.

    

Resultados

Objetivos asociados a metabolitos candidatos

Los nueve metabolitos produjeron conjuntos heterogéneos de dianas predichas mediante un programa de predicción de interacciones químico-proteína y un programa de acoplamiento molecular. El propionato, la triptamina, los ácidos biliares y la urolitina A arrojaron varias dianas con relevancia conocida en la señalización gastrointestinal. El paisaje de dianas predichas incluyó receptores de membrana canónicos, receptores nucleares, enzimas intracelulares, proteínas de señalización y proteínas relacionadas con hormonas peptídicas. Por lo tanto, los resultados posteriores se describen como genes asociados a metabolitos (MAGs), en lugar de hallazgos exclusivos de receptores (Tabla 1).

Comparación con interacciones metabolito-proteína reportadas

Para establecer un punto de referencia de la salida de predicción de objetivos frente al conocimiento experimental existente, las relaciones predichas entre metabolitos y objetivos se clasificaron en tres niveles de evidencia: (i) interacciones directas o a nivel de clase cercana entre metabolitos y proteínas respaldadas experimentalmente, en las que se ha informado que el metabolito o un metabolito endógeno estrechamente relacionado se une, activa, inhibe o regula funcionalmente la proteína codificada; (ii) interacciones respaldadas por vías o clases de objetivos, en las que el objetivo predicho pertenece a una vía sensible al metabolito o a una familia de receptores establecida, pero la evidencia directa para el par específico de metabolito-proteína es limitada; y (iii) asociaciones exclusivamente computacionales para las cuales no se identificó ninguna interacción experimental directa en la literatura revisada. Esta evaluación de referencia se utilizó para contextualizar, no para validar, los MAG predichos.

Varias predicciones recapitularon biología previamente reportada. La interacción propionato-FFAR2 se consideró respaldada experimentalmente porque FFAR2/GPR43 es un receptor canónico de ácidos grasos de cadena corta. La interacción butirato-HDAC3 se clasificó como respaldada experimental o conceptualmente porque el butirato es un inhibidor reconocido de la histona desacetilasa, y la superposición predicha implicaba a un miembro de la familia HDAC. Las predicciones asociadas a ácidos biliares que involucraban a NR1H4/FXR y VDR se consideraron respaldadas por la biología establecida de los receptores nucleares de ácidos biliares, particularmente para ácidos biliares hidrófobos como el LCA; las predicciones de FXR asociadas al ácido ursodesoxicólico (UDCA) se interpretaron con cautela porque el UDCA generalmente es un ligando de FXR más débil o dependiente del contexto. Las predicciones asociadas a la triptamina en HTR1B, HTR2A, HTR2B y HTR6 se clasificaron como respaldadas por la vía serotoninérgica, más que como interacciones directas y específicas de receptores confirmadas, porque la triptamina es una monoamina derivada del triptófano microbiano y los receptores de serotonina son reguladores establecidos de la motilidad y secreción gastrointestinales. La interacción urolitina A-CASP3 se consideró respaldada por vías debido a vínculos publicados entre la urolitina A y respuestas apoptóticas o relacionadas con caspasas, pero no por evidencia directa de unión a CASP3. Las interacciones ácido indol-3-láctico-KYAT1 y ácido indol-3-propiónico-KYAT1 se mantuvieron como hipótesis exclusivamente computacionales porque la literatura más amplia apoya la señalización en el huésped por derivados microbianos del indol, pero no la unión directa de estos metabolitos específicos a KYAT1^7,8,38,39,40.

Por consiguiente, Tabla 1 distingue la nominación computacional de dianas del nivel de respaldo experimental previo o de vía metabólica. También proporciona, para cada diana, la fuente de predicción (una predicción de interacción químico-proteína, un programa de acoplamiento molecular, o ambos), la puntuación combinada de interacción para la predicción de interacción químico-proteína y la probabilidad del programa de acoplamiento molecular cuando la diana fue identificada mediante un programa de acoplamiento molecular. Las dianas predichas sin evidencia experimental previa directa se describen como genes candidatos asociados a metabolitos que requieren validación independiente a nivel proteico y de respuesta al ligando.

Solapamiento entre los objetivos predichos y los genes diferencialmente expresados en IBS-C

La intersección entre las listas combinadas de uniones de objetivos predichos y los resultados de expresión diferencial a nivel génico identificó 17 genes únicos asociados a metabolitos predichos que presentaron una expresión diferencial significativa en la comparación entre pacientes con IBS-C y voluntarios sanos. Los 17 genes estaban regulados a la baja. Este conjunto incluyó receptores de membrana y nucleares (CASR, FFAR2, GPR68, HTR1B, HTR2A, HTR2B, HTR6, NR1H4, TBXA2R, VDR) y proteínas no receptoras (CASP3, GCG, GNAQ, GPHN, HDAC3, KYAT1, MLN) (Tabla 1, Figura 2A,B).

Todos los 17 MAG cumplieron un umbral de tasa de descubrimiento falso (FDR) inferior a 0,05; 16 de los 17 cumplieron el criterio más estricto de FDR < 0,001, mientras que el gen restante (HTR1B) fue significativo con FDR < 0,05. Siete de los 17 blancos (CASP3, GCG, GNAQ, GPHN, GPR68, HDAC3, TBXA2R) satisficieron ambos criterios: FDR < 0,001 y un cambio absoluto en el log₂ superior a 1,0 (rango de logFC de −1,34 a −1,10), lo que indica una regulación descendente fuerte y consistente para este subconjunto. Los blancos restantes mostraron una regulación descendente moderada pero estadísticamente significativa (|logFC| entre 0,45 y 0,97). Este patrón descriptivo uniforme se interpretó con cautela, considerando las características de expresión a nivel del genoma del conjunto de datos (véase la evaluación estadística más abajo).

Evaluación estadística de la superposición entre el blanco y los genes diferencialmente expresados

Para evaluar formalmente la significancia estadística del solapamiento de 17 genes, se aplicó una prueba exacta de Fisher unilateral utilizando los 17 genes diana predichos como conjunto de consulta y las 18.296 entradas únicas colapsadas por gen detectadas en GSE36701 como fondo genómico. De este fondo, 17.296 genes (94,5 %) se expresaron diferencialmente con FDR < 0,05, lo que refleja una supresión transcripcional casi universal en la comparación del tejido rectal de IBS-C. Los 17 genes diana predichos se encontraron entre los genes con expresión diferencial (solapamiento observado 17/17, 100 %). Dado el porcentaje de expresión diferencial del fondo (94,5 %), el solapamiento esperado para cualquier conjunto aleatorio de 17 genes es de 16,1 genes. La prueba exacta de Fisher arrojó un valor de p = 0,384 con una razón de momios corregida por continuidad de 2,03 (intervalo de confianza del 95 %: 0,12–33,73), lo cual no fue estadísticamente significativo a α = 0,05 (Figura 3A–C).

Este resultado indica que el solapamiento observado de 17/17 no supera el solapamiento esperado por azar según el perfil de expresión genómico de este conjunto de datos. Por consiguiente, estos hallazgos se interpretan como un patrón descriptivo direccional, en el que todos los 17 objetivos predichos se regulaban a la baja de manera consistente y significativa en el tejido de mucosa rectal de IBS-C, más que como evidencia de enriquecimiento estadístico o validación independiente frente a un fondo genómico. La prueba formal de enriquecimiento requeriría replicación en conjuntos de datos transcriptómicos con perfiles de expresión diferencial más selectivos, en los que significativamente menos de la mitad de todos los genes alcancen significancia. Debe enfatizarse que la regulación a la baja uniforme de los 17 genes solapados es una observación descriptiva y no un resultado estadístico validado por separado, ya que el fondo de expresión diferencial de este conjunto de datos también está predominantemente regulado a la baja; por tanto, una dirección descendente compartida entre los genes solapados era esperada y no se sometió a una prueba formal de direccionalidad. Esta dirección uniforme no debe interpretarse, por lo tanto, como evidencia estadística independiente de una regulación coordinada y específica del metabolito.

Patrones específicos de metabolitos

El propionato tuvo el mayor número de genes superpuestos, incluidos CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN y TBXA2R, lo que sugiere una posible participación de la señalización asociada a ácidos grasos de cadena corta y a la proteína Gq. El butirato mostró superposición con HDAC3, lo cual es consistente con la biología del butirato relacionada con la desacetilasa de histonas, aunque la disminución del ARNm por sí sola no establece una respuesta alterada al butirato. Las superposiciones asociadas con ácidos biliares incluyeron los receptores nucleares VDR y NR1H4, ambos efectores reconocidos de la señalización de ácidos biliares en el intestino^38,39. La triptamina mostró superposición con HTR1B, HTR2A, HTR2B y HTR6, implicando a la señalización serotoninérgica como un módulo candidato, un sistema con funciones bien establecidas en la motilidad y secreción gastrointestinales⁴⁰. El ácido indol-3-láctico y el ácido indol-3-propiónico mostraron superposición con KYAT1, y la urolitina A mostró superposición con CASP3.

Enriquecimiento de vías

El análisis de enriquecimiento funcional de los 17 genes superpuestos identificó vías relacionadas con la señalización downstream del receptor acoplado a proteína G (GPCR), la señalización de Gαq, la unión de ligandos al GPCR, la sinapsis serotoninérgica, la interacción entre ligandos neuroactivos y receptores, la transducción de señales de calcio, la señalización de AMPc y la secreción de hormonas peptídicas. Estos resultados son coherentes con la composición del conjunto de genes y respaldan su coherencia biológica, pero reflejan la anotación funcional de los genes enviados más que evidencia independiente de actividad a nivel de vías.

Estructura de la red de interacciones proteína-proteína

La construcción de la red de interacciones proteína-proteína y el análisis de enriquecimiento de vías se interpretaron en tres redes complementarias. En la red meta combinada de 17 genes (Red 1), la estructura con mayor respaldo por anotaciones fue un componente de señalización GPCR/Gαq centrado en GNAQ que conecta GNAQ con genes asociados a receptores, incluyendo TBXA2R, CASR, HTR2A y HTR2B. También se mantuvo una conectividad limitada entre receptores de serotonina, más notablemente entre HTR2A y HTR2B, mientras que varios otros genes permanecieron aislados o débilmente conectados en el umbral de confianza seleccionado. La red específica de propionato (Red 2) mostró una topología más restringida, con GNAQ conservando vínculos respaldados por anotaciones con CASR y TBXA2R, mientras que FFAR2, GPR68, GCG, GPHN y MLN estuvieron aislados o débilmente conectados. La red de triptamina/serotonina (Red 3) incluyó HTR1B, HTR2A, HTR2B y HTR6; dentro de este subconjunto, HTR2A y HTR2B mostraron la conexión principal respaldada por anotaciones, mientras que HTR1B y HTR6 no estuvieron directamente conectados en el umbral elegido (Figura 4A–C).

Acoplamiento molecular

Se realizó un acoplamiento molecular en cinco complejos metabolito-proteína seleccionados. Los pares ácido biliar-receptor nuclear mostraron puntuaciones Vina más favorables que urolitina A-CASP3 y triptamina-HTR2A. LCA-VDR obtuvo la mejor puntuación con −10,0 kcal/mol, seguido por LCA-NR1H4/FXR (−9,9 kcal/mol) y UDCA-NR1H4/FXR (−9,4 kcal/mol). La urolitina A-CASP3 y la triptamina-HTR2A tuvieron puntuaciones más bajas pero aún razonables de −7,1 kcal/mol (Tabla 2).

Para el complejo LCA-VDR (ID de PDB: 1DB1), la posición predicha fue respaldada por un enlace de hidrógeno convencional entre el oxígeno carboxilato de LCA y Ser278 (4,29 Å), junto con extensos contactos hidrofóbicos que involucran a Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233 y His397, y contactos adicionales de van der Waals con Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275 y Phe150. La posición mejor clasificada tuvo una puntuación de Vina de −10,0 kcal/mol, un tamaño de cavidad de 2055 ų y un centro de cuadrícula de (10, 19, 33) (Tabla 3, Figura 5A,B).

Para el complejo LCA-NR1H4/FXR (PDB ID: 3DCT), la puntuación de acoplamiento de −9,9 kcal/mol estuvo acompañada por enlaces de hidrógeno predichos que involucran a His294 e Ile335, una interacción π-Sigma con His294, y contactos hidrofóbicos Alquilo o π-Alquilo que involucran a Met290, Met328, Ala291, Leu287, Ile352 y His447, con contactos adicionales de van der Waals que favorecen la acomodación del armazón esteroideo en el bolsillo de FXR (Tabla 4, Figura 6A,B).

La conformación predicha del complejo UDCA-NR1H4/FXR (PDB ID: 3DCT) mostró un enlace de hidrógeno convencional con His447 (3,66 Å), otro enlace de hidrógeno con Gly322 (3,46 Å), una interacción π-anión con Val325 (4,96 Å) y un enlace carbono-hidrógeno con Trp469 (4,51 Å). El mapa de interacciones también identificó contactos donador-donador desfavorables con Arg395 (3,89 Å) y Gln396 (3,40 Å), lo que sugiere que la puntuación más baja de Vina para UDCA en comparación con LCA en el mismo sitio de unión puede deberse a una geometría local o electrostática menos favorable (Tabla 5, Figura 7A,B).

En el complejo urolitina A-CASP3 (ID de PDB: 2DKO), el modo de unión predicho presentó enlaces de hidrógeno convencionales con Gln161 (3,78 y 4,19 Å), Ser120 (3,95 Å) y Arg207 (3,05 y 3,77 Å), y se estabilizó adicionalmente mediante interacciones π-catión con Arg207, un enlace de hidrógeno π-donador con Cys163, y contactos adicionales π-alquilo y de Van der Waals que involucraron a Arg64, Ala162, His121, Ser205 y Trp206 (Tabla 6, Figura 8A,B).

Para el complejo triptamina-HTR2A (ID PDB: 6A93), la conformación predicha fue estabilizada por un puente salino electrostático entre la amina protonada de la triptamina y Asp155, el aspartato conservado en la hélice transmembrana 3 (D3.32 según la numeración de Ballesteros-Weinstein) que ancla la amina protonada de los ligandos aminérgicos en receptores de serotonina y relacionados^41,42,43, junto con interacciones por enlace de hidrógeno con Thr160 y Ser159, contactos aromáticos con Phe340 y Trp336, e interacciones π-alquílicas con Val156 e Ile163. Contactos adicionales de van der Waals con Tyr370, Phe339, Ser242, Phe243, Phe332 y Leu123 respaldaron un patrón de unión al sitio ortostérico (Tabla 7, Figura 9A,B).

Validación del protocolo de acoplamiento

Para evaluar la fiabilidad del protocolo de acoplamiento, se realizaron dos experimentos de control complementarios. Para los controles de reacoplamiento (positivos), se extrajeron los ligandos co-cristalizados de sus estructuras de referencia por cristalografía de rayos X y se volvieron a acoplar en sus sitios de unión nativos. La pose predicha con mayor puntuación para el análogo de vitamina D, VDX, en VDR/1DB1 presentó una desviación de 0,87 Å respecto a la posición cristalográfica, y el ligando co-cristalizado WAY-362450 en FXR/3DCT mostró una desviación de 1,79 Å; ambos valores se encuentran por debajo del umbral convencional de aceptación de 2,0 Å, lo que respalda la validez geométrica del protocolo de acoplamiento para estos sistemas receptores (Figura 10A,B). Para los controles de acoplamiento cruzado (negativos), se acopló el ácido litocólico en la caspasa-3 (2DKO), una proteasa de cisteína para la cual no es un ligando conocido, obteniéndose una puntuación predicha (−8,3 kcal/mol) 1,7 kcal/mol más débil que en su diana cognata VDR (−10,0 kcal/mol), lo cual es consistente con una selectividad predicha del sitio de unión. El triptamina acoplado en VDR produjo una puntuación predicha de −6,4 kcal/mol en comparación con −7,1 kcal/mol en su diana cognata HTR2A, una diferencia de 0,7 kcal/mol que se encuentra dentro de la incertidumbre informada para las puntuaciones de acoplamiento molecular de ligandos metabolitos a proteínas diana, lo que indica únicamente una selectividad predicha moderada para este ligando más pequeño (Figura 10C). En conjunto, estos controles indican que el protocolo de acoplamiento reproduce geometrías de unión conocidas y discrimina entre pares cognatos y no cognatos bajo las condiciones evaluadas, aunque las predicciones computacionales no sustituyen las mediciones experimentales de afinidad (Tabla 8).

Simulación de dinámica molecular

Se realizaron simulaciones de dinámica molecular para los cinco complejos priorizados durante trayectorias de producción de 200 ns. Los cuatro complejos solubles y de receptores nucleares se simularon en un solvente acuoso explícito, mientras que el complejo triptamina-HTR2A se simuló en una bicapa lipídica POPC explícita para proporcionar un entorno de membrana fisiológicamente adecuado para este receptor acoplado a proteína G. Los análisis evaluaron la estabilidad dinámica de las posiciones acopladas bajo condiciones dependientes del tiempo y permitieron la comparación del comportamiento estructural relativo entre los complejos (Tabla 9).

El perfil de RMSD del complejo LCA-VDR/1DB1 mostró un breve período de equilibrio durante los primeros 10 ns, seguido de una meseta estable, con fluctuaciones principalmente en el rango de 0.20–0.28 nm (Figura 11A). Los valores de RMSF fueron bajos, y las fluctuaciones del armazón fueron < 0.15 nm para la mayoría de los residuos (Figura 11B). El análisis de enlaces de hidrógeno mostró una red persistente de 2–5 enlaces de hidrógeno, con aumentos ocasionales hasta 7 (Figura 11C). El radio de giro (Rg) se mantuvo dentro del rango de 1.25–1.75 nm, y el área de superficie accesible al solvente (SASA) se mantuvo alrededor de 130 nm² (Figura 11D,E).

El complejo urolitina A-CASP3/2DKO mostró una mayor actividad dinámica. La DCR aumentó inicialmente y luego osciló entre 0,4 y 0,7 nm, con un breve evento de alta desviación alrededor de los 165 ns (Figura 12A). El análisis de FCR mostró una alta movilidad a nivel de residuos, con las fluctuaciones más grandes en la región del bucle flexible alrededor del residuo 175 (Figura 12B). El análisis de enlaces de hidrógeno reveló una red inicial extensa de aproximadamente 2 a 5 enlaces durante los primeros 30-40 ns, seguida principalmente por 0 a 2 enlaces intermitentes (Figura 12C). Los perfiles correspondientes del radio de giro y del SASA se muestran en la Figura 12D,E.

Para los sistemas de ácidos biliares NR1H4/FXR (3DCT), el perfil de RMSD del esqueleto permaneció dentro de un rango relativamente estrecho durante la mayor parte de la trayectoria (Figura 13A), mientras que el perfil de RMSF mostró menor movilidad en las regiones centrales y fluctuaciones más altas en las regiones flexibles (Figura 13B). El complejo LCA-3DCT mantuvo aproximadamente de tres a cuatro enlaces de hidrógeno persistentes a lo largo de toda la trayectoria, mientras que el complejo UDCA-3DCT exhibió una mayor fluctuación en los enlaces de hidrógeno y una reducción en la formación de estos después de aproximadamente 125 ns. Los perfiles del radio de giro para los sistemas unidos a LCA y UDCA se muestran en la Figura 13C,D, respectivamente, y los perfiles correspondientes de SASA se muestran en la Figura 13E,F.

Dinámica molecular de la membrana del complejo triptamina-HTR2A

Se simuló el complejo triptamina-HTR2A/6A93 durante 200 ns en un bicapa lipídica explícita de POPC compuesta por 258 moléculas lipídicas, un modelo de agua explícito de tres sitios y NaCl 0,15 M, para un tamaño total del sistema de aproximadamente 100.925 átomos^33,44,45. El receptor permaneció establemente embebido en la bicapa a lo largo de toda la trayectoria (Figura 14). La RMSD del esqueleto aumentó desde aproximadamente 0,10 nm hasta una meseta estable cercana a 0,15–0,20 nm dentro de los primeros 100 ns y se mantuvo estable posteriormente, con todos los valores por debajo de 0,25 nm, lo que indica que el receptor conservó una conformación estable en el entorno de la membrana sin desnaturalización global (Figura 15A). La RMSF por residuo mostró fluctuaciones bajas en el núcleo helicoidal transmembrana con una movilidad mayor esperada en las regiones de bucles y terminales, consistente con la flexibilidad típica de los receptores acoplados a proteínas G (Figura 15B). El radio de giro estuvo estrechamente confinado entre aproximadamente 2,06 y 2,12 nm, y el SASA fluctuó dentro de un intervalo estrecho sin deriva progresiva, lo que confirma la conservación del conjunto transmembrana compacto (Figura 15C,D).

El enlace de hidrógeno entre la proteína y el ligando se mantuvo durante toda la trayectoria (Figura 15E), con fluctuaciones importantes en el número de enlaces de hidrógeno, que variaron entre 1 y 3. Para evaluar específicamente la persistencia de la interacción iónica clave, se monitoreó la distancia mínima entre el nitrógeno del grupo amonio protonado de la triptamina y los átomos de oxígeno carboxilato de Asp155 (D3.32) durante toda la trayectoria. Esta distancia se mantuvo estrechamente distribuida alrededor de un valor medio de 0,270 nm (mínimo 0,247 nm, máximo 0,424 nm), y el contacto del puente salino (< 0,4 nm) se mantuvo durante el 99,9 % de la simulación, con solo dos breves excursiones transitorias y sin ningún evento de disociación sostenida (Figura 16). Estos resultados sugieren que la interacción iónica conservada de Asp155 fue suficiente para estabilizar la triptamina dentro del bolsillo ortostérico de HTR2A durante toda la simulación en la membrana.

Energía libre de unión MM-PBSA y descomposición por residuo

Se realizó un análisis MM-PBSA para añadir una capa adicional de priorización energética a los cinco complejos (Tabla 10). Para los cuatro complejos acuosos, la descomposición por residuo identificó los principales contribuyentes energéticos de cada modo de unión predicho. En el complejo LCA-VDR/1DB1, el ligando y Gln317 tuvieron una contribución favorable, mientras que Trp286 presentó una contribución desfavorable. En el complejo urolitina A-CASP3/2DKO, Arg64 y Arg207 mostraron contribuciones por residuo fuertemente negativas, lo que indica una estabilización polar o electrostática considerable; sin embargo, la trayectoria correspondiente permaneció altamente dinámica, lo que demuestra que una termodinámica favorable a nivel de residuos por sí sola no garantiza una estabilidad sostenida del complejo. Para los sistemas 3DCT, la unión de LCA estuvo impulsada principalmente por Arg331, mientras que la unión de UDCA implicó una red energética más distribuida que comprende Glu326, Asp394, Arg395, Arg441 y Asp470. En los cuatro sistemas acuosos, la descomposición MM-PBSA respaldó la priorización relativa de los complejos basados en LCA.

Para el complejo triptamina-HTR2A/6A93 embebido en la membrana, se realizó un análisis MM-PBSA en el subsistema proteína-ligando extraído de la trayectoria en la bicapa^46,47. Se observaron contribuciones favorables para el ligando y Asp155 (D3.32), que fue con diferencia el residuo que más contribuyó al estabilización, lo cual es consistente con la interacción de puente salino identificada tanto en los análisis de acoplamiento molecular como en los análisis de distancia a lo largo de la trayectoria. Trp137 mostró la contribución desfavorable más grande a nivel individual entre los residuos circundantes del bolsillo ortostérico (Ser86, Phe87, Phe133, Phe140, Phe141, Val156, Ser159, Thr160, Ile163, Val167, Tyr171), que en conjunto forman la red de contactos aromáticos y polares que reviste el bolsillo de unión. Estos valores representan estimaciones computacionales relativas para la priorización estructural y no corresponden a afinidades de unión experimentales.

Figura 1: Flujo de trabajo computacional para la priorización de genes hospedadores asociados a metabolitos en el IBS-C. Representación esquemática del flujo de trabajo de ocho etapas que integra la selección de metabolitos, predicción de dianas, expresión diferencial transcriptómica, análisis de superposición, enriquecimiento de redes y vías, acoplamiento molecular, simulación de dinámica molecular y análisis de energía libre de unión MM-PBSA. Haga clic aquí para ver una versión más grande de esta figura.

Figura 2: Análisis de expresión diferencial y superposición de metabolitos y dianas en la mucosa de IBS-C. (A) Gráfico de volcan de la expresión diferencial a nivel génico en GSE36701. Puntos azules, genes significativamente regulados a la baja; puntos rojos, genes significativamente regulados al alza; puntos grises, genes no significativos. Se etiquetan los genes asociados a metabolitos superpuestos seleccionados. (B) Diagrama de Venn que muestra la superposición entre 330 dianas de metabolitos predichas únicas y los genes regulados a la baja en GSE36701; se compartieron 17 genes. Haga clic aquí para ver una versión más grande de esta figura.

Figura 3: Evaluación estadística de los 17 genes diana de metabolitos predichos frente a GSE36701. (A) Cambio logarítmico en la expresión (log2) por gen para los 17 genes, coloreado según el nivel de significancia. (B) Tasa de expresión diferencial de los genes de fondo frente a los genes diana predichos, con la prueba exacta de Fisher. (C) Una tabla de contingencia de dos por dos se utiliza para la prueba exacta de Fisher. Los 17 genes diana mostraron una regulación descendente significativa; el solapamiento se interpreta como un patrón direccional descriptivo más que como un enriquecimiento estadístico. Haga clic aquí para ver una versión ampliada de esta figura.

Figura 4: Construcción de la red compuesta de interacciones proteína-proteína y redes de enriquecimiento de vías de interacciones proteína-proteína de genes asociados con metabolitos superpuestos. (A) Red 1: meta-red combinada de todos los 17 genes. (B) Red 2: red específica de propionato de ocho genes (CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN, TBXA2R). (C) Red 3: red de triptamina/serotonina de cuatro genes (HTR1B, HTR2A, HTR2B, HTR6). Las redes se generaron para Homo sapiens como mínimo, utilizando una confianza en la construcción de redes de interacciones proteína-proteína y enriquecimiento de vías ≥ 0,700. Las aristas representan asociaciones funcionales respaldadas por anotaciones Haga clic aquí para ver una versión más grande de esta figura.

Figura 5: Representación estructural tridimensional y bidimensional del ácido litocólico en complejo con VDR (ID PDB: 1DB1). (A) Representación tridimensional de superficie y en forma de cinta, con el ácido litocólico mostrado como esferas. (B) Mapa de interacciones bidimensional que muestra el enlace de hidrógeno con Ser278 y los contactos hidrofóbicos y de van der Waals circundantes. Haga clic aquí para ver una versión más grande de esta figura.

Figura 6. Representación estructural tridimensional y bidimensional del ácido litocólico en complejo con NR1H4/FXR (ID PDB: 3DCT). (A) Representación tridimensional de superficie y en forma de cinta. (B) Mapa bidimensional de interacciones que muestra enlaces de hidrógeno con His294 e Ile335, una interacción π-Sigma y contactos circundantes. Haga clic aquí para ver una versión más grande de esta figura.

Figura 7: Representación estructural tridimensional y bidimensional del ácido ursodesoxicólico en complejo con NR1H4/FXR (ID PDB: 3DCT). (A) Representación tridimensional de superficie y en forma de cinta. (B) Mapa de interacciones bidimensional que muestra enlaces de hidrógeno con His447 y Gly322, una interacción π-anión con Val325, un enlace carbono-hidrógeno con Trp469, y contactos donador-donador desfavorables con Arg395 y Gln396. Haga clic aquí para ver una versión más grande de esta figura.

Figura 8: Representación estructural tridimensional y bidimensional de la urolitina A en complejo con CASP3 (ID PDB: 2DKO). (A) Representación tridimensional de superficie y en forma de cinta. (B) Mapa bidimensional de interacciones que muestra enlaces de hidrógeno con Gln161, Ser120 y Arg207, interacciones π-catión con Arg207, un enlace de hidrógeno donador π con Cys163 y contactos circundantes. Haga clic aquí para ver una versión más grande de esta figura.

Figura 9: Representación estructural tridimensional y bidimensional de la triptamina en complejo con HTR2A (ID PDB: 6A93). (A) Representación tridimensional de superficie y tipo cinta generada mediante un programa de visualización molecular tridimensional. (B) Un mapa de interacciones bidimensional generado con una herramienta de visualización molecular y diagrama de interacciones bidimensional, que muestra el puente salino de Asp155 y otras interacciones en el sitio de unión. Haga clic aquí para ver una versión más grande de esta figura.

Figura 10: Validación del protocolo de acoplamiento molecular. (A,B) Reacoplamiento de ligandos co-cristalizados en VDR/1DB1 (RMSD 0,87 Å) y FXR/3DCT (RMSD 1,79 Å); las conformaciones cristalográficas y reacopladas se superponen, ambas por debajo del umbral de aceptación de 2,0 Å. (C) Selectividad en el acoplamiento cruzado: puntuaciones Vina de ácido litocólico y triptamina en su proteína diana frente a proteínas no diana. Haga clic aquí para ver una versión más grande de esta figura.

Figura 11. Análisis de la trayectoria de dinámica molecular del complejo LCA-VDR/1DB1 durante 200 ns. (A) Perfil de RMSD. (B) Perfil de RMSF. (C) Recuento de enlaces de hidrógeno. (D) Perfil de radio de giro. (E) Perfil de SASA. Haga clic aquí para ver una versión más grande de esta figura.

Figura 12: Análisis de la trayectoria de dinámica molecular del complejo urolitina A-CASP3/2DKO durante 200 ns. (A) Perfil de RMSD que muestra fluctuaciones conformacionales amplias y un evento transitorio de alta desviación cerca de los 165 ns. (B) Perfil de RMSF que muestra una flexibilidad pronunciada a nivel de residuos cerca del residuo 175. (C) Recuento de enlaces de hidrógeno. (D) Perfil del radio de giro. (E) Perfil de SASA. Haga clic aquí para ver una versión más grande de esta figura.

Figura 13: Análisis de la trayectoria de dinámica molecular de los sistemas de ácidos biliares NR1H4/FXR (3DCT) durante 200 ns. (A) Perfil de RMSD del esqueleto principal para el complejo 3DCT. (B) Perfil de RMSF del esqueleto principal. (C) Perfil del radio de giro para 3DCT-LCA. (D) Perfil del radio de giro para 3DCT-UDCA. (E) Perfil de SASA para 3DCT-LCA. (F) Perfil de SASA para 3DCT-UDCA. Haga clic aquí para ver una versión más grande de esta figura.

Figura 14: El complejo triptamina-HTR2A incrustado en una bicapa lipídica POPC explícita. El receptor se muestra como un modelo esquemático que atraviesa la bicapa, los lípidos POPC como líneas con los grupos fosfato de las cabezas polares resaltados, y la triptamina dentro del bolsillo ortostérico. Se muestra agua por encima y por debajo de la membrana. Haga clic aquí para ver una versión más grande de esta figura.

Figura 15: Análisis de la trayectoria de dinámica molecular del complejo triptamina-HTR2A/6A93 durante 200 ns en una bicapa lipídica POPC explícita. (A) Perfil de RMSD del esqueleto principal. (B) Perfil de RMSF por residuo. (C) Perfil del radio de giro. (D) Perfil de SASA. (E) Número de enlaces de hidrógeno entre proteína y ligando. Haga clic aquí para ver una versión más grande de esta figura.

Figura 16: Persistencia de la interacción iónica entre la triptamina y Asp155 (D3.32) a lo largo de la trayectoria de 200 ns en la membrana. Se representa gráficamente la distancia mínima entre el nitrógeno amónico de la triptamina y los átomos de oxígeno carboxilato de Asp155 frente al tiempo; la línea discontinua indica el umbral de contacto del puente salino de 0,4 nm. El contacto se mantuvo durante el 99,9 % de la simulación. Haga clic aquí para ver una versión más grande de esta figura.

Símbolo del gen	Metabolito(s) de origen	Categoría funcional	log2FC	FDR (valor P ajustado)	Nivel de significancia
GCG	Propionato	Proteína relacionada con hormona peptídica	−1.342	1.97e−7	FDR <0.001 & |logFC > 1
HDAC3	Butirato	Enzima	−1.234	2.44e−6	FDR <0.001 & |logFC| > 1
CASP3	Urolitina A	Enzima	−1.198	6.66e−7	FDR <0.001 & |logFC| > 1
GPR68	Propionato	Receptor de membrana	−1.137	4.35e−6	FDR <0.001 & |logFC| > 1
GNAQ	Propionato	Proteína de señalización intracelular	−1.122	1.05e−6	FDR <0.001 & |logFC| > 1
GPHN	Propionato	Otra proteína intracelular	−1.109	1.13e−6	FDR <0.001 & |logFC| > 1
TBXA2R	Propionato	Receptor de membrana	−1.104	4.04e−7	FDR <0.001 & |logFC| > 1
HTR6	Triptamina	Receptor de membrana	−0.967	2.17e−5	FDR <0.001
VDR	Ácido litocólico	Receptor nuclear	−0.942	5.73e−7	FDR <0.001
HTR2A	Triptamina	Receptor de membrana	−0.937	4.99e−6	FDR <0.001
FFAR2	Propionato	Receptor de membrana	−0.889	1.44e−4	FDR <0.001
NR1H4	Ácido litocólico / Ácido ursodesoxicólico	Receptor nuclear	−0.861	3.68e−6	FDR <0.001
HTR2B	Triptamina	Receptor de membrana	−0.702	1.29e−4	FDR <0.001
MLN	Propionato	Proteína relacionada con hormona peptídica	−0.605	7.39e−5	FDR <0.001
KYAT1	Ácido indol-3-láctico / Ácido indol-3-propiónico	Enzima	−0.530	3.61e−4	FDR <0.001
CASR	Propionato	Receptor de membrana	−0.483	4.05e−4	FDR <0.001
HTR1B	Triptamina	Receptor de membrana	−0.455	3.18e−2	FDR <0.05

Tabla 1: Genes diana asociados a metabolitos predichos que se solapan con genes diferencialmente expresados en el conjunto de datos del tejido rectal de IBS-C. Todos los genes solapados enumerados estaban regulados a la baja. La Tabla 1 se presenta por separado como una hoja de cálculo e incluye, para cada gen diana, el(los) metabolito(s) de origen, la categoría funcional, la fuente de predicción del gen diana (una predicción de interacción químico-proteína, un programa de acoplamiento molecular o ambos), una puntuación combinada de interacción de la predicción de interacción químico-proteína, y cuando está disponible, la probabilidad del programa de acoplamiento molecular, el nivel de predicción, el cambio en el log2 de la razón, y el FDR con el nivel de significancia de la expresión. Fuente: los valores de expresión génica se obtuvieron de la tabla de expresión diferencial GSE36701 colapsada por gen (sonda con menor FDR por gen). Los valores de la fuente de predicción del gen diana y de confianza se compiló a partir de la salida de una predicción de interacción químico-proteína y de un programa de acoplamiento molecular, utilizando umbrales de una puntuación combinada de interacción químico-proteína ≥ 0,700 y una probabilidad de acoplamiento molecular ≥ 0,70. Las puntuaciones de predicción de interacción químico-proteína son puntuaciones combinadas en una escala de 0 a 1; STP denota la probabilidad del programa de acoplamiento molecular. Nivel 1 = apoyo estricto de predicción de interacción químico-proteína; Nivel 1+ = apoyo estricto de predicción de interacción químico-proteína corroborado adicionalmente por un programa de acoplamiento molecular.

Complejo	Proteína (ID PDB)	Ligando	Puntuación Vina (kcal/mol)	Tamaño de la cavidad (A^3)	Centro de la cuadrícula X,Y,Z (A)	Caja de búsqueda (A)
LCA-VDR	VDR (1DB1)	Ácido litocólico	−10.0	2055	10, 19, 33	25 x 25 x 25
LCA-NR1H4/FXR	NR1H4/FXR (3DCT)	Ácido litocólico	−9.9	3395	137, 31, 78	25 x 25 x 25
UDCA-NR1H4/FXR	NR1H4/FXR (3DCT)	Ácido ursodesoxicólico	−9.4	3395	137, 31, 78	25 x 25 x 25
Urolitina A-CASP3	CASP3 (2DKO)	Urolitina A	−7.1	233	37, 34, 32	25 x 25 x 25
Triptamina-HTR2A	HTR2A (6A93)	Triptamina	−7.1	3238	12, −1, 61	25 x 25 x 25

Tabla 2: Resultados de acoplamiento molecular: mejores puntuaciones de acoplamiento molecular de ligandos metabolitos a proteínas diana y parámetros de cavidad para los cinco complejos proteína-ligando prioritarios. El tamaño de la cavidad se reporta en ų. Fuente: Docking_Validation/Results/Docking_Validation_Results.xlsx, hoja 'Original_Docking_Scores'. Acoplamiento molecular de ligandos metabolitos a proteínas diana; exhaustividad = 8, semilla = 42 (fija), num_modos = 9 para todos los complejos; se reporta la conformación mejor clasificada (modo 1).

Tipo de interacción	Residuo(s)	Distancia (A)	Notas
Enlace de hidrógeno convencional	Ser278	4.29	Oxígeno carboxilato de LCA
Contacto hidrófobo / Pi-alquilo	Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233, His397	-
Contacto de Van der Waals	Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275, Phe150	-

Tabla 3: Modos de unión generados para el acoplamiento del ácido litocólico con VDR (ID de PDB: 1DB1). Fuente: herramienta de visualización molecular y diagramas bidimensionales de interacción ligando-residuo, como se informa en los Resultados del manuscrito (Acoplamiento molecular). '-' indica que no se informó individualmente el valor de distancia para ese contacto.

Tipo de interacción	Residuo(s)	Distancia (A)	Notas
Enlace de hidrógeno	His294	-
Enlace de hidrógeno	Ile335	-
Interacción Pi-Sigma	His294	-
Alquilo / Pi-Alquilo (hidrófobo)	Met290, Met328, Ala291, Leu287, Ile352, His447	-
Contacto de Van der Waals	Residuos adicionales del bolsillo (no especificados individualmente en la fuente)	-	Apoya la acomodación del esqueleto esteroideo

Tabla 4: Modos de unión generados para el acoplamiento del ácido litocólico con NR1H4/FXR (ID PDB: 3DCT).

Fuente: visualización molecular y herramienta de diagramas de interacción bidimensionales; diagramas de interacción ligando-residuo 2D, según se informa en los Resultados del manuscrito (Docking molecular). '-' indica que no se informó individualmente el valor de distancia para ese contacto.

Tipo de interacción	Residuo(s)	Distancia (A)	Notas
Enlace de hidrógeno convencional	His447	3.66
Enlace de hidrógeno	Gly322	3.46
Interacción Pi-Anión	Val325	4.96
Enlace carbono-hidrógeno	Trp469	4.51
Contacto donador-donador desfavorable	Arg395	3.89
Contacto donador-donador desfavorable	Gln396	3.40

Tabla 5: Modos de unión generados para el acoplamiento molecular del ácido ursodesoxicólico con NR1H4/FXR (PDB ID: 3DCT). Fuente: herramienta de visualización molecular y diagramas bidimensionales de interacción ligando-residuo, tal como se informa en los Resultados del manuscrito (Acoplamiento molecular). '-' indica que no se informó individualmente el valor de distancia para ese contacto.

Tipo de interacción	Residuo(s)	Distancia (A)	Notas
Enlace de hidrógeno convencional	Gln161	3.78
Enlace de hidrógeno convencional	Gln161	4.19	segundo contacto
Enlace de hidrógeno convencional	Ser120	3.95
Enlace de hidrógeno convencional	Arg207	3.05
Enlace de hidrógeno convencional	Arg207	3.77	segundo contacto
Interacción Pi-Catión	Arg207	-
Enlace de hidrógeno Pi-Donador	Cys163	-
Contacto Pi-Alquilo / de van der Waals	Arg64, Ala162, His121, Ser205, Trp206	-

Tabla 6: Modos de unión generados para el acoplamiento de urolitina A con CASP3 (ID de PDB: 2DKO). Fuente: herramienta de visualización molecular y diagramas bidimensionales de interacción ligando-residuo, según se informa en los Resultados del manuscrito (Acoplamiento molecular). '-' indica que no se informó individualmente el valor de distancia para ese contacto.

Tipo de interacción	Residuo(s)	Distancia (A)	Notas
Puente salino electrostático	Asp155 (D3.32)	-	amina protonada de la triptamina
Enlace de hidrógeno	Thr160	-
Enlace de hidrógeno	Ser159	-
Contacto aromático	Phe340, Trp336	-
Interacción Pi-Alquilo	Val156, Ile163	-
Contacto de Van der Waals	Tyr370, Phe339, Ser242, Phe243, Phe332, Leu123	-

Tabla 7: Modos de unión generados para el acoplamiento molecular de la triptamina con HTR2A (PDB ID: 6A93). Fuente: herramienta de visualización molecular y diagramas bidimensionales de interacción ligando-residuo, según se informa en los Resultados del manuscrito (Acoplamiento molecular). '-' indica que no se informó individualmente el valor de distancia para ese contacto.

(A) Validación por redocking (controles positivos)
ID PDB	Proteína	Ligando de co-cristalización	Puntuación de Vina (kcal/mol)	RMSD (Å)	Umbral (Å)	Resultado
1DB1	VDR	VDX (análogo de vitamina D)	−13.0	0.87	2.0	PASA
3DCT	FXR	WAY-362450 (064)	−11.9	1.79	2.0	PASA
(B) Validación por cross-docking (controles negativos)
Ligando	Objetivo cognado (PDB)	Puntuación cognada (kcal/mol)	Objetivo no cognado (PDB)	Puntuación no cognada (kcal/mol)	Delta (kcal/mol)	Selectividad
Ácido litocólico	VDR (1DB1)	−10.0	CASP3 (2DKO)	−8.3	1.7	Confirmada
Triptamina	HTR2A (6A93)	−7.1	VDR (1DB1)	−6.4	0.7	Modesta (dentro de la incertidumbre de Vina +/−0.5–1.0)

Tabla 8: Resultados de la validación del protocolo de acoplamiento molecular: valores de RMSD en redocking (controles positivos) y puntuaciones de cross-docking (controles negativos). Fuente: Docking_Validation/Results/Docking_Validation_Results.xlsx y Docking_Validation/Logs/*.log (un acoplamiento molecular de ligandos metabolitos a proteínas diana, exhaustividad = 8, semilla = 42, caja de 25 Å × 25 Å × 25 Å). El RMSD se calculó mediante coincidencia de átomos pesados por nombre de átomo (sin superposición).

Complejo	RMSD (nm), media + / – DE (rango)	Rg (nm), media + / – DE (rango)	SASA (nm^2), media + / – DE (rango)	Puentes de hidrógeno, media + / – DE (rango)	RMSF (nm), media (máx)
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 (máx 0.600 en el residuo 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 (máx 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 (máx 0.298)
Urolitina 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 (máx 2.532 en el residuo 175)
Triptamina-HTR2A/6A93 (membrana)	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 (máx 0.319)

Tabla 9: Resumen del comportamiento de la simulación de dinámica molecular de 200 ns para los cinco complejos proteína-ligando prioritarios, incluido el sistema de triptamina-HTR2A embebido en la membrana. Fuente: archivos de utilidades de análisis de trayectorias de dinámica molecular (.xvg) — gmx rms, gmx gyrate, gmx sasa, gmx hbond, gmx rmsf — calculados durante los últimos 150 ns (50–200 ns) de cada corrida de producción de 200 ns, según el paso 8.8 del protocolo. RMSD/Rg ajustados al esqueleto principal; radio de sonda SASA de 0,14 nm; corte para donador-aceptor de enlaces de hidrógeno de 0,35 nm / 30 °. LCA-3DCT y UDCA-3DCT comparten una trayectoria del esqueleto principal de la proteína (RMSD, RMSF) con Rg/SASA/enlaces de hidrógeno específicos del ligando.

Triptamina-HTR2A/6A93 (membrana) — descomposición cuantitativa por residuo
Residuo	Contribución total de ddG (kcal/mol), media + / − DE	Dirección
Asp155 (D3.32)	−89.94 + / − 6.81	Estabilizante (dominante)
Triptamina (ligando)	−13.01 + / − 6.22	Estabilizante
Tyr171	13.62 + / − 4.54	Destabilizante
Val167	23.32 + / − 3.96	Destabilizante
Val156	20.03 + / − 3.81	Destabilizante
Thr160	4.86 + / − 3.64	Destabilizante
Ser159	24.16 + / − 3.48	Destabilizante
Ser86	24.48 + / − 3.65	Destabilizante
Phe87	35.18 + / − 4.04	Destabilizante
Phe133	32.80 + / −3.70	Destabilizante
Phe140	30.63 + / − 3.84	Destabilizante
Phe141	35.25 + / − 3.55	Destabilizante
Ile163	27.64 + / − 3.71	Destabilizante
Trp137	53.77 + / − 4.32	Destabilizante (mayor desfavorable)
Otros cuatro complejos — residuos identificados en la descomposición por residuo (cualitativo)
Complejo	Residuo	Dirección
LCA-VDR/1DB1	Ligando (LCA)	Favorable
LCA-VDR/1DB1	Gln317	Favorable
LCA-VDR/1DB1	Trp286	Desfavorable
LCA-NR1H4/FXR/3DCT	Arg331	Favorable (dominante)
UDCA-NR1H4/FXR/3DCT	Glu326	Red mixta/distribuida
UDCA-NR1H4/FXR/3DCT	Asp394	Red mixta/distribuida
UDCA-NR1H4/FXR/3DCT	Arg395	Red mixta/distribuida
UDCA-NR1H4/FXR/3DCT	Arg441	Red mixta/distribuida
UDCA-NR1H4/FXR/3DCT	Asp470	Red mixta/distribuida
Urolitina A-CASP3/2DKO	Arg64	Fuertemente favorable (polar/electrostática)
Urolitina A-CASP3/2DKO	Arg207	Fuertemente favorable (polar/electrostática)

Tabla 10: Descomposición por residuo del MM-PBSA RESUMEN CORTO: residuos estabilizantes e inestabilizantes (contribución absoluta ≥ 0,5 kcal mol⁻¹) para cada uno de los cinco complejos proteína-ligando prioritarios, incluido el sistema de triptamina-HTR2A embebido en la membrana. Fuente: Membrane Simulation/03_MMPBSA/results/FINAL_DECOMP_MMPBSA.dat (herramienta de cálculo de energía de unión mediante mecánica molecular/disolvente continuo Generalized Born (GB), descomposición por residuo, 'Complex: Total Energy Decomposition'). Los números de residuo se convirtieron de la numeración interna del sistema construido con CHARMM-GUI (desplazamiento +68) a la numeración original del PDB 6A93 utilizada en otras partes de este manuscrito.

Fuente: datos de simulación previa de dinámica molecular/1DB1,2KD0, LCA & UDCA_3DCT}/mmpbsa_*/Descomposición_NORMAL_GB_Complejo_TDC*.svg y resultados del manuscrito (energía libre de unión MM-PBSA y descomposición por residuo). Estos cuatro complejos no tienen salida numérica por residuo en formato .dat/.csv en el directorio del proyecto (solo gráficos SVG renderizados con texto en trazados vectoriales que no son extraíbles mediante máquina); únicamente se informa la identidad de los residuos y la dirección favorable o desfavorable, según se indica en el texto del manuscrito. Las contribuciones exactas en kcal/mol para estos cuatro complejos no están disponibles en el repositorio fuente.

Discusión

Este estudio computacional exploratorio demuestra un flujo de trabajo integrado y reproducible para priorizar genes del huésped y complejos proteína-ligando asociados a metabolitos microbianos, aplicado aquí a un conjunto de datos transcriptómicos públicos de la mucosa rectal en IBS-C. Mediante este flujo de trabajo, un subconjunto de genes predichos asociados a metabolitos microbianos coincidió con genes consistentemente regulados a la baja en el conjunto de datos, agrupándose dentro de vías asociadas a receptores acoplados a proteínas G, sistema serotoninérgico, señalización por calcio, ligandos-receptores neuroactivos y receptores nucleares, sistemas que se implican cada vez más en la comunicación entre microbiota y huésped^48,49. Estos hallazgos deben interpretarse estrictamente como generadores de hipótesis: el análisis no mide concentraciones de metabolitos microbianos, abundancia de proteínas receptoras, unión del ligando, activación del receptor, señalización descendente, motilidad, secreción, respuestas al dolor ni resultados clínicos. La conclusión más sólida que puede extraerse es que los genes y vías identificados son candidatos para validación experimental, más que mecanismos de enfermedad confirmados.

La importancia fundamental de este protocolo, en comparación con trabajos previos que examinaron pares individuales de metabolitos y receptores de forma aislada, radica en su integración de la predicción de objetivos, la transcriptómica pública, el análisis de redes, el acoplamiento molecular con controles de validación, la dinámica molecular y el MM-PBSA en una única tubería secuencial de priorización. Cada etapa reduce y contextualiza el conjunto de candidatos producido por la etapa anterior, y es el filtrado secuencial lo que hace que la lista final de candidatos sea manejable experimentalmente. El módulo de RCPG centrado en GNAQ y el módulo asociado al receptor de serotonina identificados aquí son biológicamente plausibles, dada la función de la señalización Gq en la activación de la fosfolipasa C, la producción de inositol 1,4,5-trisfosfato, la movilización de calcio, la secreción y la función enteroendocrina, así como los roles establecidos de la señalización derivada de ácidos grasos de cadena corta y del triptófano en la homeostasis de la mucosa y de la señalización serotoninérgica en la motilidad gastrointestinal, la secreción, la sensibilidad visceral y la comunicación intestino-cerebro^10,18,50,51.

Una característica metodológica clave de este estudio es el tratamiento del receptor HTR2A embebido en la membrana. Dado que una simulación en fase soluble no puede reproducir el entorno lipídico que rige el comportamiento conformacional de un receptor acoplado a proteína G, el complejo triptamina-HTR2A se simuló en una bicapa POPC explícita. En este entorno membranal, el receptor permaneció estructuralmente estable a lo largo de toda la trayectoria de 200 ns, y el puente de sal entre el grupo amonio de la triptamina y Asp155 (D3.32) se mantuvo prácticamente durante toda la simulación. El hecho de que tres líneas independientes de evidencia —la posición de acoplamiento, la distancia de contacto persistente a lo largo de la trayectoria y la contribución dominante por residuo en el análisis MM-PBSA— converjan en la misma interacción conservada con D3.32, aporta coherencia interna al modo de unión de la triptamina predicho, que recupera la geometría de unión canónica de los ligandos aminérgicos en los receptores de serotonina.

Existen ciertas cuestiones metodológicas que deben considerarse al reproducir este flujo de trabajo. Los errores en la estructura canónica o en los compuestos con interferencia panensayo (compuestos marcados como PAINS) se propagan mediante la predicción de dianas y el acoplamiento molecular, lo que exige una selección precisa de metabolitos y una curación quimioinformática rigurosa. Minimice conjuntos de dianas impulsados por ruido mediante la aplicación consistente de criterios de confianza (predicción de dianas mediante interacción químico-proteica: ≥0,700; programa de acoplamiento molecular: ≥0,70; construcción de redes de interacción proteína-proteína y enriquecimiento de vías: ≥0,700). Las dianas predichas deben agruparse por categoría funcional para evitar la incorrecta caracterización de todos los genes asociados a metabolitos como receptores. Un preprocesamiento preciso de las estructuras PDB, la minimización de la energía del ligando y la colocación de la cuadrícula alrededor de los residuos de unión conocidos son aspectos esenciales del acoplamiento molecular, y los controles de reacoplamiento y acoplamiento cruzado introducidos aquí proporcionan una medida objetiva de la corrección del método de acoplamiento. El margen de reproducibilidad en dinámica molecular está definido por una combinación de parametrización del campo de fuerzas, solvatación adecuada o construcción de membranas, equilibrio por etapas y un muestreo suficiente durante la fase de producción.

Los ajustes y pasos habituales de solución de problemas incluyen relajar los umbrales si la predicción del objetivo no devuelve resultados, verificar la consistencia direccional a nivel de sonda para genes con múltiples sondas e interpretar los nodos aislados de la construcción de redes de interacción proteína-proteína y enriquecimiento de vías como dependientes del umbral en lugar de biológicamente irrelevantes. Para los receptores de membrana, debe utilizarse una simulación explícita de bicapa lipídica en lugar de una simulación acuosa, como se ejemplifica con el método HTR2A descrito aquí. Cuando sea necesario un desglose energético por residuo, el cálculo debe realizarse con un motor capaz de realizar descomposiciones, y la numeración de residuos informada debe conciliarse con la numeración nativa del receptor para evitar ambigüedades. Sugerimos que los resultados de enriquecimiento de vías deben considerarse mejor como un contexto organizativo para la lista de candidatos, en lugar de una validación a nivel de vía. Mecánicamente, se producirá enriquecimiento de términos relacionados con receptores acoplados a proteínas G (GPCR), señalización serotoninérgica o señalización de calcio siempre que la lista de genes contenga múltiples genes de receptores de serotonina, independientemente de la coregulación a nivel de proteína. Los valores de RMSD, Rg y RMSF para el sistema HTR2A insertado en la membrana deben interpretarse considerando la bicapa lipídica: una disminución de Rg en la trayectoria posterior puede reflejar una adaptación conformacional inducida por la bicapa del conjunto de dominios transmembrana, en lugar de una desnaturalización global, y los enlaces de hidrógeno persistentes entre ligando y proteína deben interpretarse junto con la estabilidad general del RMSD.

Las limitaciones de este estudio son considerables y restringen la interpretación. La investigación se basó en un conjunto de datos público relativamente pequeño, y una búsqueda en los principales repositorios transcriptómicos públicos (herramienta web de análisis de expresión génica diferencial y ArrayExpress) no identificó un conjunto de datos transcriptómicos independiente del tejido rectal de pacientes con IBS-C de diseño y plataforma comparables que pudiera servir como cohorte de replicación en el momento del análisis. La ausencia de replicación transcriptómica independiente constituye una limitación importante, y ninguna afirmación de este manuscrito debe interpretarse como una validación externa de los hallazgos derivados de un único conjunto de datos. El conjunto de datos muestra una expresión diferencial casi universal (aproximadamente el 94,5 % de los genes son significativos, la gran mayoría de los cuales están regulados a la baja), una característica que hace que las estadísticas convencionales de enriquecimiento sean poco informativas y que impida conclusiones sobre la especificidad de la regulación a la baja de los genes diana en relación con el fondo genómico; por lo tanto, la superposición se informa como un patrón direccional descriptivo en lugar de un enriquecimiento estadístico. La transcriptómica masiva del tejido mucoso no puede distinguir la regulación génica real de los cambios en la composición celular. La expresión de ARNm no determina la abundancia de proteínas ni la respuesta funcional. Las bases de datos de predicción de dianas presentan sesgos de anotación, y los resultados del acoplamiento molecular, la dinámica molecular (MD) y el MM-PBSA dependen de la elección del campo de fuerza, la parametrización del ligando, la posición inicial, el tiempo de simulación y la suficiencia del muestreo. Los identificadores exactos de versiones menores o compilaciones de algunos componentes de servidores web y paquetes, incluido el servicio de parametrización de ligandos compatible con CHARMM, CHARMM-GUI y el entorno de computación estadística

Las versiones de los paquetes generados y las subversiones de la herramienta de cálculo de energía de unión mediante mecánica molecular/disolvente continuo no pudieron recuperarse completamente a partir del registro archivado del proyecto y deben indicarse como disponibles en la Tabla de Materiales aparte. Los valores de MM-PBSA son estimaciones relativas, no incluyen un término explícito de entropía conformacional y no deben interpretarse como afinidades experimentales. El estudio carece de datos metabolómicos y no puede determinar si la disponibilidad del ligando se altera en el SII-C, ni si las alteraciones observadas en la expresión son causas, consecuencias, respuestas compensatorias o correlaciones irrelevantes.

Las aplicaciones futuras de este método deberían incluir replicación transcriptómica independiente, validación mediante reacción en cadena de la polimerasa cuantitativa (qPCR) y a nivel de proteínas, localización en tipos celulares mediante transcriptómica unicelular o espacial, perfilado metabolómico de las clases de metabolitos relevantes, y ensayos funcionales de respuesta a ligandos en colonoides derivados de pacientes, explantes de mucosa u otros modelos comparables. Las comparaciones con cohortes de síndrome del intestino irritable con predominio de diarrea, síndrome del intestino irritable mixto, enfermedad inflamatoria intestinal y constipación no asociada al síndrome del intestino irritable^1,2 ayudarían a establecer la especificidad de la enfermedad. Para el componente estructural, la repetición de las trayectorias de dinámica molecular, el análisis de sensibilidad con distintas conformaciones iniciales y la documentación completa del depósito de los archivos de topología, trayectorias y entradas y salidas de MM-PBSA reforzarían aún más la reproducibilidad. Siguen siendo necesarios ensayos experimentales de respuesta a ligandos para determinar si los complejos prioritarios son funcionalmente relevantes; los resultados actuales no respaldan afirmaciones clínicas ni terapéuticas.

Materiales

Lista de materiales utilizados en este artículo

Nombre	Empresa	Número de catálogo	Comentarios
AutoDock Vina	Scripps Research / código abierto	v1.2.7; https://vina.scripps.edu/	Acoplamiento molecular de ligandos metabolitos a proteínas diana.
CGenFF/ParamChem	SilcsBio / Universidad de Maryland	v4.6; https://cgenff.com/	Parametrización del campo de fuerza del ligando para dinámica molecular.
Campo de fuerza CHARMM36m	Desarrolladores de CHARMM / código abierto	CHARMM36m; https://www.charmm.org/charmm/resources/charmm-force-fields/	Campo de fuerza proteico utilizado en simulaciones de dinámica molecular.
Constructor de membranas CHARMM-GUI	CHARMM-GUI / Universidad Lehigh	Servidor web; versión exacta no recuperable; https://www.charmmgui.org/?doc=input/membrane	Configuración de construcción y equilibrio del sistema de membrana POPC explícita.
Discovery Studio Visualizer	BIOVIA (Dassault Systèmes)	2021; https://discover.3ds.com/discovery-studio-visualizer-download	Análisis bidimensional de interacciones ligando-residuo.
GEO2R	NCBI Gene Expression Omnibus	Herramienta web; acceso entre enero y mayo de 2026; https://www.ncbi.nlm.nih.gov/geo/geo2r/	Análisis de expresión diferencial de GSE36701.
GeneCards	Instituto Weizmann de Ciencia	Base de datos web; acceso entre enero y mayo de 2026; https://www.genecards.org/	Verificación de símbolos génicos e información genética durante la estandarización de dianas.
gmx_MMPBSA	Código abierto (Valdés-Tresanco et al.)	1.5.x; https://valdes-tresanco-ms.github.io/gmx_MMPBSA/	Estimación de energía libre de unión mediante MM-PBSA y descomposición por residuo.
GROMACS	Equipo de desarrollo de GROMACS / código abierto	2024.2; https://www.gromacs.org/	Motor de simulación de dinámica molecular.
Conjunto de datos transcriptómicos GSE36701	NCBI Gene Expression Omnibus	GSE36701; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE36701	Conjunto de datos públicos de expresión en mucosa rectal de IBS-C.
Open Babel	Código abierto	3.2.0; https://openbabel.org/	Conversión de formatos de archivos químicos, generación tridimensional de ligandos y preparación de ligandos.
Base de datos OPM	Universidad de Michigan	Base de datos web; acceso entre enero y mayo de 2026; https://opm.phar.umich.edu/	Coordenadas de Orientación de Proteínas en Membranas utilizadas para alinear HTR2A.
ParmEd	Desarrolladores de ParmEd / código abierto	4.x; https://parmed.github.io/ParmEd/html/index.html	Reparticionamiento de masa de hidrógeno y procesamiento de topologías para simulaciones moleculares.
PyMOL	Schrödinger / código abierto	2.x; https://www.pymol.org/	Visualización estructural tridimensional y generación de figuras receptor-ligando.
Banco de Datos de Proteínas RCSB	RCSB PDB	Base de datos web; acceso entre enero y mayo de 2026; https://www.rcsb.org/	Fuente de estructuras proteicas experimentales y coordenadas PDB.
STITCH	Consortium STITCH (EMBL)	v5.0; https://stitch.embl.de/	Predicción de dianas en interacciones químico-proteína.
STRING	Consortium STRING / ELIXIR	v12.0; https://version-12-0.string-db.org/	Construcción de redes de interacción proteína-proteína y enriquecimiento de vías.
SwissADME	Instituto SIB de Bioinformática / Universidad de Lausana	Herramienta web; acceso entre enero y mayo de 2026; https://www.swissadme.ch/	Descriptores quimioinformáticos, predicciones farmacocinéticas y evaluación de PAINS.
SwissTargetPrediction	Instituto SIB de Bioinformática / Universidad de Lausana	Herramienta web; acceso entre enero y mayo de 2026; https://www.swisstargetprediction.ch/	Predicción basada en ligandos de dianas proteicas humanas.
Mapeo de identificadores UniProt	Consortium UniProt	Servicio web; acceso entre enero y mayo de 2026; https://www.uniprot.org/id-mapping	Mapeo de identificadores proteicos a símbolos génicos estandarizados aprobados por HGNC.
NVIDIA RTX 3080	Corporación NVIDIA	RTX 3080; ≥8 GB de VRAM; versión de CUDA/controlador no especificada en el manuscrito	Unidad de procesamiento gráfico compatible con CUDA utilizada en simulaciones de dinámica molecular.
GPU compatible con CUDA	Corporación NVIDIA	Versión del kit de herramientas CUDA no especificada en el manuscrito; ≥8 GB de VRAM	GPU compatible con CUDA con ≥8 GB de VRAM; la estación de trabajo también requería ≥32 GB de RAM y una CPU de 6 núcleos.
Ubuntu Linux	Canonical Ltd. / código abierto	22.04 LTS	Sistema operativo Linux de 64 bits.
Python 3.9	Fundación Python Software	3.9	Entorno de programación de propósito general utilizado para scripts de flujo de trabajo y análisis.
Gene Expression Omnibus (GEO)	NCBI / Biblioteca Nacional de Medicina de EE. UU.	Repositorio web público; no se especificó versión de software en el manuscrito	Repositorio público de datos de genómica funcional.
AutoDockTools/MGLTools	Laboratorio de Gráficos Moleculares, Scripps Research	1.5.7	Kit de herramientas para preparación de estructuras moleculares y entradas de acoplamiento.
Herramientas de análisis de GROMACS	Equipo de desarrollo de GROMACS / código abierto	2024.2	Utilidades para análisis de trayectorias de dinámica molecular.
Protocolo de seis pasos CHARMM-GUI	CHARMM-GUI / Universidad Lehigh	Protocolo web; versión exacta no recuperable	Flujo de trabajo web multietapa para preparación y equilibrio de sistemas moleculares.
cgenff_charmm2gmx_py3.py	Script de conversión de código abierto; fuente no especificada en el manuscrito	Versión no especificada en el manuscrito	Script de conversión de topologías de campo de fuerza.
Python	Fundación Python Software	3.9	Entorno de programación de propósito general.
SciPy	Comunidad SciPy / código abierto	Versión no especificada en el manuscrito	Biblioteca de computación científica.
scipy.stats.fisher_exact	Comunidad SciPy / código abierto	Versión de SciPy no especificada en el manuscrito	Implementación de la prueba exacta de Fisher’s.
R	Fundación R para Computación Estadística	4.3.x	Entorno de computación estadística.
Bioconductor	Proyecto Bioconductor / código abierto	3.18	Marco de software de bioinformática.
limma	Proyecto Bioconductor / código abierto	Versión no especificada en el manuscrito	Paquete para análisis de expresión génica diferencial.
Procedimiento de Benjamini–Hochberg	Método estadístico	No aplicable (procedimiento estadístico)	Método de ajuste de tasa de falsos descubrimientos.
NVIDIA RTX 3080	Corporación NVIDIA	RTX 3080; ≥8 GB de VRAM; versión de CUDA/controlador no especificada en el manuscrito	Unidad de procesamiento gráfico con al menos 8 GB de memoria de video.
GPU compatible con CUDA	Corporación NVIDIA	Versión del kit de herramientas CUDA no especificada en el manuscrito; ≥8 GB de VRAM	Unidad de procesamiento gráfico que soporta computación paralela de propósito general.
Ubuntu Linux 22.04 LTS	Canonical Ltd. / código abierto	22.04 LTS	Sistema operativo Linux de 64 bits.
TIP3P	Desarrolladores del campo de fuerza CHARMM / código abierto	TIP3P; no aplica versión de software	Modelo explícito de agua de tres sitios.
MM/PBSA	Desarrolladores de gmx_MMPBSA / código abierto	gmx_MMPBSA 1.5.x	Método de energía de unión por mecánica molecular/Poisson–Boltzmann y área superficial.