Um Fluxo de Trabalho Computacional para a Priorização de Genes do Hospedeiro Associados a Metabólitos Microbianos na Síndrome do Intestino Irritável com Predomínio de Constipação

Resumo

Este protocolo integra a predição de alvos de metabólitos microbianos, transcritômica da mucosa retal, interação proteína-proteína e enriquecimento de vias, acoplamento molecular, simulação de dinâmica molecular e estimativa da energia livre de ligação por mecânica molecular/área de superfície de Poisson-Boltzmann (MM-PBSA) para gerar uma lista curta classificada, com hipóteses geradoras, de genes hospedeiros candidatos associados a metabólitos e complexos proteína-ligante prioritizados estruturalmente para acompanhamento experimental.

Resumo

Não existe um pipeline computacional padronizado para priorizar de forma sistemática genes do hospedeiro associados a metabólitos microbianos e complexos proteína-ligante a partir de bases de dados públicas disponíveis de natureza química, genômica e estrutural. Este artigo descreve um fluxo de trabalho em oito etapas que aceita um conjunto definido pelo usuário de metabólitos derivados da microbiota intestinal e produz uma lista curta classificada de genes candidatos do hospedeiro associados a metabólitos, vias biológicas enriquecidas e complexos proteína-ligante priorizados estruturalmente para acompanhamento experimental. O pipeline integra (i) perfilagem quimioinformática de metabólitos; (ii) predição de alvos candidatos em múltiplas bases de dados, utilizando ferramenta de interação proteína-químico e predição de alvos baseada em ligante, além de um programa de docking molecular; (iii) análise de expressão gênica diferencial de dados transcriptômicos disponíveis publicamente; (iv) sobreposição entre alvos e genes diferencialmente expressos; (v) construção de rede de interação proteína-proteína e enriquecimento de vias; (vi) docking molecular com um programa de docking molecular; (vii) simulação de dinâmica molecular de 200 ns utilizando um motor de dinâmica molecular com um campo de força proteico empregado em simulações de dinâmica molecular; e (viii) estimativa da energia livre de ligação MM-PBSA. Como exemplo prático, nove metabólitos derivados ou modificados pela microbiota intestinal, representando ácidos graxos de cadeia curta, ácidos biliares, metabólitos derivados de triptofano e urolitina A, foram processados utilizando o conjunto de dados transcriptômicos públicos de mucosa retal de IBS-C GSE36701. O fluxo de trabalho classificou 17 genes únicos preditos associados a metabólitos que foram diferencialmente expressos nesse conjunto de dados. As análises de docking, simulação de dinâmica molecular e MM-PBSA priorizaram estruturalmente cinco complexos metabólito-proteína: ácido litocólico-VDR, ácido litocólico-NR1H4/FXR, ácido ursodesoxicólico-NR1H4/FXR, triptamina-HTR2A (simulada em uma bicamada lipídica explícita de 1-Palmitoil-2-oleoil-sn-glicerol-3-fosfocolina (POPC)) e urolitina A-CASP3. O protocolo foi projetado para ser adaptável a outros conjuntos de metabólitos, conjuntos de dados transcriptômicos de doenças e classes de alvos; todas as saídas são predições computacionais geradoras de hipóteses que exigem replicação transcriptômica independente, validação ao nível proteico e ensaios funcionais de resposta ao ligante antes que conclusões causais ou terapêuticas possam ser estabelecidas.

Introdução

A síndrome do intestino irritável com constipação (IBS-C) é uma desordem gastrointestinal funcional prevalente, caracterizada por dor abdominal recorrente, alterações nos hábitos intestinais, distensão abdominal e constipação, com uma prevalência global estimada em aproximadamente 10–15% da população geral^1,2. As terapias farmacológicas atuais, incluindo secretagogos, pró-cinéticos e antiespasmódicos, podem melhorar sintomas individuais em um subconjunto de pacientes; no entanto, a resposta ao tratamento permanece heterogênea e a remissão duradoura é raramente alcançada, refletindo a patobiologia complexa e multifatorial da condição^1,3,4. É, portanto, necessária uma compreensão mecanicista mais completa de como os sinais microbianos intestinais são transduzidos ao nível da mucosa para gerar hipóteses testáveis sobre novos alvos terapêuticos.

A microbiota intestinal contribui para a homeostase gastrointestinal inferior por meio da produção e biotransformação de metabólitos quimicamente diversos, incluindo ácidos graxos de cadeia curta (AGCC), ácidos biliares secundários, compostos derivados de triptofano e metabólitos derivados de polifenóis, como as urolitinas^5,6,7,8. Essas moléculas comunicam-se com células do hospedeiro por meio de um amplo e ainda parcialmente caracterizado conjunto de alvos moleculares, que vai muito além dos receptores de membrana canônicos sensíveis a metabólitos, abrangendo receptores nucleares, enzimas citosólicas, proteínas modificadoras de histonas, precursores de hormônios peptídicos e proteínas intracelulares de sinalização⁹. Alterações na composição da comunidade microbiana intestinal e nos perfis de metabólitos foram documentadas em pacientes com síndrome do intestino irritável (SII), fornecendo uma base biológica para investigar se genes do hospedeiro associados à resposta a metabólitos microbianos estão transcricionalmente perturbados na mucosa retal de pacientes com SII-C¹⁰.

O painel de nove metabólitos foi definido a priori para fornecer um conjunto compacto, quimicamente diverso e biologicamente interpretável de pequenas moléculas derivadas da microbiota intestinal ou modificadas por ela. A seleção baseou-se em cinco critérios: representação das principais classes de metabólitos microbianos envolvidos na sinalização entre hospedeiro e microbiota; exposição plausível ou conhecida na mucosa intestinal distal; disponibilidade de identificadores PubChem inequívocos e estruturas canônicas; tamanho molecular e tratabilidade estrutural para previsão e encaixe de alvos baseados em ligantes; e plausibilidade prévia de atuação em sinalizações epiteliais, neuroimunes, enteroendócrinas, mediadas por receptores nucleares ou relacionadas à motilidade no IBS-C. O painel selecionado incluiu butirato e propionato como ácidos graxos de cadeia curta (SCFAs); ácido quenodesoxicólico, ácido litocólico e ácido ursodesoxicólico como ácidos biliares; triptamina, ácido indol-3-propiônico e ácido indol-3-lático como metabólitos derivados de triptofano; e urolitina A como um metabólito de polifenóis derivado da microbiota intestinal^5,6,7,8,9,10.

A maioria das investigações computacionais e experimentais anteriores examinou pares isolados de metabólito–receptor ou metabólito–enzima separadamente, uma abordagem que não capta a natureza distribuída e convergente da sinalização por metabólitos microbianos ao longo das vias do hospedeiro^9,11. A integração de múltiplas etapas analíticas confere um poder de filtragem mutuamente reforçador que nenhuma etapa individual pode fornecer isoladamente. A predição computacional de alvos contra bancos de dados curados gera um amplo conjunto de proteínas candidatas do hospedeiro para cada metabólito. A intersecção com dados transcriptômicos relevantes para a doença filtra substancialmente esse conjunto, mantendo apenas os candidatos cujos transcritos são alterados no contexto da doença. As análises de enriquecimento de vias e de redes de interação proteína–proteína, então, mapeiam a lista reduzida de candidatos para módulos biológicos conhecidos. O encaixe molecular fornece uma avaliação computacional inicial da complementaridade do sítio de ligação para cada complexo candidato, e uma simulação adicional de dinâmica molecular (MD) de 200 ns com decomposição da energia livre de ligação por MM-PBSA acrescenta uma dimensão termodinâmica e dependente do tempo à priorização estrutural, que não está disponível apenas a partir das pontuações de encaixe. Realizar cada etapa independentemente, sem integração sistemática e filtragem sequencial, resultaria em listas de candidatos excessivamente amplas para serem viáveis experimentalmente e falharia em detectar a arquitetura convergente das vias.

Dentro do contexto deste protocolo completo, o termo “gene associado a metabólitos” (MAG) refere-se a um gene humano cujo produto proteico foi indicado como um alvo molecular putativo de um ou mais metabólitos derivados da microbiota intestinal por pelo menos uma base de dados computacionais preditivas curada, e cujo transcrito é diferencialmente expresso no conjunto de dados transcriptômicos relevante para a doença, utilizado para demonstrar o fluxo de trabalho. Essa definição operacional inclui deliberadamente receptores de membrana e receptores nucleares, enzimas citosólicas, proteínas de sinalização, precursores de hormônios peptídicos e outras proteínas intracelulares. A designação MAG não constitui evidência experimental de que um metabólito se ligue, forme um complexo proteína-ligante, ative um receptor, altere a abundância proteica ou cause doença, mas sim uma indicação hipotética derivada computacionalmente, que exige validação experimental.

Este protocolo descreve o fluxo de trabalho computacional completo em oito etapas (Figura 1) com detalhes operacionais suficientes para permitir a replicação independente, adaptação a outros painéis de metabólitos ou conjuntos de dados de doenças e extensão a outros contextos de interação hospedeiro-microbiota. O escopo do fluxo de trabalho é explicitamente definido como um quadro para geração de hipóteses e priorização estrutural que opera exclusivamente com recursos ômicos e estruturais disponíveis publicamente, não inferindo concentrações alteradas de metabólitos, estados de ativação de receptores, alterações na expressão proteica, atividade de sinalização downstream ou significância clínica a partir dos resultados computacionais isoladamente. Aqui, demonstramos o protocolo como um exemplo prático utilizando nove metabólitos derivados ou modificados pela microbiota intestinal e o conjunto de dados transcriptômicos públicos do tecido retal de pacientes com IBS-C GSE36701, com o objetivo de identificar MAGs e priorizar complexos metabólito-proteína para acompanhamento experimental subsequente.

Protocolo

A análise utilizou apenas dados transcriptômicos anonimizados e de acesso público do GSE36701 e bancos de dados químicos, proteicos e estruturais de acesso público. Os bancos de dados foram acessados entre janeiro e maio de 2026. Qualquer data de acesso posterior foi documentada na Tabela de Materiais separada.

1. Planejamento do estudo, requisitos de hardware e software

1. Defina o fluxo de trabalho antes de iniciar a análise. Utilize oito etapas: seleção de metabólitos, predição de alvos, análise de expressão diferencial, sobreposição de genes diferencialmente expressos (DEG), enriquecimento de interação proteína-proteína (PPI)/via metabólica, encaixe molecular, simulação de dinâmica molecular (MD) e estimativa de MM-PBSA.
2. Registre que o encaixe molecular, a dinâmica molecular (MD) e a estimativa de MM-PBSA são ferramentas apenas para priorização estrutural. Não interprete essas saídas como evidência experimental de ligação, ativação do receptor, alteração na abundância proteica, eficácia terapêutica ou causalidade da doença.
3. Confirme o hardware computacional antes de executar simulações de dinâmica molecular (MD). Utilize um sistema operacional Linux de 64 bits, uma CPU com 6 núcleos ou superior, uma plataforma com aceleração por GPU e ≥8 GB de memória de vídeo (VRAM), ou uma plataforma equivalente com aceleração por GPU e ≥8 GB de VRAM, com pelo menos 8 GB de VRAM, mínimo de 32 GB de RAM e pelo menos 200 GB de armazenamento livre por sistema de MD.
4. Registre o software principal: um mecanismo de dinâmica molecular, um programa de encaixe molecular de ligantes metabólitos a proteínas-alvo, conversão de formatos de arquivos químicos¹², geração tridimensional de ligantes, preparação de ligantes, um kit de ferramentas para preparação de entradas para encaixe molecular, um programa de encaixe molecular de ligantes metabólitos a proteínas-alvo, um ambiente de programação de uso geral, um ambiente de computação estatística com estrutura de software bioinformático e um pacote de análise de expressão gênica diferencial.
5. Registre as ferramentas de análise estrutural: ferramenta web para construção de sistemas de membrana, serviço de parametrização de ligantes compatível com CHARMM, ferramenta de cálculo de energia de ligação por mecânica molecular/solvente contínuo, biblioteca de conversão de topologia e parâmetros moleculares, programa de visualização molecular tridimensional e ferramenta de visualização molecular e diagrama de interações bidimensional 2021 (consulte a Tabela de Materiais para links de download e informações de versão).
6. Registre os identificadores exatos dos campos de força proteicos utilizados nas simulações de dinâmica molecular, CGenFF, CHARMM-GUI, estrutura de software R/bioinformática e versões da ferramenta de cálculo de energia de ligação por mecânica molecular/solvente contínuo na Tabela de Materiais/arquivo de ambiente. Marque os identificadores ausentes como “não recuperável”; não os infira.

2. Seleção de metabólitos e caracterização quimioinformática

1. Defina o painel de metabólitos antes da predição de alvos. Inclua 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ático (CID: 92904) e urolitina A (CID: 5488186).
2. Obtenha o Simplified Molecular Input Line Entry System (SMILES) canônico e o CID do PubChem para cada metabólito. Verifique sinônimos e estruturas duplicadas antes da predição de alvos. Armazene os identificadores finais na planilha mestra de metabólitos.
3. Envie as sequências SMILES canônicas para a ferramenta web de predição de propriedades físico-químicas e ADME¹³ (veja Tabela de Materiais). Registre o peso molecular, área topológica da superfície polar (TPSA), logP consenso, doadores de ligação de hidrogênio, aceptores de ligação de hidrogênio, ligações rotatóveis, absorção gastrointestinal predita, predição de glicoproteína P e alertas de Lipinski, Veber, Ghose, Egan, Muegge e PAINS.
4. Mantenha os metabólitos com reconhecimento estrutural bem-sucedido, peso molecular ≤500 Da e sem alertas de PAINS. Registre qualquer critério não atendido e a decisão de manter ou excluir o metabólito.
5. Atribua estados de ionização antes da predição de alvos e do encaixe molecular (docking). Utilize carboxilatos desprotonados para butirato e propionato, formas ácidas carboxílicas neutras para os ácidos biliares, forma protonada de amônio para triptamina e formas neutras para os demais metabólitos.

3. Previsão de alvos humanos candidatos

1. Abra uma previsão de alvo de interação químico-proteína¹⁴ (consulte a Tabela de Materiais). Insira o nome de cada metabólito ou o CID do PubChem, selecione Homo sapiens (ID taxonômico: 9606) e defina a pontuação mínima combinada de interação como ≥0,700.
2. Priorize canais de evidência experimentais e de bases de dados curadas em uma previsão de alvo de interação químico-proteína. Baixe a tabela completa de associações proteicas para cada metabólito.
3. Abra um programa de encaixe molecular¹⁵ (consulte a Tabela de Materiais). Envie cada string SMILES canônica com Homo sapiens selecionado e mantenha os alvos com probabilidade ≥0,70.
4. Combine as saídas da previsão de alvo de interação químico-proteína e do programa de encaixe molecular como um conjunto união para cada metabólito. Mantenha qualquer alvo que atenda a qualquer um dos limites da base de dados e remova entradas duplicadas exatas de símbolos gênicos.
5. Padronize as entradas proteicas para símbolos gênicos aprovados pelo Comitê de Nomenclatura Genética HUGO (HGNC) usando o mapeamento de identificadores proteicos para símbolos gênicos padronizados e aprovados pelo HGNC ou a base de dados integrada de informações genéticas humanas (consulte a Tabela de Materiais). Resolva sinônimos, símbolos obsoletos e anotações de isoformas para um único símbolo gênico por proteína.
6. Classifique cada alvo como receptor de membrana, receptor nuclear, enzima, proteína de sinalização intracelular, hormônio peptídico, proteína relacionada a hormônio ou outra proteína intracelular. Registre a classe na tabela de alvos.

4. Conjunto de dados transcriptômico e análise de expressão gênica diferencial

1. Acesse o GSE36701 por meio da ferramenta baseada na web do NCBI para análise de expressão gênica diferencial^16,17 (consulte a Tabela de Materiais). Registre que o conjunto de dados inclui dados de expressão de biópsias da mucosa retal de grupos com síndrome do intestino irritável com predomínio de constipação (IBS-C), síndrome do intestino irritável com predomínio de diarreia (IBS-D), síndrome do intestino irritável pós-infecciosa e voluntários saudáveis¹⁸.
2. Busque no GEO e no ArrayExpress um coorte independente para validação. Utilize combinações dos termos síndrome do intestino irritável com predomínio de constipação (IBS-C), síndrome do intestino irritável com predomínio de constipação, mucosa retal, mucosa cólica, biópsia, transcriptoma, microarranjo e RNA-seq. Registre os repositórios, os termos de busca, a data da busca e se um conjunto de dados de validação comparável foi identificado.
3. Inicie a ferramenta baseada na web para análise de expressão gênica diferencial a partir do registro do GSE36701 (consulte a Tabela de Materiais). Atribua as 18 amostras de IBS-C ao grupo IBS-C, atribua os 40 voluntários saudáveis ao grupo controle e deixe as amostras de IBS-D e de síndrome do intestino irritável pós-infecciosa não atribuídas.
4. Execute a análise de expressão diferencial utilizando o framework do pacote de análise de expressão gênica diferencial com correção da taxa de falsa descoberta (FDR) de Benjamini-Hochberg¹⁹. Baixe a tabela completa de resultados com ID da sonda, símbolo do gene, título do gene, logFC, AveExpr, estatística t moderada, valor P bruto e valor P ajustado.
5. Consolide as sondas em entradas ao nível gênico. Remova as sondas sem símbolos gênicos; mantenha a sonda com o menor FDR para símbolos duplicados; e utilize o maior valor absoluto de logFC como critério de desempate.

5. Análise de sobreposição alvo-deg e avaliação estatística

1. Intersecte cada lista predita específica de metabólito com a lista de DEGs no nível gênico com FDR < 0,05. Registre os genes sobrepostos, o metabólito de origem, logFC, valor-P ajustado e a direção da expressão.
2. Una as listas de sobreposição específicas por metabólito em uma lista não redundante de MAGs. Conte o número total de alvos preditos, sobreposições específicas por metabólito e o total de MAGs únicos.
3. Avalie a consistência direcional no nível dos sondas para genes com múltiplas sondas. Destaque qualquer gene cujas sondas discordem quanto à direção da expressão.
4. Construa a tabela de contingência do teste exato de Fisher utilizando o número total de entradas agrupadas por gene, número total de DEGs, número total de alvos preditos únicos e o número de MAGs observados. Calcule o valor-P unicaudal, a razão de chances e o intervalo de confiança de 95% usando a implementação do teste exato de Fisher.
5. Se a taxa de DEGs na base exceder 50%, relate a sobreposição como descritiva, e não como enriquecimento validado independentemente. Trate a downregulação uniforme como um padrão direcional descritivo, a menos que um teste de direcionalidade separado seja realizado.

6. Análise da Rede de Interações Proteína-Proteína e Enriquecimento de Vias

1. Envie a lista completa e única de MAGs para a construção da rede de interação proteína-proteína e enriquecimento de vias²⁰ (veja Tabela de Materiais). Selecione Homo sapiens e defina a pontuação mínima de interação como 0,700.
2. Exporte a rede combinada de construção de interação proteína-proteína e enriquecimento de vias e a tabela completa de interações. Se a mineração de texto produzir uma topologia artificialmente densa, desmarque a mineração de texto e mantenha os canais experimentais, de coexpressão e de banco de dados.
3. Gere subredes por classe de metabólito para MAGs associados a AGCC, associados a ácidos biliares e a triptamina/serotoninérgicos. Utilize as mesmas configurações de organismo e confiança da construção de rede de interação proteína-proteína e enriquecimento de vias.
4. Execute a construção de rede de interação proteína-proteína e enriquecimento de vias contra o Kyoto Encyclopedia of Genes and Genomes (KEGG)²¹, Reactome²² e Gene Ontology (GO) Biological Process^23,24. Aplique Benjamini–Hochberg BH FDR <0,05 e exporte todas as tabelas de enriquecimento.

7. Encaixe molecular

1. Recupere estruturas de receptores determinadas experimentalmente a partir do Research Collaboratory for Structural Bioinformatics Protein Data Bank RCSB PDB²⁵ (consulte a Tabela de Materiais). Utilize VDR/1DB1, NR1H4/FXR/3DCT, CASP3/2DKO e HTR2A/6A93 para os cinco complexos proteína-ligante prioritários.
2. Prepare cada receptor mantendo a cadeia A e removendo moléculas de água, ligantes co-cristalizados, cofatores, íons e registros HETATM não proteicos. Para o 6A93, remova o segmento de fusão da lisozima T4 antes da preparação do receptor.
3. Adicione hidrogênios polares, atribua cargas de Gasteiger e salve cada receptor no formato PDBQT utilizando o kit de ferramentas de estrutura molecular. Inspeccione os estados de protonação das histidinas no sítio de ligação antes da conversão para PDBQT e documente os estados selecionados.
4. Gere a estrutura tridimensional de cada ligante utilizando o kit de ferramentas de conversão de arquivos de estrutura química. Minimize a energia com o Universal Force Field (UFF) por 500 passos, atribua o estado de ionização em pH 7,4, atribua cargas de Gasteiger e salve no formato PDBQT.
5. Defina uma caixa de acoplamento de 25 Å x 25 Å x 25 Å centrada no centróide do ligante co-cristalizado. Utilize os centros (10, 19, 33) para VDR, (137, 31, 78) para FXR, (37, 34, 32) para CASP3 e (12, −1, 61) para HTR2A.
6. Realize um acoplamento molecular de ligantes metabólitos às proteínas-alvo^26,27 com exaustividade = 8, semente = 42, num_modes = 9 e energy_range = 3 kcal/mol. Registre a pontuação Vina mais alta e os valores de desvio quadrático médio (RMSD) para todas as conformações.
7. Selecione o modo 1 para cada complexo prioritário. Gere diagramas bidimensionais de ligante-resíduo em uma ferramenta de visualização molecular e de diagrama de interação bidimensional, e vistas tridimensionais do receptor-ligante em um programa de visualização molecular tridimensional.
8. Realize controles de redocking para VDR/1DB1 e FXR/3DCT. Aceite a configuração de acoplamento do receptor quando o RMSD de átomos pesados for <2,0 Å em relação à conformação cristalográfica.
9. Realize controles de cross-docking acoplando o ácido litocólico (LCA) ao CASP3 e a triptamina ao VDR. Compare as pontuações cognatas e não cognatas e registre os casos em que a diferença de pontuação for <1,0 kcal/mol.

8. Simulação de dinâmica molecular

1. Gere parâmetros de ligante com o serviço de parametrização de ligantes compatível com CHARMM²⁸ (consulte a Tabela de Materiais). Inspecione todas as pontuações de penalidade e marque qualquer parâmetro com penalidade >50.
2. Converta arquivos de fluxo de ligante em arquivos .itp e .prm compatíveis com o motor de dinâmica molecular usando o script de conversão de topologia de campo de força. Combine os arquivos de topologia do ligante e da proteína para cada complexo.
3. Aplicar a redistribuição da massa do hidrogênio com a biblioteca de conversão de parâmetros e topologia molecular. Gere topologias aquosas utilizando o campo de força da proteína empregado em simulações de dinâmica molecular²⁹ e um modelo explícito de água com três sítios³⁰.
4. Solvate os complexos aquosos em uma caixa dodecaédrica com pelo menos 1,2 nm de distância entre o soluto e a borda. Neutralize os sistemas e adicione NaCl até 0,15 M.
5. Construa o sistema de membrana tryptamina-HTR2A com uma ferramenta baseada na web para construção de sistemas de membrana^31,32,33 (consulte a Tabela de Materiais). Utilize coordenadas do receptor alinhadas ao banco de dados de orientação de proteínas de membrana³⁴ (consulte a Tabela de Materiais), uma bicamada de POPC pura, camadas de água de 22,5 Å e NaCl 0,15 M.
6. Minimize energeticamente todos os sistemas por descida mais íngreme por até 50.000 passos. Confirme a convergência com Fmáx <1000 kJ·mol^-1·nm^-1 antes da equilibração.
7. Equilibre os sistemas aquosos com os conjuntos de Número constante de partículas, volume e temperatura (NVT) e Número constante de partículas, pressão e temperatura (NPT). Equilibre o sistema de membrana utilizando o fluxo de trabalho multietapa baseado na web com seis etapas e restrições liberadas gradualmente.
8. Execute simulações de dinâmica molecular produtiva de 200 ns para os cinco complexos. Utilize um passo de tempo de 4 fs com redistribuição de massa do hidrogênio (HMR), termostato V-rescale a 310 K, barostato de Parrinello-Rahman a 1 bar, eletrostática por Ewald em malha de partículas (PME)³⁵ e vínculos LINCS³⁶.
9. Analise as trajetórias finais com utilitários de análise de trajetórias de dinâmica molecular. Calcule o RMSD da cadeia principal, flutuação quadrática média dos átomos Cα (RMSF), raio de giração, área de superfície acessível ao solvente (SASA) e ligações de hidrogênio entre proteína e ligante, utilizando os últimos 150 ns como janela principal de análise.

9. Estimativa da energia livre de ligação por MM-PBSA

1. Extraia instantâneos da trajetória para análise MM-PBSA. Utilize 2.001 quadros para cada complexo aquoso e 201 quadros processados para o subsistema HTR2A embutido na membrana.
2. Execute a ferramenta de cálculo da energia de ligação em mecânica molecular/solvente contínuo³⁷ com solvatação de Poisson-Boltzmann, constante dielétrica interna = 1, constante dielétrica externa = 80, solvatação não polar baseada em SASA e sem correção de entropia. Informe a energia livre média de ligação e o desvio padrão.
3. Realize a decomposição por resíduo para os cinco complexos. Relate os resíduos estabilizadores e desestabilizadores com contribuições absolutas ≥0,5 kcalmol⁻¹.

    

Resultados

Alvos associados a metabólitos candidatos

Os nove metabólitos produziram conjuntos heterogêneos de alvos previstos em um programa de predição de alvos de interação químico-proteína e em um programa de docking molecular. O propionato, a triptamina, os ácidos biliares e a urolitina A geraram diversos alvos com relevância conhecida para a sinalização gastrointestinal. A paisagem de alvos previstos incluiu receptores de membrana canônicos, receptores nucleares, enzimas intracelulares, proteínas de sinalização e proteínas relacionadas a hormônios peptídicos. Os resultados subsequentes são, portanto, descritos como genes associados a metabólitos (GAM), em vez de achados exclusivos de receptores (Tabela 1).

Avaliação comparativa com interações metabólito-proteína relatadas

Para comparar a saída da predição de alvos com o conhecimento experimental existente, as relações preditas entre metabolitos e alvos foram classificadas em três níveis de evidência: (i) interações diretas ou em nível de classe próximas entre metabolito e proteína com suporte experimental, nas quais o metabolito ou um metabolito endógeno estreitamente relacionado foi relatado como ligante, ativador, inibidor ou regulador funcional da proteína codificada; (ii) interações apoiadas por via metabólica ou classe de alvo, nas quais o alvo predito pertence a uma via metabólica ou família de receptores sensível ao metabolito bem estabelecida, mas há evidência direta limitada para o par exato de metabolito-proteína; e (iii) associações baseadas apenas em computação, para as quais nenhuma interação experimental direta foi identificada na literatura revisada. Essa comparação foi utilizada para contextualizar, e não validar, os MAGs preditos.

Várias previsões recapitularam biologias previamente relatadas. A interação propionato-FFAR2 foi considerada apoiada experimentalmente porque FFAR2/GPR43 é um receptor canônico de ácidos graxos de cadeia curta. A interação butirato-HDAC3 foi classificada como apoiada experimentalmente ou por evidência indireta, porque o butirato é um inibidor reconhecido de histona desacetilase, e a sobreposição prevista envolvia um membro da família HDAC. As previsões associadas aos ácidos biliares envolvendo NR1H4/FXR e VDR foram consideradas apoiadas pela biologia estabelecida dos receptores nucleares de ácidos biliares, especialmente para ácidos biliares hidrofóbicos como o LCA; as previsões de FXR associadas ao ácido ursodesoxicólico (UDCA) foram interpretadas com cautela porque o UDCA é geralmente um ligante mais fraco ou dependente do contexto para FXR. As previsões associadas à triptamina envolvendo HTR1B, HTR2A, HTR2B e HTR6 foram classificadas como apoiadas pela via serotoninérgica, e não por interações diretas e específicas confirmadas, porque a triptamina é uma monoamina derivada microbiana do triptofano, e os receptores de serotonina são reguladores estabelecidos da motilidade e secreção gastrointestinais. A interação urolitina A-CASP3 foi considerada apoiada por vias descritas na literatura entre a urolitina A e respostas apoptóticas ou relacionadas a caspases, mas não por evidência direta de ligação com CASP3. As interações ácido indol-3-lático-KYAT1 e ácido indol-3-propiônico-KYAT1 foram mantidas como hipóteses exclusivamente computacionais porque a literatura mais ampla apoia a sinalização no hospedeiro por derivados microbianos de índole, mas não a ligação direta desses metabólitos exatos com KYAT1^7,8,38,39,40.

Assim, a Tabela 1 distingue a indicação computacional de alvos do nível de evidência experimental prévia ou de suporte por vias conhecidas. Também fornece, para cada alvo, a fonte da previsão (uma previsão de alvo baseada em interação químico-proteína, um programa de encaixe molecular ou ambos), a pontuação combinada de interação para a previsão de alvo por interação químico-proteína e a probabilidade do programa de encaixe molecular quando o alvo foi identificado por esse método. Os alvos previstos sem evidência experimental direta prévia são descritos como genes candidatos associados a metabólitos que requerem validação independente ao nível proteico e de resposta ao ligante.

Sobreposição entre os alvos previstos e os genes diferencialmente expressos em IBS-C

A intersecção entre as listas unidas de alvos preditos e os resultados de expressão diferencial ao nível gênico identificou 17 genes únicos preditos como associados a metabólitos que apresentaram expressão diferencial significativa na comparação entre IBS-C e voluntários saudáveis. Os 17 genes estavam todos regulados negativamente. O conjunto incluiu receptores de membrana e nucleares (CASR, FFAR2, GPR68, HTR1B, HTR2A, HTR2B, HTR6, NR1H4, TBXA2R, VDR) e proteínas não receptoras (CASP3, GCG, GNAQ, GPHN, HDAC3, KYAT1, MLN) (Tabela 1, Figura 2A,B).

Todos os 17 MAGs atenderam a um limiar de taxa de falsa descoberta (FDR) abaixo de 0,05; 16 dos 17 atenderam ao FDR mais rigoroso < 0,001, com o gene restante (HTR1B) significativo em FDR < 0,05. Sete dos 17 alvos (CASP3, GCG, GNAQ, GPHN, GPR68, HDAC3, TBXA2R) satisfizeram tanto FDR < 0,001 quanto uma mudança absoluta no log₂ superior a 1,0 (intervalo de logFC de −1,34 a −1,10), indicando uma subregulação forte e consistente para este subconjunto. Os alvos restantes exibiram subregulação moderada, mas estatisticamente significativa (|logFC| variando de 0,45 a 0,97). Esse padrão descritivo uniforme foi interpretado com cautela, considerando as características de expressão em todo o genoma do conjunto de dados (ver avaliação estatística abaixo).

Avaliação estatística da sobreposição entre alvo e DEG

Para avaliar formalmente a significância estatística da sobreposição de 17 genes, um teste exato de Fisher unicaudal foi aplicado utilizando os 17 genes-alvo preditos como conjunto de consulta e todas as 18.296 entradas únicas colapsadas por gene detectadas em GSE36701 como fundo genômico. Deste fundo, 17.296 genes (94,5%) foram diferencialmente expressos com FDR < 0,05, refletindo uma supressão transcricional quase universal na comparação da mucosa retal de IBS-C. Todos os 17 genes-alvo preditos estavam entre os genes diferencialmente expressos (sobreposição observada 17/17, 100%). Considerando a taxa de 94,5% de expressão diferencial no fundo, a sobreposição esperada para qualquer conjunto aleatório de 17 genes é de 16,1 genes. O teste exato de Fisher resultou em p = 0,384 com uma razão de chances corrigida para continuidade de 2,03 (intervalo de confiança de 95% 0,12–33,73), o que não foi estatisticamente significativo em α = 0,05 (Figura 3A–C).

Esse resultado indica que a sobreposição observada de 17/17 não excede a sobreposição esperada por acaso segundo o perfil de expressão genômica deste conjunto de dados. Assim, esses achados são interpretados como um padrão descritivo direcional, no qual todos os 17 alvos previstos foram consistentemente e significativamente subregulados no tecido da mucosa retal de IBS-C, ao invés de constituírem evidência de enriquecimento estatístico ou validação independente em relação a um fundo genômico. Um teste formal de enriquecimento exigiria replicação em conjuntos de dados transcriptômicos com perfis de expressão diferencial mais seletivos, nos quais significativamente menos da metade de todos os genes alcance significância. Deve-se enfatizar que a subregulação uniforme de todos os 17 genes sobrepostos é uma observação descritiva, e não um resultado estatístico validado separadamente, pois o fundo de expressão diferencial deste conjunto de dados é ele próprio predominantemente subregulado, sendo esperada uma direção descendente compartilhada entre os genes sobrepostos, a qual não foi submetida a um teste formal de direcionalidade. Essa direção uniforme não deve, portanto, ser interpretada como evidência estatística independente de regulação coordenada e específica ao metabólito.

Padrões específicos de metabólitos

O propionato apresentou o maior número de genes sobrepostos, incluindo CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN e TBXA2R, sugerindo possível envolvimento da sinalização mediada por ácidos graxos de cadeia curta e associada à proteína Gq. O butirato apresentou sobreposição com HDAC3, compatível com a biologia do butirato relacionada à desacetilase de histonas, embora a simples redução do mRNA não estabeleça alteração na resposta ao butirato. As sobreposições associadas aos ácidos biliares incluíram os receptores nucleares VDR e NR1H4, ambos reconhecidos como efetores da sinalização dos ácidos biliares no intestino^38,39. A triptamina apresentou sobreposição com HTR1B, HTR2A, HTR2B e HTR6, implicando a sinalização serotoninérgica como um módulo candidato, um sistema com papéis bem estabelecidos na motilidade e secreção gastrointestinais⁴⁰. O ácido indol-3-lático e o ácido indol-3-propiônico apresentaram sobreposição com KYAT1, e a urolitina A apresentou sobreposição com CASP3.

Enriquecimento de vias

A análise de enriquecimento funcional dos 17 genes sobrepostos identificou vias relacionadas à sinalização downstream do receptor acoplado à proteína G (GPCR), sinalização Gαq, ligação de ligantes ao GPCR, sinapse serotoninérgica, interação entre ligantes neuroativos e receptores, transdução do sinal de cálcio, sinalização de cAMP e secreção de hormônios peptídicos. Esses resultados são consistentes com a composição do conjunto de genes e apoiam sua coerência biológica, mas refletem a anotação funcional dos genes submetidos, e não evidência independente de atividade em nível de via.

Estrutura da rede de interação proteína-proteína

A construção da rede de interação entre proteínas e a análise de enriquecimento de vias foram interpretadas em três redes complementares. Na meta-rede combinada de 17 genes (Rede 1), a estrutura mais evidente com suporte de anotação foi um componente centralizado em GNAQ do sinal GPCR/Gαq, ligando GNAQ a genes associados a receptores, incluindo TBXA2R, CASR, HTR2A e HTR2B. A conectividade limitada entre receptores de serotonina também foi mantida, mais proeminentemente entre HTR2A e HTR2B, enquanto vários outros genes permaneceram isolados ou fracamente conectados no limiar de confiança selecionado. A rede específica ao propionato (Rede 2) apresentou uma topologia mais restrita, com GNAQ mantendo ligações com suporte de anotação para CASR e TBXA2R, enquanto FFAR2, GPR68, GCG, GPHN e MLN estavam isolados ou fracamente conectados. A rede de triptamina/serotonina (Rede 3) incluiu HTR1B, HTR2A, HTR2B e HTR6; nesse subconjunto, HTR2A e HTR2B mostraram a principal conexão com suporte de anotação, enquanto HTR1B e HTR6 não estavam diretamente conectados no limiar escolhido (Figura 4A–C).

Encaixe molecular

O encaixe molecular foi realizado em cinco complexos metabólito-proteína selecionados. Os pares ácido biliar-receptor nuclear apresentaram pontuações Vina mais favoráveis do que urolitina A-CASP3 e triptamina-HTR2A. LCA-VDR obteve a melhor pontuação, de −10,0 kcal/mol, seguida por LCA-NR1H4/FXR (−9,9 kcal/mol) e UDCA-NR1H4/FXR (−9,4 kcal/mol). A urolitina A-CASP3 e a triptamina-HTR2A apresentaram pontuações menores, mas ainda razoáveis, de −7,1 kcal/mol (Tabela 2).

Para o complexo LCA-VDR (PDB ID: 1DB1), a pose prevista foi sustentada por uma ligação de hidrogênio convencional entre o oxigênio carboxilato da LCA e Ser278 (4,29 Å), juntamente com extensos contatos hidrofóbicos envolvendo Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233 e His397, além de contatos adicionais de van der Waals com Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275 e Phe150. A pose com maior classificação apresentou um escore Vina de −10,0 kcal/mol, um tamanho de cavidade de 2055 ų e um centro da grade em (10, 19, 33) (Tabela 3, Figura 5A,B).

Para o complexo LCA-NR1H4/FXR (ID PDB: 3DCT), a pontuação de acoplamento de −9,9 kcal/mol foi acompanhada por ligações de hidrogênio previstas envolvendo His294 e Ile335, uma interação π-Sigma com His294 e contatos hidrofóbicos Alquila ou π-Alquila envolvendo Met290, Met328, Ala291, Leu287, Ile352 e His447, com contatos adicionais de van der Waals que apoiam acomodação do esqueleto esteroide na bolsa do FXR (Tabela 4, Figura 6A,B).

A pose prevista do complexo UDCA-NR1H4/FXR (PDB ID: 3DCT) mostrou uma ligação de hidrogênio convencional com His447 (3,66 Å), outra ligação de hidrogênio com Gly322 (3,46 Å), uma interação π-ânion com Val325 (4,96 Å) e uma ligação carbono-hidrogênio com Trp469 (4,51 Å). O mapa de interações também identificou contatos doadores-doador desfavoráveis com Arg395 (3,89 Å) e Gln396 (3,40 Å), sugerindo que a pontuação mais baixa do Vina para a UDCA em comparação com a LCA no mesmo sítio receptor pode ser devida a uma geometria local ou eletrostática menos favorável (Tabela 5, Figura 7A,B).

No complexo urolitina A-CASP3 (ID PDB: 2DKO), o modo de ligação previsto apresentou ligações de hidrogênio convencionais com Gln161 (3,78 e 4,19 Å), Ser120 (3,95 Å) e Arg207 (3,05 e 3,77 Å), sendo ainda estabilizado por interações π-cátion com Arg207, uma ligação de hidrogênio doadora π com Cys163, e contatos adicionais π-alquila e de van der Waals envolvendo Arg64, Ala162, His121, Ser205 e Trp206 (Tabela 6, Figura 8A,B).

Para o complexo triptamina-HTR2A (ID PDB: 6A93), a pose prevista foi estabilizada por uma ponte salina eletrostática entre a amina protonada da triptamina e Asp155, o aspartato conservado na hélice transmembrana 3 (D3.32 na numeração de Ballesteros-Weinstein) que fixa a amina protonada de ligantes aminérgicos em receptores de serotonina e relacionados^41,42,43, juntamente com ligações de hidrogênio com Thr160 e Ser159, contatos aromáticos com Phe340 e Trp336, e interações π-alquila com Val156 e Ile163. Contatos adicionais de van der Waals com Tyr370, Phe339, Ser242, Phe243, Phe332 e Leu123 sustentaram um padrão de ligação ao sítio ortostérico (Tabela 7, Figura 9A,B).

Validação do protocolo de acoplamento molecular

Para avaliar a confiabilidade do protocolo de docking, dois experimentos controle complementares foram realizados. Para os controles de redocking (positivos), ligantes co-cristalizados foram extraídos de suas estruturas de referência por raios X e recolocados por docking em seus sítios de ligação nativos. A conformação prevista com maior classificação para o análogo da vitamina D, VDX, no complexo VDR/1DB1 apresentou um desvio de 0,87 Å em relação à posição cristalográfica, e o ligante do co-cristal WAY-362450 no complexo FXR/3DCT apresentou um desvio de 1,79 Å; ambos os valores ficaram abaixo do limite convencional de aceitação de 2,0 Å, o que sustenta a validade geométrica do protocolo de docking para esses sistemas de receptores (Figura 10A,B). Para os controles de cross-docking (negativos), o ácido litocólico foi submetido a docking na caspase-3 (2DKO), uma protease de cisteína para a qual ele não é um ligante conhecido, resultando em uma pontuação prevista (−8,3 kcal/mol) 1,7 kcal/mol mais fraca do que no seu alvo cognato VDR (−10,0 kcal/mol), o que é consistente com a seletividade prevista para o sítio de ligação. A triptamina submetida a docking no VDR resultou em uma pontuação prevista de −6,4 kcal/mol em comparação com −7,1 kcal/mol no seu alvo cognato HTR2A, uma diferença de 0,7 kcal/mol que está dentro da incerteza relatada para pontuações de docking molecular de ligantes metabólitos com proteínas-alvo, indicando, portanto, apenas uma seletividade prevista modesta para este ligante menor (Figura 10C). Em conjunto, esses controles indicam que o protocolo de docking reproduz geometrias de ligação conhecidas e discrimina pares cognatos de não cognatos nas condições testadas, embora as previsões computacionais restantes não substituam as medidas experimentais de afinidade (Tabela 8).

Simulação de dinâmica molecular

Simulações de dinâmica molecular foram realizadas para os cinco complexos priorizados ao longo de trajetórias de produção de 200 ns. Os quatro complexos solúveis e de receptores nucleares foram simulados em solvente aquoso explícito, enquanto o complexo triptamina-HTR2A foi simulado em uma bicamada lipídica POPC explícita para fornecer um ambiente de membrana fisiologicamente apropriado para este receptor acoplado à proteína G. As análises testaram a estabilidade dinâmica das poses encaixadas sob condições dependentes do tempo e permitiram a comparação do comportamento estrutural relativo entre os complexos (Table 9).

O perfil de RMSD do complexo LCA-VDR/1DB1 mostrou um curto período de equilibração durante os primeiros 10 ns, seguido por uma fase estável, com flutuações principalmente na faixa de 0,20–0,28 nm (Figura 11A). Os valores de RMSF foram baixos, e as flutuações da cadeia principal foram < 0,15 nm para a maioria dos resíduos (Figura 11B). A análise de ligações de hidrogênio revelou uma rede persistente de 2–5 ligações de hidrogênio, com aumentos ocasionais para 7 (Figura 11C). O raio de giração (Rg) manteve-se na faixa de 1,25–1,75 nm, e a área de superfície acessível ao solvente (SASA) permaneceu em torno de 130 nm² (Figura 11D,E).

O complexo urolitina A-CASP3/2DKO exibiu maior atividade dinâmica. O RMSD aumentou inicialmente e depois oscilou entre 0,4 e 0,7 nm, com um breve evento de alta variação por volta de 165 ns (Figura 12A). A análise de RMSF mostrou alta mobilidade ao nível dos resíduos, com as maiores flutuações na região do loop flexível ao redor do resíduo 175 (Figura 12B). A análise de ligação de hidrogênio revelou uma rede inicial extensa de cerca de 2 a 5 ligações nos primeiros 30 a 40 ns, seguida principalmente por 0 a 2 ligações intermitentes (Figura 12C). Os perfis correspondentes do raio de giração e da SASA são apresentados na Figura 12D,E.

Para os sistemas de ácidos biliares NR1H4/FXR (3DCT), o perfil de RMSD do esqueleto permaneceu dentro de uma faixa relativamente estreita durante a maior parte da trajetória (Figura 13A), enquanto o perfil de RMSF mostrou menor mobilidade nas regiões centrais e flutuações mais elevadas nas regiões flexíveis (Figura 13B). O complexo LCA-3DCT manteve aproximadamente de três a quatro ligações de hidrogênio persistentes ao longo de toda a trajetória, enquanto o complexo UDCA-3DCT exibiu maior flutuação nas ligações de hidrogênio e uma redução no número de ligações após aproximadamente 125 ns. Os perfis de raio de giração para os sistemas ligados ao LCA e ao UDCA são mostrados na Figura 13C,D, respectivamente, e os perfis correspondentes de SASA são apresentados na Figura 13E,F.

Dinâmica molecular da membrana do complexo triptamina-HTR2A

O complexo triptamina-HTR2A/6A93 foi simulado por 200 ns em uma bicamada lipídica explícita de POPC composta por 258 moléculas de lipídios, um modelo de água explícito de três sítios e 0,15 M de NaCl, totalizando aproximadamente 100.925 átomos^33,44,45. O receptor permaneceu estável e incorporado na bicamada ao longo de toda a trajetória (Figura 14). O RMSD da cadeia principal aumentou de aproximadamente 0,10 nm para uma platô estável próximo a 0,15–0,20 nm nos primeiros 100 ns e permaneceu estável posteriormente, com todos os valores abaixo de 0,25 nm, indicando que o receptor manteve uma conformação estável no ambiente da membrana sem desdobramento global (Figura 15A). O RMSF por resíduo mostrou baixas flutuações no núcleo helicoidal transmembrana, com mobilidade esperada maior nas regiões de alças e terminais, consistente com a flexibilidade típica de GPCR (Figura 15B). O raio de giração foi rigidamente confinado entre aproximadamente 2,06 e 2,12 nm, e a SASA flutuou dentro de uma faixa estreita sem deriva progressiva, confirmando ambos a preservação do feixe transmembrana compacto (Figura 15C,D).

As ligações de hidrogênio entre a proteína e o ligante foram mantidas ao longo de toda a trajetória (Figura 15E), com grandes flutuações no número de ligações de hidrogênio, variando de 1 a 3. Para avaliar especificamente a persistência da interação iônica chave, monitorou-se ao longo de toda a trajetória a distância mínima entre o nitrogênio do grupo amônio protonado da triptamina e os átomos de oxigênio carboxilato da Asp155 (D3.32). Essa distância permaneceu estritamente distribuída em torno de uma média de 0,270 nm (mínimo de 0,247 nm, máximo de 0,424 nm), e o contato da ponte salina (< 0,4 nm) foi mantido durante 99,9% da simulação, com apenas duas breves excursões transitórias e nenhum evento sustentado de dissociação (Figura 16). Esses resultados sugerem que a interação iônica conservada com a Asp155 foi suficiente para estabilizar a triptamina dentro do sítio ortostérico da HTR2A ao longo de toda a simulação na membrana.

Energia livre de ligação MM-PBSA e decomposição por resíduo

A análise MM-PBSA foi realizada para adicionar uma camada adicional de priorização energética aos cinco complexos (Tabela 10). Para os quatro complexos aquosos, a decomposição por resíduo identificou os principais contribuintes energéticos para cada modo de ligação previsto. No complexo LCA-VDR/1DB1, o ligante e Gln317 apresentaram contribuição favorável, enquanto Trp286 apresentou contribuição desfavorável. No complexo urolitina A-CASP3/2DKO, Arg64 e Arg207 exibiram contribuições por resíduo fortemente negativas, indicando uma estabilização polar ou eletrostática substancial; contudo, a trajetória correspondente permaneceu altamente dinâmica, demonstrando que a termodinâmica favorável ao nível dos resíduos isoladamente não garante estabilidade sustentada do complexo. Para os sistemas 3DCT, a ligação de LCA foi principalmente impulsionada por Arg331, enquanto a ligação de UDCA envolveu uma rede energética mais distribuída, composta por Glu326, Asp394, Arg395, Arg441 e Asp470. Nos quatro sistemas aquosos, a decomposição MM-PBSA apoiou a priorização relativa dos complexos baseados em LCA.

Para o complexo triptamina-HTR2A/6A93 embutido na membrana, a análise MM-PBSA foi realizada no subsistema proteína-ligante extraído da trajetória na bicamada^46,47. Foram observadas contribuições favoráveis para o ligante e Asp155 (D3.32), que foi de longe o resíduo com maior contribuição estabilizadora, consistente com a interação por ponte de sal identificada nas análises de docking e de distância ao longo da trajetória. Trp137 exibiu a maior contribuição desfavorável por resíduo entre os resíduos circundantes do bolso ortostérico (Ser86, Phe87, Phe133, Phe140, Phe141, Val156, Ser159, Thr160, Ile163, Val167, Tyr171), que juntos formam a rede de contatos aromáticos e polares que reveste o bolso de ligação. Esses valores representam estimativas computacionais relativas para priorização estrutural e não são afinidades de ligação experimentais.

Figura 1: Fluxo de trabalho computacional para a priorização de genes hospedeiros associados a metabólitos na IBS-C. Representação esquemática do fluxo de trabalho em oito etapas, integrando seleção de metabólitos, predição de alvos, expressão diferencial transcriptômica, análise de sobreposição, enriquecimento de rede e de vias, encaixe molecular, simulação de dinâmica molecular e análise da energia livre de ligação MM-PBSA. Clique aqui para visualizar uma versão maior desta figura.

Figura 2: Análise da expressão diferencial e sobreposição de alvos metabólitos na mucosa de IBS-C. (A) Gráfico de vulcão da expressão diferencial em nível gênico em GSE36701. Pontos azuis, genes significativamente subregulados; pontos vermelhos, genes significativamente superregulados; pontos cinzas, genes não significativos. Genes associados a metabólitos sobrepostos selecionados estão indicados. (B) Diagrama de Venn mostrando a sobreposição entre 330 alvos de metabólitos únicos previstos e genes subregulados em GSE36701; 17 genes foram compartilhados. Clique aqui para visualizar uma versão maior desta figura.

Figura 3: Avaliação estatística dos 17 genes-alvo de metabólitos preditos em relação ao GSE36701. (A) Variação do log2 por gene para todos os 17 genes, coloridos por categoria de significância. (B) Taxa de expressão diferencial de genes de fundo versus alvos preditos, com teste exato de Fisher. (C) Uma tabela de contingência dois por dois é utilizada para o teste exato de Fisher. Todos os 17 alvos foram significativamente regulados negativamente; a sobreposição é interpretada como um padrão descritivo direcional, e não como enriquecimento estatístico. Clique aqui para visualizar uma versão maior desta figura.

Figura 4: Construção da rede de interação proteína-proteína composta e redes de enriquecimento de vias de interação proteína-proteína de genes associados a metabólitos sobrepostos. (A) Rede 1: meta-rede combinada de todos os 17 genes. (B) Rede 2: rede específica de propionato de oito genes (CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN, TBXA2R). (C) Rede 3: rede de triptamina/serotonina de quatro genes (HTR1B, HTR2A, HTR2B, HTR6). As redes foram geradas para Homo sapiens no mínimo, utilizando confiança de construção de rede de interação proteína-proteína e enriquecimento de vias ≥ 0,700. As arestas representam associação funcional apoiada por anotação Clique aqui para visualizar uma versão maior desta figura.

Figura 5: Representação estrutural tridimensional e bidimensional do ácido litocólico em complexo com o VDR (ID PDB: 1DB1). (A) Representação tridimensional em superfície e esquemática, com o ácido litocólico mostrado como esferas. (B) Mapa de interações bidimensional mostrando a ligação de hidrogênio com a Ser278 e as interações hidrofóbicas e de van der Waals ao redor. Clique aqui para visualizar uma versão maior desta figura.

Figura 6. Representação estrutural tridimensional e bidimensional do ácido litocólico em complexo com NR1H4/FXR (ID PDB: 3DCT). (A) Representação tridimensional em superfície e esquemática. (B) Mapa de interações bidimensional mostrando ligações de hidrogênio com His294 e Ile335, uma interação π-Sigma e contatos adjacentes. Clique aqui para visualizar uma versão maior desta figura.

Figura 7: Representação estrutural tridimensional e bidimensional do ácido ursodesoxicólico em complexo com NR1H4/FXR (ID PDB: 3DCT). (A) Representação tridimensional em superfície e esquemática. (B) Mapa bidimensional de interações mostrando ligações de hidrogênio com His447 e Gly322, uma interação π-ânion com Val325, uma ligação carbono-hidrogênio com Trp469 e contatos doador-doador desfavoráveis com Arg395 e Gln396. Clique aqui para visualizar uma versão maior desta figura.

Figura 8: Representação estrutural tridimensional e bidimensional da urolitina A em complexo com CASP3 (ID PDB: 2DKO). (A) Representação tridimensional em superfície e no formato de fita. (B) Mapa de interações bidimensional mostrando ligações de hidrogênio com Gln161, Ser120 e Arg207, interações π-cátion com Arg207, uma ligação de hidrogênio doadora π com Cys163 e contatos circundantes. Clique aqui para visualizar uma versão maior desta figura.

Figura 9: Representação estrutural tridimensional e bidimensional da triptamina em complexo com HTR2A (ID PDB: 6A93). (A) Representação tridimensional em superfície e esquemática gerada em um programa de visualização molecular tridimensional. (B) Mapa de interações bidimensional gerado em uma ferramenta de visualização molecular e diagrama de interações bidimensional, ilustrando a ponte salina de Asp155 e interações adicionais no sítio de ligação. Clique aqui para visualizar uma versão maior desta figura.

Figura 10: Validação do protocolo de docking. (A,B) Redocking de ligantes co-cristalizados na VDR/1DB1 (RMSD 0,87 Å) e na FXR/3DCT (RMSD 1,79 Å); as conformações cristalográficas e redockeadas estão sobrepostas, ambas abaixo do limite de aceitação de 2,0 Å. (C) Seletividade em cross-docking: pontuações Vina de ligantes cognatos versus não cognatos para o ácido litocólico e a triptamina. Clique aqui para visualizar uma versão maior desta figura.

Figura 11. Análise da trajetória de dinâmica molecular do complexo LCA-VDR/1DB1 ao longo de 200 ns. (A) Perfil de RMSD. (B) Perfil de RMSF. (C) Contagem de ligações de hidrogênio. (D) Perfil de raio de giração. (E) Perfil de SASA. Clique aqui para visualizar uma versão maior desta figura.

Figura 12: Análise da trajetória de dinâmica molecular do complexo urolitina A-CASP3/2DKO ao longo de 200 ns. (A) Perfil de RMSD mostrando flutuações conformacionais amplas e um evento transitório de alta desvio próximo a 165 ns. (B) Perfil de RMSF mostrando flexibilidade acentuada ao nível dos resíduos próximo ao resíduo 175. (C) Contagem de ligações de hidrogênio. (D) Perfil de raio de giração. (E) Perfil de SASA. Clique aqui para visualizar uma versão maior desta figura.

Figura 13: Análise da trajetória de dinâmica molecular dos sistemas de ácidos biliares NR1H4/FXR (3DCT) ao longo de 200 ns. (A) Perfil de RMSD da cadeia principal para o complexo 3DCT. (B) Perfil de RMSF da cadeia principal. (C) Perfil de raio de giração para 3DCT-LCA. (D) Perfil de raio de giração para 3DCT-UDCA. (E) Perfil de SASA para 3DCT-LCA. (F) Perfil de SASA para 3DCT-UDCA. Clique aqui para visualizar uma versão maior desta figura.

Figura 14: O complexo triptamina-HTR2A embutido em uma bicamada lipídica POPC explícita. O receptor é mostrado como um modelo esquemático atravessando a bicamada, os lipídios POPC como linhas com os grupos fosfato das cabeças destacados e a triptamina dentro do sítio ortostérico. A água é mostrada acima e abaixo da membrana. Clique aqui para visualizar uma versão maior desta figura.

Figura 15: Análise da trajetória de dinâmica molecular do complexo triptamina-HTR2A/6A93 ao longo de 200 ns em uma bicamada lipídica POPC explícita. (A) Perfil de RMSD da cadeia principal. (B) Perfil de RMSF por resíduo. (C) Perfil de raio de giração. (D) Perfil de SASA. (E) Contagem de ligações de hidrogênio entre proteína e ligante. Clique aqui para visualizar uma versão maior desta figura.

Figura 16: Persistência da interação iônica entre a triptamina e Asp155 (D3.32) ao longo da trajetória de 200 ns na membrana. A distância mínima entre o nitrogênio amônio da triptamina e os átomos de oxigênio carboxilato do Asp155 é representada em função do tempo; a linha tracejada indica o limite de contato da ponte de sal em 0,4 nm. O contato foi mantido durante 99,9% da simulação. Clique aqui para visualizar uma versão maior desta figura.

Símbolo do Gene	Metabólito(s) de Origem	Categoria Funcional	log2FC	FDR (valor-P ajustado)	Nível de Significância
GCG	Propionato	Proteína relacionada a hormônio peptídico	−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 sinalização intracelular	−1.122	1.05e−6	FDR <0.001 & |logFC| > 1
GPHN	Propionato	Outra 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 a hormônio peptídico	−0.605	7.39e−5	FDR <0.001
KYAT1	Ácido indol-3-lático / Á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

Tabela 1: Genes-alvo associados a metabólitos preditos que se sobrepõem aos genes diferencialmente expressos no conjunto de dados do mucosa retal de IBS-C. Todos os genes sobrepostos listados foram regulados negativamente. A Tabela 1 é enviada separadamente como uma tabela de planilha e lista, para cada alvo, o(s) metabólito(s) de origem, categoria funcional, fonte de predição do alvo (uma predição de alvo por interação químico-proteína, um programa de encaixe molecular ou ambos), um escore combinado de interação da predição de alvo por interação químico-proteína e a probabilidade do programa de encaixe molecular, quando disponível, o nível de predição, a variação de expressão em log2, e o FDR com o nível de significância da expressão. Fonte: os valores de expressão gênica foram obtidos da tabela de expressão diferencial GSE36701 com genes agrupados (sonda com menor FDR por gene). Os valores de fonte de predição de alvo e de confiança foram compilados a partir da saída de uma predição de alvo por interação químico-proteína e de um programa de encaixe molecular, utilizando limiares de escore combinado de interação químico-proteína ≥ 0,700 e probabilidade do programa de encaixe molecular ≥ 0,70. Os escores da predição de alvo por interação químico-proteína são escores combinados numa escala de 0 a 1; STP denota a probabilidade do programa de encaixe molecular. Nível 1 = suporte rigoroso da predição de alvo por interação químico-proteína; Nível 1+ = suporte rigoroso da predição de alvo por interação químico-proteína corroborado por um programa de encaixe molecular.

Complexo	Proteína (ID PDB)	Ligante	Pontuação Vina (kcal/mol)	Tamanho da Cavidade (A^3)	Centro da Grade X,Y,Z (A)	Caixa de Busca (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

Tabela 2: Resultados de docking molecular: melhores classificações de docking molecular de ligantes metabólitos a proteínas-alvo, escores e parâmetros de cavidade para os cinco complexos proteína-ligante prioritários. O tamanho da cavidade é informado em ų. Fonte: Docking_Validation/Results/Docking_Validation_Results.xlsx, planilha 'Original_Docking_Scores'. Um docking molecular de ligantes metabólitos a proteínas-alvo; exaustividade = 8, semente = 42 (fixa), num_modes = 9 para todos os complexos; pose com melhor classificação (modo 1) relatada.

Tipo de Interação	Resíduo(s)	Distância (A)	Notas
Ligação de hidrogênio convencional	Ser278	4.29	Oxigênio carboxilato de LCA
Contato hidrofóbico / Pi-Alquila	Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233, His397	-
Contato de Van der Waals	Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275, Phe150	-

Tabela 3: Modos de ligação gerados para o encaixe do ácido litocólico com VDR (ID PDB: 1DB1). Fonte: ferramenta de visualização molecular e diagramas bidimensionais de interação ligante-resíduo, conforme relatado nos Resultados do manuscrito (Encaixe molecular). '-' indica que o valor da distância não foi relatado individualmente para esse contato.

Tipo de Interação	Resíduo(s)	Distância (A)	Notas
Ligação de hidrogênio	His294	-
Ligação de hidrogênio	Ile335	-
Interação Pi-Sigma	His294	-
Alquila / Pi-Alquila (hidrofóbica)	Met290, Met328, Ala291, Leu287, Ile352, His447	-
Contato de Van der Waals	Resíduos adicionais do bolso (não especificados individualmente na fonte)	-	Contribui para acomodação do esqueleto esteroide

Tabela 4: Modos de ligação gerados para o encaixe do ácido litocólico com NR1H4/FXR (ID PDB: 3DCT).

Fonte: visualização molecular e ferramenta de diagrama de interação bidimensional diagramas de interação ligante-resíduo 2D, conforme relatado nos Resultados do manuscrito (Docking molecular). '-' indica que o valor da distância não foi individualmente relatado para esse contato.

Tipo de Interação	Resíduo(s)	Distância (A)	Observações
Ligação de hidrogênio convencional	His447	3.66
Ligação de hidrogênio	Gly322	3.46
Interação Pi-Ânion	Val325	4.96
Ligação carbono-hidrogênio	Trp469	4.51
Contato doador-doador desfavorável	Arg395	3.89
Contato doador-doador desfavorável	Gln396	3.40

Tabela 5: Modos de ligação gerados para o encaixe do ácido ursodesoxicólico com NR1H4/FXR (ID PDB: 3DCT). Fonte: ferramenta de visualização molecular e diagramas bidimensionais de interação ligante-resíduo, conforme relatado nos Resultados do manuscrito (Encaixe molecular). '-' indica que o valor da distância não foi individualmente relatado para esse contato.

Tipo de Interação	Resíduo(s)	Distância (A)	Notas
Ligação de hidrogênio convencional	Gln161	3.78
Ligação de hidrogênio convencional	Gln161	4.19	segundo contato
Ligação de hidrogênio convencional	Ser120	3.95
Ligação de hidrogênio convencional	Arg207	3.05
Ligação de hidrogênio convencional	Arg207	3.77	segundo contato
Interação Pi-Cátion	Arg207	-
Ligação de hidrogênio doadora Pi	Cys163	-
Contato Pi-Alquila / van der Waals	Arg64, Ala162, His121, Ser205, Trp206	-

Tabela 6: Modos de ligação gerados para o encaixe molecular da urolitina A com CASP3 (ID PDB: 2DKO). Fonte: ferramenta de visualização molecular e diagramas bidimensionais de interação ligante-resíduo, conforme relatado nos Resultados do manuscrito (Encaixe molecular). '-' indica que o valor da distância não foi relatado individualmente para esse contato.

Tipo de Interação	Resíduo(s)	Distância (A)	Observações
Ponte salina eletrostática	Asp155 (D3.32)	-	amina protonada da triptamina
Ligação de hidrogênio	Thr160	-
Ligação de hidrogênio	Ser159	-
Contato aromático	Phe340, Trp336	-
Interação Pi-Alquila	Val156, Ile163	-
Contato de Van der Waals	Tyr370, Phe339, Ser242, Phe243, Phe332, Leu123	-

Tabela 7: Modos de ligação gerados para o encaixe molecular da triptamina com HTR2A (ID PDB: 6A93). Fonte: ferramenta de visualização molecular e diagramas bidimensionais de interação ligante-resíduo, conforme relatado nos Resultados do manuscrito (Encaixe molecular). '-' indica que o valor da distância não foi individualmente informado para esse contato.

(A) Validação por redocking (controles positivos)
ID PDB	Proteína	Ligante Cocrystalizado	Pontuação Vina (kcal/mol)	RMSD (Å)	Limiar (Å)	Resultado
1DB1	VDR	VDX (análogo da vitamina D)	−13.0	0.87	2.0	PASSOU
3DCT	FXR	WAY-362450 (064)	−11.9	1.79	2.0	PASSOU
(B) Validação por cross-docking (controles negativos)
Ligante	Alvo Cognato (PDB)	Pontuação Cognata (kcal/mol)	Alvo Não Cognato (PDB)	Pontuação Não Cognata (kcal/mol)	Delta (kcal/mol)	Seletividade
Á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 da incerteza da Vina +/−0,5–1,0)

Tabela 8: Resultados da validação do protocolo de docking: valores de RMSD de redocking (controles positivos) e pontuações de cross-docking (controles negativos). Fonte: Docking_Validation/Results/Docking_Validation_Results.xlsx e Docking_Validation/Logs/*.log (um docking molecular de ligantes metabólitos a proteínas-alvo, exaustividade = 8, semente = 42, caixa de 25 Å × 25 Å × 25 Å). RMSD calculado por correspondência de átomos pesados e nomes de átomos (sem sobreposição).

Complexo	RMSD (nm), média + / – DP (intervalo)	Rg (nm), média + / – DP (intervalo)	SASA (nm^2), média + / – DP (intervalo)	Ligações de hidrogênio, média + / – DP (intervalo)	RMSF (nm), média (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 no resíduo 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 no resíduo 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)

Tabela 9: Resumo do comportamento da simulação de dinâmica molecular de 200 ns para os cinco complexos proteína-ligante prioritários, incluindo o sistema triptamina-HTR2A embutido na membrana. Fonte: arquivos de utilitários de análise de trajetórias de dinâmica molecular (.xvg) — gmx rms, gmx gyrate, gmx sasa, gmx hbond, gmx rmsf — calculados ao longo dos últimos 150 ns (50–200 ns) de cada simulação produtiva de 200 ns, conforme o passo 8.8 do protocolo. RMSD/Rg ajustados ao esqueleto principal; raio da sonda SASA de 0,14 nm; corte para doador-aceptor de ligação de hidrogênio de 0,35 nm / 30 °. LCA-3DCT e UDCA-3DCT compartilham uma trajetória do esqueleto principal da proteína (RMSD, RMSF) com Rg/SASA/ligações de hidrogênio específicas ao ligante.

Triptamina-HTR2A/6A93 (membrana) — decomposição quantitativa por resíduo
Resíduo	Contribuição total de ddG (kcal/mol), média + / − DP	Direção
Asp155 (D3.32)	−89.94 + / − 6.81	Estabilizante (dominante)
Triptamina (ligante)	−13.01 + / − 6.22	Estabilizante
Tyr171	13.62 + / − 4.54	Desestabilizante
Val167	23.32 + / − 3.96	Desestabilizante
Val156	20.03 + / − 3.81	Desestabilizante
Thr160	4.86 + / − 3.64	Desestabilizante
Ser159	24.16 + / − 3.48	Desestabilizante
Ser86	24.48 + / − 3.65	Desestabilizante
Phe87	35.18 + / − 4.04	Desestabilizante
Phe133	32.80 + / −3.70	Desestabilizante
Phe140	30.63 + / − 3.84	Desestabilizante
Phe141	35.25 + / − 3.55	Desestabilizante
Ile163	27.64 + / − 3.71	Desestabilizante
Trp137	53.77 + / − 4.32	Desestabilizante (maior desfavorável)
Outros quatro complexos — resíduos identificados na decomposição por resíduo (qualitativo)
Complexo	Resíduo	Direção
LCA-VDR/1DB1	Ligante (LCA)	Favorável
LCA-VDR/1DB1	Gln317	Favorável
LCA-VDR/1DB1	Trp286	Desfavorável
LCA-NR1H4/FXR/3DCT	Arg331	Favorável (dominante)
UDCA-NR1H4/FXR/3DCT	Glu326	Rede mista/distribuída
UDCA-NR1H4/FXR/3DCT	Asp394	Rede mista/distribuída
UDCA-NR1H4/FXR/3DCT	Arg395	Rede mista/distribuída
UDCA-NR1H4/FXR/3DCT	Arg441	Rede mista/distribuída
UDCA-NR1H4/FXR/3DCT	Asp470	Rede mista/distribuída
Urolitina A-CASP3/2DKO	Arg64	Forte favorável (polar/eletrostático)
Urolitina A-CASP3/2DKO	Arg207	Forte favorável (polar/eletrostático)

Tabela 10: Decomposição por resíduo do MM-PBSA – RESUMO BREVE: resíduos estabilizadores e desestabilizadores (contribuição absoluta ≥ 0,5 kcal mol⁻¹) para cada um dos cinco complexos prioritários proteína-ligante, incluindo o sistema de triptamina-HTR2A embutido na membrana. Fonte: Membrane Simulation/03_MMPBSA/results/FINAL_DECOMP_MMPBSA.dat (ferramenta de cálculo da energia de ligação por mecânica molecular/solvente contínuo Generalized Born (GB), decomposição por resíduo, 'Complex: Total Energy Decomposition'). Os números dos resíduos foram convertidos da numeração interna do sistema construído pelo CHARMM-GUI (deslocamento +68) para a numeração original do PDB 6A93 utilizada em outras partes deste manuscrito.

Fonte: dados de simulação de dinâmica molecular anteriores/1DB1,2KD0, LCA & UDCA_3DCT}/mmpbsa_*/Decomposition_NORMAL_GB_Complex_TDC*.svg e resultados do manuscrito (energia livre de ligação MM-PBSA e decomposição por resíduo). Esses quatro complexos não possuem saída numérica por resíduo em arquivos .dat/.csv no diretório do projeto (apenas gráficos SVG renderizados com texto em caminhos vetoriais que não são extraíveis por máquina); apenas a identidade dos resíduos e a direção favorável/desfavorável, conforme declarado no texto do manuscrito, são informadas. As contribuições exatas em kcal/mol para esses quatro complexos não estão disponíveis no repositório de origem.

Discussão

Este estudo computacional exploratório demonstra um fluxo de trabalho integrado e reproduzível para priorizar genes do hospedeiro e complexos proteína-ligante associados a metabólitos microbianos, aplicado aqui a um conjunto de dados transcriptômicos públicos de mucosa retal de IBS-C. Utilizando este fluxo de trabalho, um subconjunto de genes preditos como associados a metabólitos microbianos apresentou sobreposição com genes consistentemente subexpressos no conjunto de dados, agrupando-se em vias associadas a GPCR, serotoninérgicas, de sinalização de cálcio, ligante-receptor neuroativo e receptores nucleares, sistemas cada vez mais implicados na comunicação entre microbiota e hospedeiro^48,49. Esses achados devem ser interpretados estritamente como geradores de hipóteses: a análise não mede concentrações de metabólitos microbianos, abundância de proteínas receptoras, ligação do ligante, ativação do receptor, sinalização downstream, motilidade, secreção, respostas à dor ou desfechos clínicos. A conclusão mais robusta sustentável é que os genes e vias identificados são candidatos para validação experimental, e não mecanismos de doença confirmados.

A importância crítica deste protocolo, em comparação com trabalhos anteriores que examinam pares isolados de metabólito-receptor, reside na integração da predição de alvos, transcriptômica pública, análise de rede, acoplamento com controles de validação, dinâmica molecular e MM-PBSA em um único fluxo sequencial de priorização. Cada etapa reduz e contextualiza o conjunto de candidatos gerado pela etapa anterior, sendo o filtro sequencial o que torna a lista final de candidatos viável experimentalmente. O módulo de GPCR centrado em GNAQ e o módulo associado ao receptor de serotonina identificados aqui são biologicamente plausíveis, considerando o papel da sinalização Gq na ativação da fosfolipase C, produção de inositol 1,4,5-trisfosfato, mobilização de cálcio, secreção e função enteroendócrina, bem como os papéis estabelecidos da sinalização derivada de ácidos graxos de cadeia curta e de triptofano na homeostase da mucosa e da sinalização serotoninérgica na motilidade gastrointestinal, secreção, sensibilidade visceral e comunicação intestino-cérebro^10,18,50,51.

Um aspecto metodológico fundamental deste estudo é o tratamento do receptor HTR2A inserido na membrana. Como uma simulação em fase solúvel não pode reproduzir o ambiente lipídico que rege o comportamento conformacional de um receptor acoplado à proteína G, o complexo triptamina-HTR2A foi simulado em uma bicamada POPC explícita. Nesse ambiente de membrana, o receptor manteve-se estruturalmente estável ao longo de toda a trajetória de 200 ns, e a ponte de sal entre o grupo amônio da triptamina e o Asp155 (D3.32) foi mantida essencialmente durante toda a simulação. O fato de três linhas independentes de evidência — a pose de docking, a distância de contato persistente ao longo da trajetória e a contribuição dominante por resíduo no cálculo MM-PBSA — convergirem para a mesma interação conservada com D3.32 confere consistência interna ao modo de ligação da triptamina previsto, o qual recapitula a geometria de ligação canônica de ligantes aminérgicos em receptores de serotonina.

Existem certas questões metodológicas a serem consideradas ao reproduzir este fluxo de trabalho. Erros na estrutura canônica ou compostos com interferência em ensaios multissérie (compostos sinalizados como PAINS) se propagam por meio da predição de alvos e do encaixe molecular, exigindo seleção precisa de metabólitos e curadoria quimioinformática. Minimize conjuntos de alvos gerados por ruído aplicando critérios de confiança de forma consistente (predição de alvos por interação químico-proteica: ≥0,700; programa de encaixe molecular: ≥0,70; construção de rede de interação proteína-proteína e enriquecimento de vias: ≥0,700). Os alvos preditos devem ser agrupados por categoria funcional para evitar a caracterização incorreta de todos os genes associados aos metabólitos como receptores. O pré-processamento preciso das estruturas PDB, a minimização da energia do ligante e o posicionamento da grade ao redor dos resíduos de ligação conhecidos são aspectos essenciais do encaixe molecular, e os controles de reencaixe e encaixe cruzado introduzidos aqui fornecem uma medida objetiva da correção da metodologia de encaixe. O limite de reprodutibilidade na dinâmica molecular é definido por uma combinação de parametrização do campo de força, solvatação adequada ou montagem de membrana, equilibração em etapas e amostragem suficiente na fase produtiva.

Adaptações e etapas de solução de problemas típicas incluem relaxar os limiares caso a previsão do alvo não retorne resultados, verificar a consistência direcional no nível da sonda para genes com múltiplas sondas e interpretar nós isolados de construção de rede de interação proteína-proteína e enriquecimento de vias como dependentes do limiar, e não como biologicamente irrelevantes. Para receptores de membrana, deve-se utilizar simulação explícita de bicamada lipídica em vez de simulação aquosa, conforme exemplificado pelo método HTR2A descrito aqui. Quando for necessária uma decomposição energética por resíduo, o cálculo deve ser realizado com um mecanismo capaz de decomposição, e a numeração dos resíduos informada deve ser compatibilizada com a numeração nativa do receptor para evitar ambiguidade. Sugerimos que os resultados de enriquecimento de vias devam ser tratados principalmente como um contexto organizacional para a lista de candidatos, e não como validação em nível de via. Mecanicamente, o enriquecimento de termos relacionados a RCPG, sinalização serotoninérgica ou sinalização de cálcio ocorrerá sempre que a lista de genes contiver múltiplos genes de receptores de serotonina, independentemente da coregulação em nível proteico. Os valores de RMSD, Rg e RMSF para o sistema HTR2A inserido na membrana devem ser interpretados considerando a bicamada lipídica: uma redução do Rg na trajetória posterior pode refletir uma adaptação conformacional induzida pela bicamada no feixe transmembrana, e não um desdobramento global; além disso, as ligações de hidrogênio persistentes entre ligante e proteína devem ser interpretadas em conjunto com a estabilidade geral do RMSD.

As limitações deste estudo são substanciais e restringem a interpretação. A pesquisa foi baseada em um único conjunto de dados público relativamente pequeno, e uma busca nos principais repositórios públicos de transcriptômica (ferramenta web de análise de expressão gênica diferencial e ArrayExpress) não identificou um conjunto de dados independente de transcriptoma da mucosa retal de pacientes com IBS-C com desenho e plataforma comparáveis que pudesse servir como coorte de replicação no momento da análise. A ausência de replicação transcriptômica independente é uma limitação importante, e nenhuma afirmação neste manuscrito deve ser interpretada como validação externa dos achados obtidos a partir de um único conjunto de dados. O conjunto de dados mostra expressão diferencial quase universal (aproximadamente 94,5% dos genes são significativos, sendo a grande maioria regulados negativamente), uma característica que torna as estatísticas convencionais de enriquecimento não informativas e impede conclusões sobre a especificidade da subexpressão dos genes-alvo em relação ao fundo genômico; portanto, a sobreposição é relatada como um padrão descritivo direcional em vez de enriquecimento estatístico. A transcriptômica de mucosa em bloco não consegue distinguir regulação gênica real de alterações na composição celular. A expressão de mRNA não determina a abundância proteica nem a resposta funcional. Os bancos de dados de predição de alvos sofrem com viés de anotação, e os resultados de docking, MD e MM-PBSA dependem da escolha do campo de força, parametrização do ligante, posição inicial, tempo de simulação e da adequação da amostragem. Identificadores exatos de versão/subversão para alguns componentes de servidores web e pacotes, incluindo o serviço de parametrização de ligantes compatível com CHARMM, CHARMM-GUI e ambiente de computação estatística

as versões dos pacotes criados e as subversões da ferramenta de cálculo de energia de ligação por mecânica molecular/solvente contínuo não puderam ser totalmente recuperadas do registro arquivado do projeto e devem ser relatadas como disponíveis na Tabela de Materiais separada. Os valores de MM-PBSA são estimativas relativas, não contêm um termo explícito de entropia conformacional e não devem ser interpretados como afinidades experimentais. O estudo não possui dados metabolômicos e não pode determinar se a disponibilidade do ligante é alterada no IBS-C, ou se as alterações na expressão observadas são causas, consequências, respostas compensatórias ou correlações irrelevantes.

Aplicações futuras deste método devem incluir replicação transcriptômica independente, validação por reação em cadeia da polimerase quantitativa (qPCR) e em nível proteico, localização por tipo celular mediante transcriptômica de célula única ou espacial, perfilhamento metabolômico das classes de metabólitos relevantes e ensaios funcionais de resposta a ligantes em colonoides derivados de pacientes, explantes da mucosa ou modelos comparáveis. Comparações com coortes de síndrome do intestino irritável com predomínio de diarreia, síndrome do intestino irritável mista, doença inflamatória intestinal e constipação não relacionada à síndrome do intestino irritável^1,2 ajudariam a estabelecer a especificidade da doença. Para o componente estrutural, a repetição das trajetórias de dinâmica molecular (MD), a realização de análises de sensibilidade com poses iniciais alternativas e a documentação completa do depósito da topologia, das trajetórias e dos arquivos de entrada e saída do MM-PBSA reforçariam ainda mais a reprodutibilidade. Ensaios experimentais de resposta a ligantes permanecem necessários para determinar se os complexos priorizados são funcionalmente relevantes; os resultados atuais não sustentam afirmações clínicas ou terapêuticas.

Materiais

Lista de materiais utilizados neste artigo

Nome	Empresa	Número de catálogo	Comentários
AutoDock Vina	Scripps Research / código aberto	v1.2.7; https://vina.scripps.edu/	Acoplamento molecular de ligantes metabólitos a proteínas-alvo.
CGenFF/ParamChem	SilcsBio / University of Maryland	v4.6; https://cgenff.com/	Parametrização do campo de forças do ligante para dinâmica molecular.
Campo de forças CHARMM36m	Desenvolvedores do CHARMM / código aberto	CHARMM36m; https://www.charmm.org/charmm/resources/charmm-force-fields/	Campo de forças proteico utilizado em simulações de dinâmica molecular.
CHARMM-GUI Membrane Builder	CHARMM-GUI / Lehigh University	Servidor web; versão exata não recuperável; https://www.charmm-gui.org/?doc=input/membrane	Configuração da construção e equilibração do sistema de membrana POPC explícita.
Discovery Studio Visualizer	BIOVIA (Dassault Systèmes)	2021; https://discover.3ds.com/discovery-studio-visualizer-download	Análise bidimensional das interações ligante-resíduo.
GEO2R	NCBI Gene Expression Omnibus	Ferramenta web; acesso de jan. a maio de 2026; https://www.ncbi.nlm.nih.gov/geo/geo2r/	Análise de expressão diferencial do GSE36701.
GeneCards	Weizmann Institute of Science	Banco de dados web; acesso de jan. a maio de 2026; https://www.genecards.org/	Verificação de símbolos gênicos e informações genéticas durante a padronização de alvos.
gmx_MMPBSA	Código aberto (Valdés-Tresanco et al.)	1.5.x; https://valdes-tresanco-ms.github.io/gmx_MMPBSA/	Estimativa da energia livre de ligação MM-PBSA e decomposição por resíduo.
GROMACS	Equipe de desenvolvimento do GROMACS / código aberto	2024.2; https://www.gromacs.org/	Mecanismo de simulação de dinâmica molecular.
Conjunto de dados transcriptômicos GSE36701	NCBI Gene Expression Omnibus	GSE36701; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE36701	Conjunto de dados públicos de expressão no epitélio retal de IBS-C.
Open Babel	Código aberto	3.2.0; https://openbabel.org/	Conversão de formatos de arquivos químicos, geração tridimensional de ligantes e preparação de ligantes.
Banco de dados OPM	University of Michigan	Banco de dados web; acesso de jan. a maio de 2026; https://opm.phar.umich.edu/	Coordenadas da Orientação de Proteínas em Membranas usadas para alinhar HTR2A.
ParmEd	Desenvolvedores do ParmEd / código aberto	4.x; https://parmed.github.io/ParmEd/html/index.html	Reparticionamento de massa de hidrogênio e processamento de topologia para simulação molecular.
PyMOL	Schrödinger / código aberto	2.x; https://www.pymol.org/	Visualização estrutural tridimensional e figuras receptor-ligante.
Banco de Dados de Proteínas RCSB	RCSB PDB	Banco de dados web; acesso de jan. a maio de 2026; https://www.rcsb.org/	Fonte de estruturas proteicas experimentais e coordenadas PDB.
STITCH	Consortium STITCH (EMBL)	v5.0; https://stitch.embl.de/	Previsão de alvos de interação químico-proteína.
STRING	Consortium STRING / ELIXIR	v12.0; https://version-12-0.string-db.org/	Construção de rede de interação proteína-proteína e enriquecimento de vias.
SwissADME	SIB Swiss Institute of Bioinformatics / University of Lausanne	Ferramenta web; acesso de jan. a maio de 2026; https://www.swissadme.ch/	Descritores quimioinformáticos, previsões farmacocinéticas e avaliação de PAINS.
SwissTargetPrediction	SIB Swiss Institute of Bioinformatics / University of Lausanne	Ferramenta web; acesso de jan. a maio de 2026; https://www.swisstargetprediction.ch/	Previsão baseada em ligante de alvos proteicos humanos.
Mapeamento de ID UniProt	Consortium UniProt	Serviço web; acesso de jan. a maio de 2026; https://www.uniprot.org/id-mapping	Mapeamento de identificadores proteicos para símbolos gênicos padronizados aprovados pelo HGNC.
NVIDIA RTX 3080	NVIDIA Corporation	RTX 3080; ≥8 GB de VRAM; versão CUDA/driver não especificada no manuscrito	Unidade de processamento gráfico compatível com CUDA usada em simulações de dinâmica molecular.
GPU compatível com CUDA	NVIDIA Corporation	Versão do toolkit CUDA não especificada no manuscrito; ≥8 GB de VRAM	GPU compatível com CUDA com ≥8 GB de VRAM; a estação de trabalho também exigia ≥32 GB de RAM e CPU de 6 núcleos.
Ubuntu Linux	Canonical Ltd. / código aberto	22.04 LTS	Sistema operacional Linux de 64 bits.
Python 3.9	Python Software Foundation	3.9	Ambiente de programação de uso geral usado para scripts de fluxo de trabalho e análise.
Gene Expression Omnibus (GEO)	NCBI / U.S. National Library of Medicine	Repositório web público; nenhuma versão de software especificada no manuscrito	Repositório público de dados de genômica funcional.
AutoDockTools/MGLTools	Molecular Graphics Laboratory, Scripps Research	1.5.7	Ferramenta de preparação de estrutura molecular e entrada para acoplamento.
Ferramentas de análise do GROMACS	Equipe de desenvolvimento do GROMACS / código aberto	2024.2	Utilitários de análise de trajetória de dinâmica molecular.
Protocolo de seis etapas do CHARMM-GUI	CHARMM-GUI / Lehigh University	Protocolo web; versão exata não recuperável	Fluxo de trabalho web baseado em múltiplas etapas para preparação e equilibração de sistema molecular.
cgenff_charmm2gmx_py3.py	Script de conversão de código aberto; fonte não especificada no manuscrito	Versão não especificada no manuscrito	Script de conversão de topologia de campo de forças.
Python	Python Software Foundation	3.9	Ambiente de programação de uso geral.
SciPy	Comunidade SciPy / código aberto	Versão não especificada no manuscrito	Biblioteca de computação científica.
scipy.stats.fisher_exact	Comunidade SciPy / código aberto	Versão SciPy não especificada no manuscrito	Implementação do teste exato de Fisher’s.
R	R Foundation for Statistical Computing	4.3.x	Ambiente de computação estatística.
Bioconductor	Projeto Bioconductor / código aberto	3.18	Plataforma de software bioinformático.
limma	Projeto Bioconductor / código aberto	Versão não especificada no manuscrito	Pacote de análise de expressão gênica diferencial.
Procedimento de Benjamini–Hochberg	Método estatístico	Não aplicável (procedimento estatístico)	Método de ajuste da taxa de falsas descobertas.
NVIDIA RTX 3080	NVIDIA Corporation	RTX 3080; ≥8 GB de VRAM; versão CUDA/driver não especificada no manuscrito	Unidade de processamento gráfico com pelo menos 8 GB de memória de vídeo.
GPU compatível com CUDA	NVIDIA Corporation	Versão do toolkit CUDA não especificada no manuscrito; ≥8 GB de VRAM	Unidade de processamento gráfico com suporte à computação paralela de propósito geral.
Ubuntu Linux 22.04 LTS	Canonical Ltd. / código aberto	22.04 LTS	Sistema operacional Linux de 64 bits.
TIP3P	Desenvolvedores do campo de forças CHARMM / código aberto	TIP3P; nenhuma versão de software aplicável	Modelo explícito de água com três sítios.
MM/PBSA	Desenvolvedores do gmx_MMPBSA / código aberto	gmx_MMPBSA 1.5.x	Método de energia de ligação por mecânica molecular/Poisson–Boltzmann com área de superfície.