Вычислительный рабочий процесс для приоритизации генов хозяина, связанных с микробными метаболитами, при синдроме раздраженного кишечника с преобладанием запоров

Краткое содержание

Данный протокол объединяет прогнозирование мишеней микробных метаболитов, транскриптомику слизистой оболочки прямой кишки, анализ белок-белковых взаимодействий и обогащения путей, молекулярный докинг, моделирование молекулярной динамики и оценку свободной энергии связывания методом молекулярной механики/поверхности Пуассона — Больцмана (MM-PBSA) для формирования ранжированного списка гипотетических генов-кандидатов хозяина, связанных с метаболитами, и структурно приоритетных белок-лигандных комплексов для последующей экспериментальной проверки.

Аннотация

На данный момент не существует стандартизированного вычислительного конвейера для систематического определения приоритетных генов хозяина, ассоциированных с микробными метаболитами, и белок-лигандных комплексов на основе общедоступных химических, геномных и структурных баз данных. В данной статье описывается восьмиэтапный рабочий процесс, который принимает заданный пользователем набор метаболитов, производных микробиоты кишечника, и формирует ранжированный сокращенный список потенциальных ассоциированных с метаболитами генов хозяина, обогащенных биологических путей и структурно приоритетных белок-лигандных комплексов для последующих экспериментальных исследований. Конвейер объединяет в себе (i) хемоинформационное профилирование метаболитов; (ii) предсказание потенциальных мишеней с использованием нескольких баз данных с помощью инструментов предсказания мишеней на основе белково-химического взаимодействия и лигандов, а также программы молекулярного докинга; (iii) анализ дифференциальной экспрессии генов на основе общедоступных транскриптомных данных; (iv) поиск пересечений между мишенями и дифференциально экспрессируемыми генами; (v) построение сети белок-белковых взаимодействий и анализ обогащения путей; (vi) молекулярный докинг с помощью программы молекулярного докинга; (vii) молекулярно-динамическое моделирование продолжительностью 200 ns с использованием движка молекулярной динамики и силового поля белков, применяемого для моделирования молекулярной динамики; и (viii) оценку свободной энергии связывания методом MM-PBSA. В качестве практического примера девять производных или модифицированных микробиотой кишечника метаболитов, представляющих короткоцепочечные жирные кислоты, желчные кислоты, метаболиты триптофана и уролитин A, были обработаны с использованием общедоступного набора транскриптомных данных слизистой оболочки прямой кишки при IBS-C GSE36701. В результате рабочего процесса были ранжированы 17 уникальных предсказанных ассоциированных с метаболитами генов, которые дифференциально экспрессировались в данном наборе данных. Анализы докинга, молекулярно-динамического моделирования и MM-PBSA позволили структурно приоритезировать пять комплексов метаболит-белок: литохолевая кислота-VDR, литохолевая кислота-NR1H4/FXR, урсодезоксихолевая кислота-NR1H4/FXR, триптамин-HTR2A (моделирование проводилось в явном липидном бислое 1-Palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC)) и уролитин A-CASP3. Данный протокол разработан с возможностью адаптации к другим наборам метаболитов, транскриптомным наборам данных заболеваний и классам мишеней; все результаты представляют собой генерирующие гипотезы вычислительные предсказания, которые требуют независимой транскриптомной репликации, валидации на уровне белка и функциональных анализов ответа на лиганд, прежде чем можно будет сделать причинно-следственные или терапевтические выводы.

Введение

Синдром раздраженного кишечника с запорами (СРК-З) представляет собой распространенное функциональное желудочно-кишечное расстройство, характеризующееся рецидивирующими болями в животе, изменением характера дефекации, вздутием и запорами; по оценкам, глобальная распространенность этого заболевания составляет примерно 10–15% от общей численности населения^1,2. Современные методы фармакологической терапии, включая секретагоги, прокинетики и спазмолитики, могут облегчить отдельные симптомы у части пациентов, однако ответ на лечение остается неоднородным, а стойкая ремиссия достигается редко, что отражает сложную многофакторную патобиологию данного состояния^1,3,4. Таким образом, для формулирования проверяемых гипотез о новых терапевтических мишенях требуется более полное механистическое понимание того, как сигналы микробиоты кишечника передаются на слизистом уровне.

Микробиота кишечника способствует поддержанию гомеостаза в нижних отделах желудочно-кишечного тракта посредством продукции и биотрансформации химически разнообразных метаболитов, включая короткоцепочечные жирные кислоты (КЦЖК), вторичные желчные кислоты, соединения, производные триптофана, и метаболиты полифенолов, такие как уролитины^5,6,7,8. Эти молекулы взаимодействуют с клетками хозяина через широкий и не полностью охарактеризованный набор молекулярных мишеней, который выходит далеко за рамки канонических мембранных рецепторов, распознающих метаболиты, и включает ядерные рецепторы, цитозольные ферменты, белки, модифицирующие гистоны, предшественники пептидных гормонов и внутриклеточные сигнальные белки⁹. Изменения в составе микробного сообщества кишечника и профилях метаболитов были зафиксированы у пациентов с синдромом раздраженного кишечника (СРК), что дает биологическое обоснование для исследования того, нарушена ли транскрипция генов хозяина, связанных с реакцией на микробные метаболиты, в слизистой оболочке прямой кишки при СРК с запорами (СРК-З)¹⁰.

Панель из девяти метаболитов была определена априори для создания компактного, химически разнообразного и биологически интерпретируемого набора малых молекул, производных от микробиоты кишечника или модифицированных ею. Отбор основывался на пяти критериях: представленность основных классов микробных метаболитов, участвующих в сигнализации между хозяином и микробиотой; доказанное или вероятное воздействие на слизистую оболочку дистального отдела кишечника; наличие однозначных идентификаторов PubChem и канонических структур; размер молекулы и структурная доступность для предсказания мишеней на основе лигандов и докинга; а также предварительная вероятность участия в эпителиальной, нейроиммунной, энтероэндокринной, ядерно-рецепторной сигнализации или процессах, связанных с моторикой при СРК с запорами (IBS-C). Выбранная панель включала бутират и пропионат в качестве короткоцепочечных жирных кислот (SCFA); хенодезоксихолевую кислоту, литохолевую кислоту и урсодезоксихолевую кислоту в качестве желчных кислот; триптамин, индол-3-пропионовую кислоту и индол-3-молочную кислоту в качестве метаболитов, производных триптофана; а также уролитин A в качестве полифенольного метаболита, производного микробиоты кишечника^5,6,7,8,9,10.

Большинство предыдущих вычислительных и экспериментальных исследований были сосредоточены на изучении отдельных пар метаболит–рецептор или метаболит–фермент в изоляции; такой подход не позволяет охватить распределенный, конвергентный характер сигналинга микробных метаболитов в путях хозяина^9,11. Интеграция нескольких аналитических этапов обеспечивает взаимно усиливающую фильтрующую способность, которую не может обеспечить ни один отдельный этап самостоятельно. Вычислительное предсказание мишеней с использованием курируемых баз данных позволяет получить широкий набор белков-кандидатов хозяина для каждого метаболита. Пересечение с релевантными заболеванию транскриптомными данными существенно фильтрует этот набор, оставляя только тех кандидатов, чьи транскрипты изменены в контексте заболевания. Затем анализ обогащения путей и сети белок-белковых взаимодействий позволяют соотнести сокращенный список кандидатов с известными биологическими модулями. Молекулярный докинг обеспечивает первоначальную вычислительную оценку комплементарности связывающего кармана для каждого кандидатного комплекса, а дополнительная симуляция молекулярной динамики (MD) продолжительностью 200 ns с разложением свободной энергии связывания по методу MM-PBSA добавляет временное разрешение и термодинамическое измерение к структурной приоритизации, что недоступно при использовании только показателей докинга. Выполнение каждого шага независимо, без систематической интеграции и последовательной фильтрации, привело бы к созданию списков кандидатов, слишком обширных для экспериментальной проверки, и не позволило бы выявить конвергентную архитектуру путей.

В рамках данного протокола под термином «ген, связанный с метаболитом» (metabolite-associated gene, MAG) понимается ген человека, белковый продукт которого был определен как предполагаемая молекулярная мишень одного или нескольких метаболитов, производных микробиоты кишечника, по крайней мере одной курируемой базой данных вычислительного прогнозирования, и транскрипт которого дифференциально экспрессируется в соответствующем заболеванию транскриптомном наборе данных, используемом для демонстрации рабочего процесса. Это операциональное определение намеренно включает мембранные и ядерные рецепторы, цитозольные ферменты, сигнальные белки, предшественники пептидных гормонов и другие внутриклеточные белки. Обозначение MAG не является экспериментальным доказательством того, что метаболит связывается, образует белок-лигандный комплекс, активирует рецептор, изменяет количество белка или вызывает заболевание; это полученная вычислительным путем гипотеза, требующая экспериментального подтверждения.

В данном протоколе подробно описан полный восьмиэтапный вычислительный рабочий процесс (Рисунок 1), содержащий достаточно операционных деталей для обеспечения независимого воспроизведения, адаптации к другим панелям метаболитов или наборам данных о заболеваниях, а также расширения на другие контексты взаимодействия между хозяином и микробиотой. Данный рабочий процесс четко определен как основа для генерации гипотез и структурной приоритизации, которая оперирует исключительно общедоступными омиксными и структурными ресурсами и не позволяет делать выводы об изменении концентрации метаболитов, состояниях активации рецепторов, изменениях экспрессии белков, активности нижележащих сигнальных путей или клинической значимости только на основании результатов вычислительного анализа. В данной работе мы демонстрируем протокол на конкретном примере, используя девять метаболитов, производных от микробиоты кишечника или модифицированных ею, и общедоступный транскриптомный набор данных слизистой оболочки прямой кишки при IBS-C GSE36701, с целью идентификации MAG и приоритизации комплексов метаболит–белок для последующего экспериментального подтверждения.

Протокол

При проведении анализа использовались только общедоступные деидентифицированные транскриптомные данные из GSE36701, а также общедоступные химические, белковые и структурные базы данных. Доступ к базам данных осуществлялся в период с января по май 2026 года. Любая более поздняя дата доступа была задокументирована в отдельной таблице материалов (Table of Materials).

1. Дизайн исследования, требования к оборудованию и программному обеспечению

1. Определите рабочий процесс перед началом анализа. Используйте восемь этапов: выбор метаболитов, предсказание мишеней, анализ дифференциальной экспрессии, перекрытие мишеней и дифференциально экспрессируемых генов (DEG), анализ белок-белковых взаимодействий (PPI) / обогащение путей, молекулярный докинг, МД-моделирование и расчет MM-PBSA.
2. Отметьте, что докинг, МД и MM-PBSA являются исключительно инструментами структурной приоритизации. Не интерпретируйте эти результаты как экспериментальные доказательства связывания, активации рецептора, изменения содержания белка, терапевтической эффективности или причинно-следственной связи с заболеванием.
3. Проверьте компьютерное оборудование перед запуском МД-моделирования. Используйте 64-разрядную ОС Linux, 6-ядерный процессор или лучше, платформу с GPU-ускорением с объемом VRAM ≥8 GB или эквивалентную платформу с GPU-ускорением с объемом VRAM ≥8 GB, минимум 32 GB оперативной памяти и не менее 200 GB свободного места на диске для каждой МД-системы.
4. Зафиксируйте основное программное обеспечение: движок молекулярной динамики, программу для молекулярного докинга метаболитов-лигандов к целевым белкам, средство конвертации химических форматов файлов¹², программу генерации трехмерных структур лигандов, инструмент подготовки лигандов, инструментарий по подготовке входных данных для докинга, программу для молекулярного докинга метаболитов-лигандов к целевым белкам, среду программирования общего назначения, среду статистических вычислений с биоинформатическим программным каркасом и пакет для анализа дифференциальной экспрессии генов.
5. Зафиксируйте инструменты структурного анализа: веб-инструмент для построения мембранных систем, сервис параметризации лигандов, совместимый с CHARMM, инструмент расчета энергии связывания методом молекулярной механики/континуальной модели растворителя, библиотеку конвертации молекулярной топологии и параметров, программу трехмерной молекулярной визуализации, а также инструмент молекулярной визуализации и построения двухмерных диаграмм взаимодействий 2021 (см. Таблицу материалов для получения ссылок на скачивание и информации о версиях).
6. Запишите в Таблицу материалов/файл окружения точные идентификаторы версий использованного силового поля белков для МД-моделирования, CGenFF, CHARMM-GUI, программного каркаса R/биоинформатики и инструмента расчета энергии связывания методом молекулярной механики/континуальной модели растворителя. Отсутствующие идентификаторы помечайте как «невосстановимые» (not recoverable); не пытайтесь их угадать.

2. Выбор метаболитов и хемоинформатическая характеристика

1. Определите панель метаболитов перед прогнозированием мишеней. Включите бутират (PubChem CID: 264), пропионат (CID: 1032), хенодезоксихолевую кислоту (CID: 10133), литохолевую кислоту (CID: 9903), урсодезоксихолевую кислоту (CID: 31401), триптамин (CID: 1150), индол-3-пропионовую кислоту (CID: 3744), индол-3-молочную кислоту (CID: 92904) и уролитин A (CID: 5488186).
2. Получите канонические строки SMILES (Simplified Molecular Input Line Entry System) и идентификаторы PubChem CID для каждого метаболита. Проверьте синонимы и дубликаты структур перед прогнозированием мишеней. Сохраните итоговые идентификаторы в основном листе метаболитов.
3. Введите канонические строки SMILES в веб-инструмент для предсказания физико-химических свойств и ADME¹³ (см. Таблицу материалов). Зафиксируйте молекулярную массу, топологическую полярную площадь поверхности (TPSA), консенсусный logP, количество доноров и акцепторов водородных связей, количество вращающихся связей, прогнозируемую абсорбцию в желудочно-кишечном тракте, прогноз по P-гликопротеину, а также соответствие правилам Липински, Вебера, Гозе, Игана, Мюгге и наличие PAINS-оповещений.
4. Оставьте метаболиты, для которых структура была успешно распознана, молекулярная масса составляет ≤500 Da и отсутствуют PAINS-оповещения. Зафиксируйте любое несоответствие критериям и решение об удержании или исключении метаболита.
5. Определите состояния ионизации перед последующим прогнозированием мишеней и докингом. Используйте депротонированные карбоксилаты для бутирата и пропионата, нейтральные формы карбоновых кислот для желчных кислот, протонированный аммоний для триптамина и нейтральные формы для остальных метаболитов.

3. Прогнозирование потенциальных мишеней человека

1. Откройте инструмент для прогнозирования мишеней химического взаимодействия с белками¹⁴ (см. Таблицу материалов). Введите название каждого метаболита или PubChem CID, выберите Homo sapiens (ID таксономии: 9606) и установите минимальное комбинированное значение показателя взаимодействия ≥0.700.
2. При прогнозировании мишеней химического взаимодействия с белками отдайте приоритет каналам данных из экспериментальных и курируемых баз данных. Загрузите полную таблицу ассоциаций с белками для каждого метаболита.
3. Откройте программу для молекулярного докинга¹⁵ (см. Таблицу материалов). Введите каждую каноническую строку SMILES, выбрав Homo sapiens, и оставьте мишени с вероятностью ≥0.70.
4. Объедините результаты прогнозирования мишеней химического взаимодействия с белками и программы молекулярного докинга в виде объединения множеств для каждого метаболита. Оставьте все мишени, соответствующие порогу любой из баз данных, и удалите повторяющиеся записи символов генов.
5. Стандартизируйте записи белков в соответствии с утвержденными символами генов Комитета по номенклатуре генов человека (HGNC), используя сопоставление идентификаторов белков со стандартизированными символами HGNC или интегрированную базу данных информации о генах человека (см. Таблицу материалов). Приведите псевдонимы, устаревшие символы и аннотации изоформ к одному символу гена для каждого белка.
6. Классифицируйте каждую мишень как мембранный рецептор, ядерный рецептор, фермент, внутриклеточный сигнальный белок, пептидный гормон, гормон-связывающий белок или другой внутриклеточный белок. Занесите класс в таблицу мишеней.

4. Транскриптомный набор данных и анализ дифференциальной экспрессии генов

1. Получите доступ к набору данных GSE36701 с помощью веб-инструмента анализа дифференциальной экспрессии генов NCBI^16,17 (см. Таблицу материалов). Отметьте, что набор данных включает данные экспрессии биоптатов слизистой оболочки прямой кишки из групп пациентов с СРК с преобладанием запоров (СРК-З), СРК с преобладанием диареи (СРК-Д), постинфекционным СРК и группы здоровых добровольцев¹⁸.
2. Выполните поиск независимой валидационной когорты в базах GEO и ArrayExpress. Используйте комбинации терминов IBS-C, constipation-predominant irritable bowel syndrome, rectal mucosa, colonic mucosa, biopsy, transcriptome, microarray и RNA-seq. Зафиксируйте репозитории, поисковые запросы, дату поиска и информацию о том, был ли найден сопоставимый валидационный набор данных.
3. Запустите веб-инструмент анализа дифференциальной экспрессии генов из записи GSE36701 (см. Таблицу материалов). Назначьте 18 образцов СРК-З группе СРК-З, 40 здоровых добровольцев — контрольной группе, а образцы СРК-Д и постинфекционного СРК оставьте нераспределенными.
4. Проведите анализ дифференциальной экспрессии с использованием программного пакета для анализа дифференциальной экспрессии генов с коррекцией ложноположительных результатов (FDR) по методу Бенджамини-Хохберга¹⁹. Загрузите полную таблицу результатов, содержащую ID зонда, символ гена, название гена, logFC, AveExpr, модерированную t-статистику, исходное значение P-value и скорректированное значение P-value.
5. Объедините данные зондов на уровне генов. Удалите зонды, не имеющие символов генов; при наличии дубликатов символов оставьте зонд с наименьшим значением FDR; в случае равенства значений FDR используйте большее абсолютное значение logFC.

5. Анализ перекрытия мишеней с дифференциально экспрессируемыми генами и статистическая оценка

1. Пересеките каждый список предсказанных мишеней для конкретного метаболита со списком дифференциально экспрессируемых генов (DEG) на уровне генов при FDR < 0.05. Зафиксируйте пересекающиеся гены, исходный метаболит, logFC, скорректированное P-значение и направление экспрессии.
2. Объедините списки пересечений по конкретным метаболитам в недублирующий список MAG. Подсчитайте общее количество предсказанных мишеней, пересечений по конкретным метаболитам и общее количество уникальных MAG.
3. Оцените согласованность направления экспрессии на уровне зондов для генов с несколькими зондами. Отметьте любой ген, зонды которого указывают на разное направление экспрессии.
4. Составьте таблицу сопряженности для точного теста Фишера, используя общее количество записей с коллапсом генов, общее число DEG, общее число уникальных предсказанных мишеней и наблюдаемые MAG. Рассчитайте одностороннее P-значение, отношение шансов и 95% доверительный интервал с помощью реализации точного теста Фишера.
5. Если фоновая частота DEG превышает 50%, интерпретируйте пересечение как описательное, а не как независимо подтвержденное обогащение. Рассматривайте единообразное снижение экспрессии как описательный направленный паттерн, если не проводится отдельный тест на направленность.

6. Анализ сети белок-белковых взаимодействий и обогащение сигнальных путей

1. Загрузите полный список уникальных MAG в программу для построения сети белок-белковых взаимодействий и анализа обогащения путей²⁰ (см. Таблицу материалов). Выберите Homo sapiens и установите минимальный показатель взаимодействия (interaction score) на уровне 0.700.
2. Экспортируйте объединенную сеть белок-белковых взаимодействий и анализа обогащения путей, а также полную таблицу взаимодействий. Если текстовый майнинг (text-mining) приводит к созданию искусственно плотной топологии, отключите этот метод и оставьте только экспериментальные данные, данные по коэкспрессии и базы данных.
3. Сгенерируйте подсети по классам метаболитов для MAG, связанных с КЦЖК (SCFA), желчными кислотами и триптамином/серотонином. Используйте те же настройки организма и уровни достоверности, что и при построении сети белок-белковых взаимодействий и анализа обогащения путей.
4. Проведите анализ обогащения путей по базам данных Kyoto Encyclopedia of Genes and Genomes (KEGG)²¹, Reactome²² и Gene Ontology (GO) Biological Process^23,24. Примените поправку Бенджамини-Хохберга (BH FDR) <0.05 и экспортируйте все таблицы обогащения.

7. Молекулярный докинг

1. Получите экспериментально определенные структуры рецепторов из Protein Data Bank RCSB PDB²⁵ Исследовательского коллаборатора по структурной биоинформатике (см. Таблицу материалов). Для пяти приоритетных белок-лигандных комплексов используйте VDR/1DB1, NR1H4/FXR/3DCT, CASP3/2DKO и HTR2A/6A93.
2. Подготовьте каждый рецептор, сохранив цепь A и удалив воду, кокристаллизованные лиганды, кофакторы, ионы и записи HETATM, не относящиеся к белку. Для 6A93 перед подготовкой рецептора удалите сегмент слияния с T4-лизоцимом.
3. Добавьте полярные атомы водорода, присвойте заряды Гастейгера и сохраните каждый рецептор в формате PDBQT с помощью инструментария для анализа молекулярных структур. Перед конвертацией в PDBQT проверьте состояния протонирования гистидинов в сайте связывания и задокументируйте выбранные состояния.
4. Сгенерируйте трехмерную структуру каждого лиганда в инструменте конвертации файлов химических структур. Проведите минимизацию энергии с использованием универсального силового поля (UFF) в течение 500 шагов, установите состояние ионизации для pH 7.4, присвойте заряды Гастейгера и сохраните в формате PDBQT.
5. Определите область докинга (док-бокс) размером 25 Å x 25 Å x 25 Å с центром в центре масс кокристаллизованного лиганда. Используйте следующие координаты центра: (10, 19, 33) для VDR, (137, 31, 78) для FXR, (37, 34, 32) для CASP3 и (12, −1, 61) для HTR2A.
6. Выполните молекулярный докинг лигандов-метаболитов к целевым белкам^26,27 со следующими параметрами: exhaustiveness = 8, seed = 42, num_modes = 9 и energy_range = 3 kcal/mol. Зафиксируйте наилучший показатель Vina score и значения среднеквадратичного отклонения (RMSD) для всех поз.
7. Выберите модус 1 для каждого приоритетного комплекса. Постройте двухмерные диаграммы взаимодействий лиганд-остатки с помощью инструментов молекулярной визуализации и построения 2D-диаграмм взаимодействий, а также создайте трехмерные виды рецептор-лиганд в программе трехмерной молекулярной визуализации.
8. Проведите контрольный редокинг для VDR/1DB1 и FXR/3DCT. Считайте настройки докинга рецептора приемлемыми, если RMSD по тяжелым атомам составляет <2.0 Å относительно кристаллографической позы.
9. Проведите контроль перекрестного докинга, выполнив докинг литохолевой кислоты (LCA) в CASP3 и триптамина в VDR. Сравните показатели для когнатных и некогнатных пар и зафиксируйте случаи, когда разница в показателях составляет <1.0 kcal/mol.

8. Моделирование молекулярной динамики

1. Сгенерируйте параметры лиганда с помощью сервиса параметризации лигандов, совместимого с CHARMM²⁸ (см. Таблицу материалов). Проверьте все показатели штрафов (penalty scores) и отметьте любой параметр со значением штрафа >50.
2. Преобразуйте стрим-файлы лигандов в файлы .itp и .prm, совместимые с движком молекулярной динамики, с помощью скрипта конвертации топологии силового поля. Объедините файлы топологии лиганда и белка для каждого комплекса.
3. Выполните перераспределение масс водорода с помощью библиотеки конвертации молекулярной топологии и параметров. Создайте водные топологии с использованием силового поля белка, применявшегося для моделирования молекулярной динамики²⁹, и явной трехточечной модели воды³⁰.
4. Проведите сольватацию водных комплексов в додекаэдрическом боксе с минимальным расстоянием от растворенного вещества до края бокса 1.2 nm. Нейтрализуйте системы и добавьте NaCl до концентрации 0.15 M.
5. Соберите мембранную систему триптамин-HTR2A с помощью веб-инструмента для конструирования мембранных систем^31,32,33 (см. Таблицу материалов). Используйте координаты рецептора, выровненные по базе данных ориентации мембранных белков³⁴ (см. Таблицу материалов), чистый бислой POPC, слои воды толщиной 22.5 Å и NaCl 0.15 M.
6. Проведите минимизацию энергии всех систем методом наискорейшего спуска до 50,000 шагов. Перед уравновешиванием подтвердите сходимость при Fmax <1000 kJmol^-1nm^-1.
7. Проведите уравновешивание водных систем в фазах с постоянным числом частиц, объемом и температурой (NVT), а также с постоянным числом частиц, давлением и температурой (NPT). Уравновешивание мембранной системы выполните с использованием шестиэтапного веб-процесса многостадийной подготовки и уравновешивания молекулярных систем с постепенным снятием ограничений.
8. Выполните продуктивное моделирование MD продолжительностью 200 ns для всех пяти комплексов. Используйте шаг интегрирования 4 fs с перераспределением масс водорода (HMR), термостат V-rescale при 310 K, баростат Парринелло-Рамана при 1 bar, электростатику по методу Particle mesh Ewald (PME)³⁵ и ограничения LINCS³⁶.
9. Проанализируйте финальные траектории с помощью утилит анализа траекторий молекулярной динамики. Рассчитайте RMSD остова, среднеквадратичные флуктуации (RMSF) атомов Cα, радиус гирации, площадь поверхности, доступную для растворителя (SASA), и водородные связи белок-лиганд, используя последние 150 ns в качестве основного окна анализа.

9. Оценка свободной энергии связывания методом MM-PBSA

1. Извлеките снимки траектории для анализа MM-PBSA. Используйте 2001 кадр для каждого комплекса в водной среде и 201 обработанный кадр для подсистемы HTR2A, встроенной в мембрану.
2. Запустите инструмент расчета энергии связывания методом молекулярной механики/континуального растворителя³⁷ с сольватацией по Пуассону-Больцману, внутренней диэлектрической проницаемостью = 1, внешней диэлектрической проницаемостью = 80, неполярной сольватацией на основе SASA и без поправки на энтропию. Укажите среднюю свободную энергию связывания и стандартное отклонение.
3. Выполните декомпозицию по остаткам для всех пяти комплексов. Укажите стабилизирующие и дестабилизирующие остатки с абсолютным вкладом ≥0.5 kcalmol⁻¹.

    

Результаты

Потенциальные мишени, связанные с метаболитами

Девять метаболитов привели к формированию гетерогенных предсказанных наборов мишеней в программе прогнозирования химических взаимодействий белок-лиганд и программе молекулярного докинга. Пропионат, триптамин, желчные кислоты и уролитин А дали несколько мишеней, имеющих известное значение для желудочно-кишечной сигнализации. Предсказанный ландшафт мишеней включал канонические мембранные рецепторы, ядерные рецепторы, внутриклеточные ферменты, сигнальные белки и белки, связанные с пептидными гормонами. Таким образом, последующие результаты описываются как метаболит-ассоциированные гены (MAG), а не только как результаты по рецепторам (Таблица 1).

Сравнение с опубликованными данными о взаимодействиях метаболитов и белков

Для сопоставления результатов прогнозирования мишеней с имеющимися экспериментальными данными предсказанные связи между метаболитами и мишенями были разделены на три уровня доказательности: (i) экспериментально подтвержденные прямые или близкие взаимодействия на уровне классов метаболит-белок, при которых было сообщено, что данный метаболит или тесно связанный с ним эндогенный метаболит связывает, активирует, ингибирует или функционально регулирует кодируемый белок; (ii) взаимодействия, подтвержденные на уровне путей или классов мишеней, при которых предсказанная мишень принадлежит к установленному метаболит-зависимому пути или семейству рецепторов, но прямые доказательства для конкретной пары метаболит-белок ограничены; и (iii) исключительно вычислительные ассоциации, для которых в просмотренной литературе не было выявлено прямых экспериментальных взаимодействий. Данное сопоставление использовалось для контекстуализации, а не для валидации предсказанных MAGs.

Несколько прогнозов подтвердили ранее описанные биологические закономерности. Взаимодействие пропионат-FFAR2 было расценено как подтвержденное экспериментально, поскольку FFAR2/GPR43 является каноническим рецептором короткоцепочечных жирных кислот. Взаимодействие бутират-HDAC3 было классифицировано как подтвержденное экспериментально или на уровне класса, так как бутират является признанным ингибитором гистондеацетилазы, а прогнозируемое перекрытие касалось представителя семейства HDAC. Прогнозы, связанные с желчными кислотами и включающие NR1H4/FXR и VDR, были признаны подтвержденными установленными данными о биологии ядерных рецепторов желчных кислот, особенно для гидрофобных желчных кислот, таких как LCA; прогнозы для FXR, связанные с урсодезоксихолевой кислотой (UDCA), интерпретировались с осторожностью, поскольку UDCA, как правило, является более слабым или контекстно-зависимым лигандом FXR. Прогнозы для HTR1B, HTR2A, HTR2B и HTR6, связанные с триптамином, были классифицированы как подтвержденные серотонинергическим путем, а не как подтвержденные прямые рецептор-специфические взаимодействия, поскольку триптамин представляет собой микробиомный моноамин, производное триптофана, а серотониновые рецепторы являются установленными регуляторами моторики и секреции желудочно-кишечного тракта. Взаимодействие уролитин A-CASP3 считалось подтвержденным на уровне пути на основании опубликованных данных о связях между уролитином A и апоптотическими/каспаза-зависимыми ответами, но не прямыми доказательствами связывания с CASP3. Прогнозы для индол-3-молочной кислоты-KYAT1 и индол-3-пропионовой кислоты-KYAT1 были сохранены как чисто вычислительные гипотезы, поскольку в более широкой литературе приводятся данные о передаче сигналов хозяину посредством микробных производных индола, но нет данных о прямом связывании этих конкретных метаболитов с KYAT1^7,8,38,39,40.

Соответственно, в Таблице 1 представлено разграничение между вычислительным определением потенциальных мишеней и уровнем имеющихся экспериментальных данных или данных о путях. Также для каждой мишени указаны источник предсказания (предсказание мишени на основе химического взаимодействия с белком, программа молекулярного докинга или и то, и другое), комбинированный показатель взаимодействия для предсказания мишени на основе химического взаимодействия с белком и вероятность в программе молекулярного докинга, если мишень была определена с помощью данной программы. Предсказанные мишени, для которых нет прямых предварительных экспериментальных доказательств, описаны как гены-кандидаты, ассоциированные с метаболитами, которые требуют независимой валидации на уровне белка и ответа на лиганд.

Пересечение между предсказанными мишенями и дифференциально экспрессируемыми генами при IBS-C

Пересечение объединенного списка предсказанных мишеней и результатов дифференциальной экспрессии на уровне генов позволило выявить 17 уникальных предсказанных генов, связанных с метаболитами, экспрессия которых была значимо изменена при сравнении группы с IBS-C и группы здоровых добровольцев. Все 17 генов были подавлены. В этот набор вошли мембранные и ядерные рецепторы (CASR, FFAR2, GPR68, HTR1B, HTR2A, HTR2B, HTR6, NR1H4, TBXA2R, VDR) и нерецепторные белки (CASP3, GCG, GNAQ, GPHN, HDAC3, KYAT1, MLN) (Таблица 1, Рисунок 2A,B).

Все 17 MAG соответствовали порогу частоты ложноположительных результатов (FDR) ниже 0,05; 16 из 17 соответствовали более строгому критерию FDR < 0,001, при этом оставшийся ген (HTR1B) был значимым при FDR < 0,05. Семь из 17 мишеней (CASP3, GCG, GNAQ, GPHN, GPR68, HDAC3, TBXA2R) удовлетворяли обоим критериям: FDR < 0,001 и абсолютному значению логарифма кратности изменения log₂ более 1,0 (диапазон logFC от −1,34 до −1,10), что указывает на выраженное и стабильное снижение экспрессии для данной подгруппы. Остальные мишени демонстрировали умеренное, но статистически значимое снижение экспрессии (|logFC| в диапазоне от 0,45 до 0,97). Данная единообразная описательная закономерность интерпретировалась с осторожностью с учетом характеристик общегеномной экспрессии набора данных (см. статистическую оценку ниже).

Статистическая оценка перекрытия целевых DEG

Для формальной оценки статистической значимости перекрытия 17 генов был применен односторонний точный тест Фишера, где в качестве запрашиваемого набора использовались 17 предсказанных генов-мишеней, а в качестве геномного фона — все 18 296 уникальных записей с объединенными генами, обнаруженных в GSE36701. Из этого фона 17 296 генов (94,5%) были дифференциально экспрессированы при FDR < 0,05, что отражает почти универсальное транскрипционное подавление при сравнении слизистой оболочки прямой кишки при IBS-C. Все 17 предсказанных генов-мишеней оказались среди дифференциально экспрессированных генов (наблюдаемое перекрытие 17/17, 100%). Учитывая, что уровень дифференциальной экспрессии в фоновом наборе составляет 94,5%, ожидаемое перекрытие для любого случайно выбранного набора из 17 генов составляет 16,1 гена. Точный тест Фишера дал значение p = 0,384 с поправкой на непрерывность отношения шансов 2,03 (95% доверительный интервал 0,12–33,73), что не было статистически значимым при α = 0,05 (Рисунок 3A–C).

Этот результат указывает на то, что наблюдаемое перекрытие 17/17 не превышает перекрытия, ожидаемого случайно при таком профиле экспрессии по всему геному в данном наборе данных. Соответственно, эти результаты интерпретируются как описательная направленная закономерность, при которой все 17 предсказанных мишеней последовательно и значимо подвергались даунрегуляции в тканях слизистой оболочки прямой кишки при IBS-C, а не как доказательство статистического обогащения или независимой валидации на фоне всего генома. Формальное тестирование на обогащение потребовало бы воспроизведения в транскриптомных наборах данных с более селективными профилями дифференциальной экспрессии, в которых уровень значимости достигают существенно меньше половины всех генов. Следует подчеркнуть, что единообразная даунрегуляция всех 17 перекрывающихся генов является описательным наблюдением, а не отдельно валидированным статистическим результатом, поскольку фон дифференциальной экспрессии в данном наборе данных сам по себе преимущественно характеризуется даунрегуляцией; общее направление снижения экспрессии среди перекрывающихся генов было ожидаемым и не подвергалось формальному тесту на направленность. Следовательно, эта единообразная направленность не должна интерпретироваться как независимое статистическое доказательство скоординированной, метаболит-специфичной регуляции.

Специфические паттерны метаболитов

Пропионат имел наибольшее количество перекрывающихся генов, включая CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN и TBXA2R, что указывает на возможную роль сигнальных путей, чувствительных к короткоцепочечным жирным кислотам и связанным с Gq. Бутират перекрывался с HDAC3, что согласуется с биологией гистон-деацетилаз, связанной с бутиратом, хотя одно лишь снижение экспрессии мРНК не доказывает изменение чувствительности к бутирату. Перекрывающиеся гены, связанные с желчными кислотами, включали ядерные рецепторы VDR и NR1H4, оба из которых являются признанными эффекторами передачи сигналов желчных кислот в кишечнике^38,39. Триптамин перекрывался с HTR1B, HTR2A, HTR2B и HTR6, что указывает на серотонинергическую передачу сигналов как на потенциальный модуль — систему с хорошо изученной ролью в моторике и секреции желудочно-кишечного тракта⁴⁰. Индол-3-молочная кислота и индол-3-пропионовая кислота перекрывались с KYAT1, а уролитин A — с CASP3.

Обогащение путей

Функциональный анализ обогащения 17 перекрывающихся генов выявил пути, связанные с нижележащим сигналингом G-белоксопряженных рецепторов (GPCR), сигналингом Gαq, связыванием лигандов GPCR, серотонинергическим синапсом, взаимодействием нейроактивных лигандов с рецепторами, трансдукцией кальциевого сигнала, сигналингом cAMP и секрецией пептидных гормонов. Эти результаты согласуются с составом набора генов и подтверждают его биологическую согласованность, однако они отражают функциональную аннотацию представленных генов, а не являются независимым доказательством активности на уровне путей.

Структура сети белок-белковых взаимодействий

Построение сети белок-белковых взаимодействий и анализ обогащения путей были интерпретированы с использованием трех дополняющих друг друга сетей. В комбинированной мета-сети из 17 генов (Сеть 1) наиболее выраженной структурой, подтвержденной аннотацией, был компонент сигналинга GPCR/Gαq с центром в GNAQ, связывающий GNAQ с рецептор-ассоциированными генами, включая TBXA2R, CASR, HTR2A и HTR2B. Также сохранилась ограниченная связность серотониновых рецепторов, наиболее заметная между HTR2A и HTR2B, в то время как несколько других генов остались изолированными или слабо связанными при выбранном пороге достоверности. Сеть, специфичная для пропионата (Сеть 2), продемонстрировала более ограниченную топологию: GNAQ сохранил подтвержденные аннотацией связи с CASR и TBXA2R, тогда как FFAR2, GPR68, GCG, GPHN и MLN были изолированы или слабо связаны. Сеть триптамина/серотонина (Сеть 3) включала HTR1B, HTR2A, HTR2B и HTR6; в этом подмножестве HTR2A и HTR2B показали основную связь, подтвержденную аннотацией, в то время как HTR1B и HTR6 не были напрямую связаны при выбранном пороге (Рисунок 4A–C).

Молекулярный докинг

Молекулярный докинг был проведен для пяти выбранных комплексов метаболит-белок. Пары «желчная кислота — ядерный рецептор» продемонстрировали более благоприятные значения по шкале Vina, чем пары уролитин A-CASP3 и триптамин-HTR2A. Наилучший показатель был получен для LCA-VDR (−10.0 kcal/mol), за которым следовали LCA-NR1H4/FXR (−9.9 kcal/mol) и UDCA-NR1H4/FXR (−9.4 kcal/mol). Уролитин A-CASP3 и триптамин-HTR2A имели более низкие, но все же приемлемые показатели, составившие −7.1 kcal/mol (Таблица 2).

Для комплекса LCA-VDR (PDB ID: 1DB1) прогнозируемая поза была подтверждена наличием обычной водородной связи между кислородом карбоксилата LCA и Ser278 (4.29 Å), а также обширными гидрофобными контактами с Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233 и His397, и дополнительными контактами Ван-дер-Ваальса с Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275 и Phe150. Поза с наивысшим рейтингом имела показатель Vina −10.0 kcal/mol, размер полости 2055 ų и центр сетки (10, 19, 33) (Таблица 3, Рисунок 5A,B).

Для комплекса LCA-NR1H4/FXR (PDB ID: 3DCT) значение скоринга докинга −9.9 kcal/mol сопровождалось предсказанными водородными связями с участием His294 и Ile335, π-Sigma взаимодействием с His294, а также гидрофобными Alkyl или π-Alkyl контактами с Met290, Met328, Ala291, Leu287, Ile352 и His447; дополнительные ван-дер-ваальсовы контакты подтверждали размещение стероидного остова в кармане FXR (Таблица 4, Рисунок 6A,B).

Предполагаемая конформация комплекса UDCA-NR1H4/FXR (PDB ID: 3DCT) продемонстрировала наличие стандартной водородной связи с His447 (3,66 Å), еще одной водородной связи с Gly322 (3,46 Å), π-анионного взаимодействия с Val325 (4,96 Å) и углерод-водородной связи с Trp469 (4,51 Å). Карта взаимодействий также выявила неблагоприятные донор-донорные контакты с Arg395 (3,89 Å) и Gln396 (3,40 Å), что позволяет предположить, что более низкий показатель Vina для UDCA по сравнению с LCA в одном и том же кармане рецептора может быть обусловлен менее благоприятной локальной геометрией или электростатикой (Таблица 5, Рисунок 7A,B).

В комплексе уролитин A-CASP3 (PDB ID: 2DKO) прогнозируемый способ связывания характеризовался образованием обычных водородных связей с Gln161 (3.78 и 4.19 Å), Ser120 (3.95 Å) и Arg207 (3.05 и 3.77 Å), а также дополнительно стабилизировался за счет π-катионных взаимодействий с Arg207, π-донорной водородной связи с Cys163 и дополнительных π-алкильных контактов и контактов Ван-дер-Ваальса с участием Arg64, Ala162, His121, Ser205 и Trp206 (Таблица 6, Рисунок 8A,B).

Для комплекса триптамин-HTR2A (PDB ID: 6A93) предсказанная поза была стабилизирована электростатическим солевым мостиком между протонированным амином триптамина и Asp155 — консервативным аспартатом 3-й трансмембранной спирали (D3.32 по нумерации Ballesteros-Weinstein), который закрепляет протонированный амин аминергических лигандов в серотониновых и родственных им рецепторах^41,42,43, а также водородными связями с Thr160 и Ser159, ароматическими контактами с Phe340 и Trp336 и π-алкильными взаимодействиями с Val156 и Ile163. Дополнительные ван-дер-ваальсовы контакты с Tyr370, Phe339, Ser242, Phe243, Phe332 и Leu123 подтверждали паттерн связывания в ортостерическом кармане (Таблица 7, Рисунок 9A,B).

Валидация протокола докинга

Для оценки надежности протокола докинга были проведены два дополнительных контрольных эксперимента. Для положительного контроля (редокинг) кокристаллизованные лиганды были извлечены из их референтных рентгеноструктурных моделей и повторно задокированы в соответствующие нативные сайты связывания. Самая высокоранжированная предсказанная поза для аналога витамина D VDX в VDR/1DB1 отклонилась на 0.87 Å от кристаллографического положения, а кокристаллизованный лиганд WAY-362450 в FXR/3DCT отклонился на 1.79 Å; оба значения оказались ниже общепринятого порога приемлемости 2.0 Å, что подтверждает геометрическую достоверность протокола докинга для данных рецепторных систем (Рисунок 10A,B). Для отрицательного контроля (кросс-докинг) литохолевая кислота была задокирована в каспазу-3 (2DKO), цистеиновую протеазу, для которой она не является известным лигандом; полученный предсказанный скор (−8.3 kcal/mol) был на 1.7 kcal/mol слабее, чем для ее когнатной мишени VDR (−10.0 kcal/mol), что согласуется с предсказанной селективностью сайта связывания. Докинг триптамина в VDR дал предсказанный скор −6.4 kcal/mol по сравнению с −7.1 kcal/mol для его когнатной мишени HTR2A; разница в 0.7 kcal/mol находится в пределах зарегистрированной погрешности скоринга молекулярного докинга метаболических лигандов к белкам-мишеням и, следовательно, указывает лишь на умеренную предсказанную селективность для этого меньшего лиганда (Рисунок 10C). В совокупности эти контрольные эксперименты показывают, что протокол докинга воспроизводит известные геометрии связывания и дифференцирует когнатные и некогнатные пары в протестированных условиях, оставаясь при этом вычислительным прогнозом, который не заменяет экспериментальные измерения аффинности (Таблица 8).

Моделирование молекулярной динамики

Молекулярно-динамическое моделирование пяти приоритетных комплексов проводилось с использованием продуктивных траекторий длительностью 200 ns. Четыре комплекса с растворимыми белками и ядерными рецепторами моделировались в явном водном растворителе, в то время как комплекс триптамина с HTR2A моделировался в явном липидном бислое POPC для обеспечения физиологически соответствующего мембранного окружения для этого рецептора, сопряженного с G-белком. Анализ позволил проверить динамическую стабильность задокканных поз в условиях зависимости от времени и сравнить относительное структурное поведение различных комплексов (Таблица 9).

Профиль RMSD комплекса LCA-VDR/1DB1 продемонстрировал короткий период уравновешивания в течение первых 10 ns, за которым следовало стабильное плато с колебаниями преимущественно в диапазоне 0,20–0,28 nm (Рисунок 11A). Значения RMSF были низкими, а колебания основной цепи составляли < 0,15 nm для большинства остатков (Рисунок 11B). Анализ водородных связей выявил устойчивую сеть из 2–5 водородных связей с эпизодическим увеличением их количества до 7 (Рисунок 11C). Радиус инерции (Rg) поддерживался в диапазоне 1,25–1,75 nm, а доступная растворителю площадь поверхности (SASA) оставалась на уровне около 130 nm² (Рисунок 11D,E).

Комплекс urolithin A-CASP3/2DKO продемонстрировал более высокую динамическую активность. Значение RMSD изначально увеличивалось, а затем колебалось в диапазоне от 0.4 до 0.7 nm, с кратковременным событием высокого отклонения примерно на 165 ns (Рисунок 12A). Анализ RMSF показал высокую мобильность на уровне остатков, при этом наибольшие флуктуации наблюдались в области гибкой петли вокруг остатка 175 (Рисунок 12B). Анализ водородных связей выявил наличие обширной сети из примерно 2–5 связей в первые 30–40 ns, после чего наблюдалось в основном от 0 до 2 перемежающихся связей (Рисунок 12C). Соответствующие профили радиуса гирации и SASA представлены на Рисунке 12D,E.

Для систем NR1H4/FXR (3DCT) с желчными кислотами профиль RMSD остова оставался в относительно узком диапазоне на протяжении большей части траектории (Рисунок 13A), в то время как профиль RMSF продемонстрировал более низкую подвижность в основных областях и более сильные флуктуации в гибких областях (Рисунок 13B). Комплекс LCA-3DCT поддерживал примерно три-четыре устойчивые водородные связи на протяжении всей траектории, тогда как комплекс UDCA-3DCT характеризовался более выраженными флуктуациями водородных связей и снижением их количества примерно после 125 ns. Профили радиуса гирации для систем с LCA и UDCA показаны на Рисунке 13C,D соответственно, а соответствующие профили SASA представлены на Рисунке 13E,F.

Молекулярная динамика мембранного комплекса триптамин-HTR2A

Моделирование комплекса триптамин-HTR2A/6A93 проводилось в течение 200 ns в явном липидном бислое POPC, состоящем из 258 молекул липидов, с использованием явной трехточечной модели воды и 0.15 M NaCl; общий размер системы составил примерно 100,925 атомов^33,44,45. Рецептор оставался стабильно встроенным в бислой на протяжении всей траектории (Рисунок 14). RMSD остова увеличилось примерно с 0.10 nm до стабильного плато в районе 0.15–0.20 nm в течение первых 100 ns и оставалось стабильным в дальнейшем, при этом все значения были ниже 0.25 nm, что указывает на сохранение рецептором стабильной конформации в мембранном окружении без глобального развертывания (Рисунок 15A). RMSF по остаткам продемонстрировала низкие флуктуации в трансмембранном спиральном ядре с ожидаемой более высокой мобильностью в петлевых и терминальных областях, что соответствует типичной гибкости GPCR (Рисунок 15B). Радиус гирации был строго ограничен диапазоном примерно от 2.06 до 2.12 nm, а SASA колебалась в узком диапазоне без прогрессирующего дрейфа, что подтверждает сохранение компактного трансмембранного пучка (Рисунок 15C,D).

Водородные связи между белком и лигандом сохранялись на протяжении всей траектории (Рисунок 15E), при этом количество водородных связей существенно колебалось в диапазоне от 1 до 3. Для специфической оценки стабильности ключевого ионного взаимодействия на протяжении всей траектории контролировалось минимальное расстояние между протонированным азотом аммония триптамина и атомами кислорода карбоксилатной группы Asp155 (D3.32). Это расстояние имело узкое распределение вокруг среднего значения 0.270 nm (минимум 0.247 nm, максимум 0.424 nm), и контакт солевого мостика (< 0.4 nm) поддерживался в течение 99.9% времени симуляции с лишь двумя краткими переходными отклонениями, без событий устойчивой диссоциации (Рисунок 16). Эти результаты позволяют предположить, что консервативное ионное взаимодействие Asp155 было достаточным для стабилизации триптамина в ортостерическом кармане HTR2A на протяжении всего моделирования мембраны.

Расчет свободной энергии связывания методом MM-PBSA и декомпозиция по остаткам

Для добавления дополнительного уровня энергетической приоритизации пяти комплексов был проведен анализ MM-PBSA (Таблица 10). Для четырех водных комплексов декомпозиция по остаткам позволила определить основные энергетические вклады для каждого предсказанного режима связывания. В комплексе LCA-VDR/1DB1 благоприятный вклад внесли лиганд и Gln317, в то время как Trp286 оказал неблагоприятное воздействие. В комплексе уролитин A-CASP3/2DKO остатки Arg64 и Arg207 продемонстрировали сильно отрицательные вклады по остаткам, что указывает на значительную полярную или электростатическую стабилизацию; тем не менее, соответствующая траектория оставалась высокодинамичной, что доказывает, что благоприятная энергетика на уровне отдельных остатков сама по себе не гарантирует устойчивой стабильности комплекса. Для систем 3DCT связывание LCA определялось преимущественно Arg331, тогда как связывание UDCA включало более распределенную энергетическую сеть, состоящую из Glu326, Asp394, Arg395, Arg441 и Asp470. В четырех водных системах декомпозиция MM-PBSA подтвердила относительную приоритизацию комплексов на основе LCA.

Для комплекса триптамина с HTR2A/6A93, встроенного в мембрану, анализ MM-PBSA проводился для подсистемы белок–лиганд, извлеченной из траектории бислоя^46,47. Благоприятный вклад был отмечен для лиганда и Asp155 (D3.32), который оказался наиболее значимым стабилизирующим остатком, что согласуется с взаимодействием посредством солевого мостика, выявленным как при докинге, так и при анализе расстояний по траектории. Trp137 продемонстрировал наибольший неблагоприятный вклад на уровне отдельного остатка среди окружающих остатков ортостерического кармана (Ser86, Phe87, Phe133, Phe140, Phe141, Val156, Ser159, Thr160, Ile163, Val167, Tyr171), которые в совокупности образуют сеть ароматических и полярных контактов, выстилающую связывающий карман. Данные значения представляют собой относительные расчетные оценки для структурной приоритизации и не являются экспериментальными значениями аффинности связывания.

Рисунок 1: Вычислительный алгоритм приоритизации генов хозяина, связанных с метаболитами, при IBS-C. Схематическое представление восьмиэтапного рабочего процесса, объединяющего отбор метаболитов, предсказание мишеней, анализ дифференциальной экспрессии транскриптома, анализ перекрытий, обогащение сетей и путей, молекулярный докинг, симуляцию молекулярной динамики и анализ свободной энергии связывания методом MM-PBSA. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 2: Анализ дифференциальной экспрессии и перекрытия мишеней метаболитов в слизистой оболочке при IBS-C. (A) График вулкана дифференциальной экспрессии генов в GSE36701. Синие точки — значимо подавленные гены; красные точки — значимо активированные гены; серые точки — гены без значимых изменений. Отмечены выбранные перекрывающиеся гены, связанные с метаболитами. (B) Диаграмма Венна, показывающая перекрытие между 330 уникальными предсказанными мишенями метаболитов и подавленными генами в GSE36701; общими оказались 17 генов. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 3: Статистическая оценка 17 предсказанных генов-мишеней метаболитов по данным GSE36701. (A) Значение log2 fold change для каждого из 17 генов, окрашенные по уровню значимости. (B) Частота дифференциальной экспрессии фоновых генов по сравнению с предсказанными мишенями, с применением точного теста Фишера. (C) Для точного теста Фишера используется таблица сопряженности 2x2. Все 17 мишеней были значимо подавлены; перекрытие интерпретируется как описательный направленный паттерн, а не как статистическое обогащение. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 4: Построение композитной сети белок-белковых взаимодействий и анализ обогащения путей в сетях белок-белковых взаимодействий перекрывающихся генов, связанных с метаболитами. (A) Сеть 1: комбинированная мета-сеть всех 17 генов. (B) Сеть 2: специфическая для пропионата сеть из восьми генов (CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN, TBXA2R). (C) Сеть 3: сеть триптамина/серотонина из четырех генов (HTR1B, HTR2A, HTR2B, HTR6). Сети были построены для Homo sapiens с использованием параметров построения сети белок-белковых взаимодействий и обогащения путей с уровнем достоверности ≥ 0.700. Ребра представляют функциональные связи, подтвержденные аннотациями. Нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 5: Трехмерное и двухмерное структурное представление литохолевой кислоты в комплексе с VDR (PDB ID: 1DB1). (A) Трехмерное представление поверхности и упрощенной структуры (cartoon), где литохолевая кислота показана в виде сфер. (B) Двумерная карта взаимодействий, демонстрирующая водородную связь с Ser278, а также окружающие гидрофобные контакты и взаимодействия Ван-дер-Ваальса. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 6. Трехмерное и двумерное структурное представление литохолевой кислоты в комплексе с NR1H4/FXR (PDB ID: 3DCT). (A) Трехмерное представление поверхности и ленточной модели. (B) Двумерная карта взаимодействий, показывающая водородные связи с His294 и Ile335, π-сигма-взаимодействие и окружающие контакты. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 7: Трехмерное и двумерное структурное представление урсодезоксихолевой кислоты в комплексе с NR1H4/FXR (PDB ID: 3DCT). (A) Трехмерное представление поверхности и ленточной модели. (B) Двумерная карта взаимодействий, показывающая водородные связи с His447 и Gly322, $\pi$-анионное взаимодействие с Val325, углерод-водородную связь с Trp469 и неблагоприятные донор-донорные контакты с Arg395 и Gln396. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 8: Трехмерное и двумерное структурное представление уролитина A в комплексе с CASP3 (PDB ID: 2DKO). (A) Трехмерное представление поверхности и ленточная модель. (B) Двумерная карта взаимодействий, демонстрирующая водородные связи с Gln161, Ser120 и Arg207, π-катионные взаимодействия с Arg207, π-донорную водородную связь с Cys163 и окружающие контакты. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 9: Трехмерное и двухмерное структурное представление триптамина в комплексе с HTR2A (PDB ID: 6A93). (A>) Трехмерное представление поверхности и структуры «лента», созданное в программе трехмерной молекулярной визуализации. (B>) Двумерная карта взаимодействий, созданная с помощью инструментов молекулярной визуализации и построения двухмерных диаграмм взаимодействий, иллюстрирующая солевой мостик Asp155 и дополнительные взаимодействия в сайте связывания. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 10: Валидация протокола докинга. (A,B) Редокинг кокристаллизованных лигандов в VDR/1DB1 (RMSD 0.87 Å) и FXR/3DCT (RMSD 1.79 Å); кристаллографические и редокированные позы наложены друг на друга, обе находятся ниже порога приемлемости 2.0 Å. (C) Селективность кросс-докинга: показатели Vina для когнатных и некогнатных связей литохолевой кислоты и триптамина. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 11. Анализ траектории молекулярной динамики комплекса LCA-VDR/1DB1 в течение 200 ns. (A) Профиль RMSD. (B) Профиль RMSF. (C) Количество водородных связей. (D) Профиль радиуса гирации. (E) Профиль SASA. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 12: Анализ траектории молекулярной динамики комплекса уролитин A-CASP3/2DKO в течение 200 ns. (A) Профиль RMSD, демонстрирующий широкие конформационные флуктуации и кратковременное событие с высоким отклонением около 165 ns. (B) Профиль RMSF, демонстрирующий выраженную гибкость на уровне остатков вблизи остатка 175. (C) Количество водородных связей. (D) Профиль радиуса гирации. (E) Профиль SASA. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 13: Анализ траекторий молекулярной динамики систем желчных кислот NR1H4/FXR (3DCT) в течение 200 ns. (A) Профиль RMSD остова для комплекса 3DCT. (B) Профиль RMSF остова. (C) Профиль радиуса гирации для 3DCT-LCA. (D) Профиль радиуса гирации для 3DCT-UDCA. (E) Профиль SASA для 3DCT-LCA. (F) Профиль SASA для 3DCT-UDCA. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 14: Комплекс триптамин-HTR2A, встроенный в эксплицитный липидный бислой POPC. Рецептор представлен в виде упрощенной схемы (cartoon), пересекающей бислой; липиды POPC показаны в виде линий с выделенными фосфатными головными группами; триптамин расположен внутри ортостерического кармана. Вода показана над и под мембраной. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 15: Анализ траектории молекулярной динамики комплекса триптамина-HTR2A/6A93 в течение 200 ns в эксплицитном липидном бислое POPC. (A) Профиль RMSD остова. (B) Поостаточный профиль RMSF. (C) Профиль радиуса гирации. (D) Профиль SASA. (E) Количество водородных связей между белком и лигандом. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 16: Сохранение ионного взаимодействия между триптамином и Asp155 (D3.32) на протяжении 200 ns мембранной траектории. На графике представлена зависимость минимального расстояния между атомом азота аммонийной группы триптамина и атомами кислорода карбоксилата Asp155 от времени; пунктирная линия отмечает порог контакта солевого мостика 0.4 nm. Контакт сохранялся в течение 99.9% времени симуляции. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Символ гена	Исходный метаболит(ы)	Функциональная категория	log2FC	FDR (скорректированное значение P)	Уровень значимости
GCG	Пропионат	Белок, связанный с пептидными гормонами	−1.342	1.97e−7	FDR <0.001 & |logFC > 1
HDAC3	Бутират	Фермент	−1.234	2.44e−6	FDR <0.001 & |logFC| > 1
CASP3	Уролитинин A	Фермент	−1.198	6.66e−7	FDR <0.001 & |logFC| > 1
GPR68	Пропионат	Мембранный рецептор	−1.137	4.35e−6	FDR <0.001 & |logFC| > 1
GNAQ	Пропионат	Белок внутриклеточного сигналинга	−1.122	1.05e−6	FDR <0.001 & |logFC| > 1
GPHN	Пропионат	Другой внутриклеточный белок	−1.109	1.13e−6	FDR <0.001 & |logFC| > 1
TBXA2R	Пропионат	Мембранный рецептор	−1.104	4.04e−7	FDR <0.001 & |logFC| > 1
HTR6	Триптамин	Мембранный рецептор	−0.967	2.17e−5	FDR <0.001
VDR	Литохолевая кислота	Ядерный рецептор	−0.942	5.73e−7	FDR <0.001
HTR2A	Триптамин	Мембранный рецептор	−0.937	4.99e−6	FDR <0.001
FFAR2	Пропионат	Мембранный рецептор	−0.889	1.44e−4	FDR <0.001
NR1H4	Литохолевая кислота / Урсодезоксихолевая кислота	Ядерный рецептор	−0.861	3.68e−6	FDR <0.001
HTR2B	Триптамин	Мембранный рецептор	−0.702	1.29e−4	FDR <0.001
MLN	Пропионат	Белок, связанный с пептидными гормонами	−0.605	7.39e−5	FDR <0.001
KYAT1	Индол-3-молочная кислота / Индол-3-пропионовая кислота	Фермент	−0.530	3.61e−4	FDR <0.001
CASR	Пропионат	Мембранный рецептор	−0.483	4.05e−4	FDR <0.001
HTR1B	Триптамин	Мембранный рецептор	−0.455	3.18e−2	FDR <0.05

Таблица 1: Прогнозируемые гены-мишени, связанные с метаболитами, которые перекрываются с дифференциально экспрессируемыми генами в наборе данных слизистой оболочки прямой кишки при IBS-C. Все перечисленные перекрывающиеся гены были подавлены (downregulated). Таблица 1 представлена отдельно в виде электронной таблицы и содержит для каждой мишени: исходный метаболит(ы), функциональную категорию, источник прогнозирования мишени (прогноз взаимодействия химическое вещество-белок, программа молекулярного докинга или и то, и другое), комбинированный показатель взаимодействия прогноза взаимодействия химическое вещество-белок и, при наличии, вероятность программы молекулярного докинга, уровень прогнозирования, log2 fold change и FDR с уровнем значимости экспрессии. Источник: Значения экспрессии генов были получены из таблицы дифференциальной экспрессии GSE36701 с объединением по генам (зонд с наименьшим FDR на ген). Источник прогнозирования мишени и значения достоверности были скомпилированы из результатов прогноза взаимодействия химическое вещество-белок и вывода программы молекулярного докинга с использованием пороговых значений комбинированного показателя взаимодействия прогноза взаимодействия химическое вещество-белок ≥ 0.700 и вероятности программы молекулярного докинга ≥ 0.70. Показатели прогноза взаимодействия химическое вещество-белок представляют собой комбинированные показатели по шкале 0–1; STP обозначает вероятность программы молекулярного докинга. Уровень 1 = строгая поддержка прогноза взаимодействия химическое вещество-белок; Уровень 1+ = строгая поддержка прогноза взаимодействия химическое вещество-белок, подтвержденная программой молекулярного докинга.

Комплекс	Белок (PDB ID)	Лиганд	Балл Vina (ккал/моль)	Размер полости (ų)	Центр сетки X,Y,Z (А)	Поле поиска (A)
LCA-VDR	VDR (1DB1)	Литохолевая кислота	−10.0	2055	10, 19, 33	25 x 25 x 25
LCA-NR1H4/FXR	NR1H4/FXR (3DCT)	Литохолевая кислота	−9.9	3395	137, 31, 78	25 x 25 x 25
УДХК-NR1H4/FXR	NR1H4/FXR (3DCT)	Урсодезоксихолевая кислота	−9.4	3395	137, 31, 78	25 x 25 x 25
Уролитин А-CASP3	CASP3 (2DKO)	Уролитин А	−7.1	233	37, 34, 32	25 x 25 x 25
Триптамин-HTR2A	HTR2A (6A93)	Триптамин	−7.1	3238	12, −1, 61	25 x 25 x 25

Таблица 2: Результаты молекулярного докинга: наивысшие баллы молекулярного докинга лиганд-метаболитов с белками-мишенями и параметры полости для пяти приоритетных белково-лигандных комплексов. Размер полости указан в ų. Источник: Docking_Validation/Results/Docking_Validation_Results.xlsx, лист 'Original_Docking_Scores'. Молекулярный докинг лиганд-метаболитов с белками-мишенями; exhaustiveness = 8, seed = 42 (фиксированный), num_modes = 9 для всех комплексов; представлены позы с наивысшим рейтингом (мода 1).

Тип взаимодействия	Остаток(и)	Расстояние (А)	Примечания
Обычная водородная связь	Ser278	4.29	кислород карбоксильной группы LCA
Гидрофобный / Пи-алкильный контакт	Leu230, Val234, Trp286, Val300, His305, Tyr295, Leu233, His397	-
Контакт Ван-дер-Ваальса	Met272, Leu313, Ile271, Ile268, Leu309, Phe422, Val418, Ala231, Ala303, Cys288, Ser275, Phe150	-

Таблица 3: Режимы связывания, полученные при докинге литохолевой кислоты с VDR (PDB ID: 1DB1). Источник: инструменты молекулярной визуализации и построения двумерных диаграмм взаимодействий; двумерные диаграммы взаимодействия лиганда с остатками, как указано в разделе «Результаты (Молекулярный докинг)» рукописи. «-» означает, что значение расстояния для данного контакта не было указано отдельно.

Тип взаимодействия	Остаток(и)	Расстояние (A)	Примечания
Водородная связь	His294	-
Водородная связь	Ile335	-
$\pi$-$\sigma$-взаимодействие	His294	-
Алкильный / Пи-алкильный (гидрофобный)	Met290, Met328, Ala291, Leu287, Ile352, His447	-
Контакт Ван-дер-Ваальса	Дополнительные остатки кармана (не указаны индивидуально в источнике)	-	Обеспечивает размещение стероидного остова

Таблица 4: Способы связывания, полученные при докинге литохолевой кислоты с NR1H4/FXR (PDB ID: 3DCT).

Источник: инструмент молекулярной визуализации и построения двумерных диаграмм взаимодействий; двумерные диаграммы взаимодействий лиганд-остаток, представленные в разделе «Результаты» (Молекулярный докинг) рукописи. «-» означает, что значение расстояния для данного контакта не было указано отдельно.

Тип взаимодействия	Остаток(и)	Расстояние (A)	Примечания
Обычная водородная связь	His447	3.66
Водородная связь	Gly322	3.46
$\pi$-анионное взаимодействие	Val325	4.96
Углерод-водородная связь	Trp469	4.51
Неблагоприятный донор-донорный контакт	Arg395	3.89
Неблагоприятный донорно-донорский контакт	Gln396	3.40

Таблица 5: Способы связывания, полученные при докинге урсодезоксихолевой кислоты с NR1H4/FXR (PDB ID: 3DCT).Источник: инструмент молекулярной визуализации и двумерные диаграммы взаимодействия лиганд-остаток, представленные в разделе «Результаты» (Молекулярный докинг) рукописи. «-» означает, что значение расстояния для данного контакта не было указано отдельно.

Тип взаимодействия	Остаток(и)	Расстояние (А)	Примечания
Обычная водородная связь	Gln161	3.78
Обычная водородная связь	Gln161	4.19	вторичный контакт
Обычная водородная связь	Ser120	3.95
Обычная водородная связь	Arg207	3.05
Классическая водородная связь	Arg207	3.77	вторичный контакт
$\pi$-катионное взаимодействие	Arg207	-
Водородная связь с $\pi$-донором	Cys163	-
Пи-алкильный / ван-дер-ваальсов контакт	Arg64, Ala162, His121, Ser205, Trp206	-

Таблица 6: Режимы связывания, полученные при докинге уролитина А с CASP3 (PDB ID: 2DKO). Источник: инструмент молекулярной визуализации и двухмерные диаграммы взаимодействий лиганда с остатками, представленные в разделе «Результаты (Молекулярный докинг)» рукописи. «-» означает, что значение расстояния для данного контакта не было указано отдельно.

Тип взаимодействия	Остаток(и)	Расстояние (A)	Примечания
Электростатический солевой мостик	Asp155 (D3.32)	-	протонированный амин триптамина
Водородная связь	Thr160	-
Водородная связь	Ser159	-
Ароматический контакт	Phe340, Trp336	-
Пи-алкильное взаимодействие	Val156, Ile163	-
Контакт Ван-дер-Ваальса	Tyr370, Phe339, Ser242, Phe243, Phe332, Leu123	-

Таблица 7: Модели связывания, полученные при докинге триптамина с HTR2A (PDB ID: 6A93).Источник: инструмент молекулярной визуализации и двумерных диаграмм взаимодействий — двумерные диаграммы взаимодействий лиганд-остаток, представленные в разделе «Результаты» (Молекулярный докинг) рукописи. «-» означает, что значение расстояния для данного контакта не было указано отдельно.

(A) Валидация передокинга (положительные контроли)
PDB ID	Белок	Кокристаллический лиганд	Vina Score (kcal/mol)	RMSD (A)	Порог (A)	Результат
1DB1	VDR	VDX (аналог витамина D)	−13.0	0.87	2.0	ПРОЙДЕНО
3DCT	FXR	WAY-362450 (064)	−11.9	1.79	2.0	ПРОЙДЕНО
(B) Валидация кросс-докинга (отрицательные контроли)
Лиганд	Когнатная мишень (PDB)	Когнатный балл (kcal/mol)	Некогнатная мишень (PDB)	Некогнатный балл (kcal/mol)	Дельта (kcal/mol)	Селективность
Литохолевая кислота	VDR (1DB1)	−10.0	CASP3 (2DKO)	−8.3	1.7	Подтверждена
Триптамин	HTR2A (6A93)	−7.1	VDR (1DB1)	−6.4	0.7	Умеренная (в пределах погрешности Vina +/−0.5–1.0)

Таблица 8: Результаты валидации протокола докинга: значения RMSD при редокинге (положительные контроли) и баллы кросс-докинга (отрицательные контроли). Источник: Docking_Validation/Results/Docking_Validation_Results.xlsx и Docking_Validation/Logs/*.log (молекулярный докинг метаболитов-лигандов к белкам-мишеням, exhaustiveness = 8, seed = 42, бокс 25 Å × 25 Å × 25 Å). RMSD вычислено путем сопоставления тяжелых атомов по их именам (без суперпозиции).

Комплекс	RMSD (нм), среднее ± SD (диапазон)	Rg (нм), среднее ± SD (диапазон)	SASA (нм²), среднее ± SD (диапазон)	H-связи, среднее ± СД (диапазон)	RMSF (нм), среднее (макс.)
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 (макс. 0,600 в остатке 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 (макс. 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 (макс. 0,298)
Urolithin A-CASP3/2DKO	0.521 + / – 0.058 (0.244–0.755)	1.892 + / − 0.024 (1.839–1.984)	134.9 + / − 3.0 (126.4–146.4)	0.6 + / − 0.7 (0–3)	1,172 (макс. 2,532 в остатке 175)
Триптамин-HTR2A/6A93 (мембранная фракция)	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 (макс. 0,319)

Таблица 9: Сводные данные поведения в ходе молекулярно-динамического моделирования продолжительностью 200 ns для пяти приоритетных белково-лигандных комплексов, включая встроенную в мембрану систему триптамин-HTR2A.Источник: файлы утилит анализа траекторий молекулярной динамики (.xvg) — gmx rms, gmx gyrate, gmx sasa, gmx hbond, gmx rmsf — рассчитанные за последние 150 ns (50–200 ns) каждого 200 ns продуктивного запуска согласно шагу 8.8 протокола. RMSD/Rg рассчитаны с подгонкой по остову; радиус зонда SASA 0.14 nm; граничный критерий донора-акцептора H-связи 0.35 nm / 30 °. Для LCA-3DCT и UDCA-3DCT используется одна траектория остова белка (RMSD, RMSF) с отдельными значениями Rg/SASA/H-связей для каждого лиганда.

Триптамин-HTR2A/6A93 (мембрана) — количественное порезидуильное разложение
Остаток	Общий вклад $\Delta\Delta\text{G}$ (ккал/моль), среднее $\pm$ SD	Направление
Asp155 (D3.32)	−89.94 + / − 6.81	Стабилизирующий (доминантный)
Триптамин (лиганд)	−13.01 + / − 6.22	Стабилизация
Tyr171	13.62 + / − 4.54	Дестабилизация
Val167	23.32 + / − 3.96	Дестабилизирующий
Val156	20.03 + / − 3.81	Дестабилизация
Thr160	4.86 + / − 3.64	Дестабилизация
Ser159	24.16 + / − 3.48	Дестабилизирующий
Ser86	24.48 + / − 3.65	Дестабилизация
Phe87	35.18 + / − 4.04	Дестабилизация
Phe133	32.80 + / −3.70	Дестабилизирующий
Phe140	30.63 + / − 3.84	Дестабилизация
Phe141	35.25 + / − 3.55	Дестабилизация
Ile163	27.64 + / − 3.71	Дестабилизация
Trp137	53.77 + / − 4.32	Дестабилизирующий (наиболее неблагоприятный)
Остальные четыре комплекса — остатки, определенные при порезидуальном разложении (качественном)
Комплекс	Осадок	Направление
LCA-VDR/1DB1	Лиганд (LCA)	Благоприятный
LCA-VDR/1DB1	Gln317	Благоприятный
LCA-VDR/1DB1	Trp286	Неблагоприятный
LCA-NR1H4/FXR/3DCT	Arg331	Благоприятный (доминантный)
UDCA-NR1H4/FXR/3DCT	Glu326	Смешанная/распределенная сеть
UDCA-NR1H4/FXR/3DCT	Asp394	Смешанная/распределенная сеть
UDCA-NR1H4/FXR/3DCT	Arg395	Смешанная/распределенная сеть
UDCA-NR1H4/FXR/3DCT	Arg441	Смешанная/распределенная сеть
UDCA-NR1H4/FXR/3DCT	Asp470	Смешанная/распределенная сеть
Уролитин А-CASP3/2DKO	Arg64	Сильно благоприятный (полярный/электростатический)
Urolithin A-CASP3/2DKO	Arg207	Сильно благоприятное (полярное/электростатическое)

Таблица 10: Поостаточный анализ разложения MM-PBSA. КРАТКАЯ АННОТАЦИЯ: стабилизирующие и дестабилизирующие остатки (абсолютный вклад ≥ 0,5 kcal mol⁻¹) для каждого из пяти приоритетных белок-лигандных комплексов, включая встроенную в мембрану систему триптамин-HTR2A.Источник: Membrane Simulation/03_MMPBSA/results/FINAL_DECOMP_MMPBSA.dat (инструмент расчета энергии связывания методом молекулярной механики/континуального растворителя — поостаточное разложение Generalized Born (GB), «Complex: Total Energy Decomposition»). Номера остатков переведены из внутренней нумерации системы, созданной в CHARMM-GUI (смещение +68), в исходную нумерацию PDB 6A93, используемую в остальной части данной рукописи.

Источник: Данные предыдущего моделирования МД Simulation data/1DB1,2KD0, LCA & UDCA_3DCT}/mmpbsa_*/Decomposition_NORMAL_GB_Complex_TDC*.svg и раздел результатов рукописи (свободная энергия связывания MM-PBSA и разложение по остаткам). Для этих четырех комплексов в директории проекта отсутствуют числовые выходные данные по остаткам в форматах .dat/.csv (доступны только отрисованные SVG-графики с текстом в виде векторных путей, который невозможно извлечь машинным способом); в тексте рукописи указаны только идентификаторы остатков и направление влияния (благоприятное/неблагоприятное). Точные значения вкладов в ккал/моль для этих четырех комплексов в исходном репозитории отсутствуют.

Обсуждение

Данное поисковое вычислительное исследование демонстрирует интегрированный воспроизводимый рабочий процесс для приоритизации генов хозяина, связанных с микробными метаболитами, и комплексов белок-лиганд, примененный в данном случае к общедоступному набору транскриптомных данных слизистой оболочки прямой кишки при СРК-З. С помощью этого рабочего процесса подмножество предсказанных генов, связанных с микробными метаболитами, совпало с последовательно сниженными генами в наборе данных, которые объединились в кластеры путей GPCR, серотонинергических, кальциевых сигнальных, нейроактивных лиганд-рецепторных и ядерных рецептор-ассоциированных путей — систем, которые все чаще вовлекаются в коммуникацию между микробиотой и хозяином^48,49. Эти результаты следует интерпретировать строго как генерацию гипотез: анализ не измеряет концентрации микробных метаболитов, содержание рецепторных белков, связывание лигандов, активацию рецепторов, нисходящий сигналинг, моторику, секрецию, болевые реакции или клинические исходы. Наиболее обоснованным выводом является то, что выявленные гены и пути являются кандидатами для экспериментальной валидации, а не подтвержденными механизмами заболевания.

Критическая значимость данного протокола по сравнению с предыдущими работами, в которых отдельные пары «метаболит-рецептор» изучались изолированно, заключается в интеграции прогнозирования мишеней, общедоступных транскриптомных данных, сетевого анализа, докинга с контролями валидации, молекулярной динамики и метода MM-PBSA в единый последовательный конвейер приоритезации. Каждый этап сужает и контекстуализирует набор кандидатов, полученный на предыдущем этапе, и именно эта последовательная фильтрация делает итоговый список кандидатов доступным для экспериментальной проверки. Выявленные здесь модуль GPCR, центрированный на GNAQ, и модуль, ассоциированный с серотониновым рецептором, биологически обоснованы, учитывая роль сигналинга Gq в активации фосфолипазы C, выработке инозитол-1,4,5-трифосфата, мобилизации кальция, секреции и энтероэндокринной функции, а также установленную роль сигналинга короткоцепочечных жирных кислот и производных триптофана в гомеостазе слизистой оболочки и серотонинергического сигналинга в моторике желудочно-кишечного тракта, секреции, висцеральной чувствительности и взаимодействии оси «кишечник-мозг»^10,18,50,51.

Ключевой методологической особенностью данного исследования является подход к работе с встроенным в мембрану рецептором HTR2A. Поскольку моделирование в растворе не позволяет воспроизвести липидное окружение, определяющее конформационное поведение рецептора, сопряженного с G-белком, комплекс триптамин-HTR2A моделировался в эксплицитном бислое POPC. В условиях данной мембраны рецептор оставался структурно стабильным на протяжении всей траектории в 200 ns, а солевой мостик между аммонийной группой триптамина и Asp155 (D3.32) сохранялся практически в течение всего моделирования. Тот факт, что три независимых доказательства — поза докинга, постоянное расстояние контакта на протяжении всей траектории и доминирующий вклад каждого остатка в расчетах MM-PBSA — указывают на одно и то же консервативное взаимодействие с D3.32, подтверждает внутреннюю согласованность предсказанного способа связывания триптамина, который повторяет каноническую геометрию связывания аминергических лигандов с серотониновыми рецепторами.

При воспроизведении данного рабочего процесса следует учитывать определенные методологические аспекты. Ошибки в канонической структуре или наличие соединений, помеченных как Pan-assay interference compounds (PAINS), распространяются на этапы предсказания мишеней и докинга, что требует точного подбора метаболитов и хемоинформатической курации. Чтобы минимизировать наборы мишеней, вызванные шумом, необходимо последовательно применять критерии достоверности (предсказание взаимодействий химическое вещество-белок: ≥0.700; программа молекулярного докинга: ≥0.70; построение сети белок-белковых взаимодействий и анализ обогащения путей: ≥0.700). Предсказанные мишени следует группировать по функциональным категориям, чтобы избежать ошибочной характеристики всех ассоциированных с метаболитом генов как рецепторов. Важными аспектами докинга являются точная предварительная обработка структур PDB, минимизация энергии лигандов и размещение сетки вокруг известных остатков связывания, а введенные здесь контроли в виде редокинга и кросс-докинга обеспечивают объективную оценку корректности методологии докинга. Границы воспроизводимости в молекулярной динамике определяются сочетанием параметризации силового поля, правильной сольватации или построения мембраны, поэтапной эквилибрации и адекватного отбора проб на этапе продуктивного моделирования.

Типичные варианты адаптации и шаги по устранению неполадок включают снижение строгости порогов, если целевое предсказание не дает результатов, проверку согласованности направлений на уровне зондов для генов с несколькими зондами, а также интерпретацию отдельных узлов при построении сетей белок-белковых взаимодействий и анализе обогащения путей как зависящих от порога, а не как биологически нерелевантных. Для мембранных рецепторов следует использовать эксплицитное моделирование липидного бислоя вместо моделирования в водной среде, как показано на примере метода HTR2A, описанного здесь. В случаях, когда необходима декомпозиция энергии по остаткам, вычисления должны выполняться с использованием движка, поддерживающего декомпозицию, а сообщаемая нумерация остатков должна быть приведена в соответствие с нативной нумерацией рецептора во избежание неоднозначности. Мы полагаем, что результаты обогащения путей лучше рассматривать как организационный контекст для списка кандидатов, а не как валидацию на уровне путей. С технической точки зрения, обогащение терминов GPCR, серотонинергических путей или сигналинга кальция будет происходить всякий раз, когда список генов содержит несколько генов серотониновых рецепторов, независимо от корегуляции на уровне белков. Значения RMSD, Rg и RMSF для встроенной в мембрану системы HTR2A следует интерпретировать с учетом липидного бислоя: снижение Rg на более поздних этапах траектории может отражать вызванную бислоем конформационную адаптацию трансмембранного пучка, а не глобальное развертывание, а устойчивые водородные связи между лигандом и белком следует интерпретировать совместно с общей стабильностью RMSD.

Ограничения данного исследования существенны и накладывают определенные рамки на интерпретацию результатов. Исследование было основано на одном относительно небольшом общедоступном наборе данных; поиск в основных публичных репозиториях транскриптомных данных (веб-инструмент для анализа дифференциальной экспрессии генов и ArrayExpress) на момент анализа не выявил независимого набора данных транскриптома слизистой оболочки прямой кишки при СРК-З сопоставимого дизайна и платформы, который мог бы послужить когортой для репликации. Отсутствие независимой транскриптомной репликации является серьезным ограничением, и ни одно утверждение в данной рукописи не должно интерпретироваться как внешняя валидация результатов, полученных на основе одного набора данных. Набор данных демонстрирует почти повсеместную дифференциальную экспрессию (примерно 94,5% генов значимы, при этом подавляющее большинство из них понижено), что делает стандартную статистику обогащения неинформативной и не позволяет сделать выводы о специфичности подавления целевых генов относительно геномного фона; поэтому перекрытие представлено как описательный направленный паттерн, а не как статистическое обогащение. Общая транскриптомная analysis слизистой оболочки не позволяет отличить фактическую регуляцию генов от изменений в клеточном составе. Экспрессия мРНК не определяет содержание белка или функциональный ответ. Базы данных для предсказания мишеней подвержены систематической ошибке аннотирования, а результаты докинга, MD и MM-PBSA зависят от выбора силового поля, параметризации лиганда, начального положения, времени симуляции и адекватности выборки. Точные идентификаторы версий/сборок для некоторых веб-серверов и компонентов пакетов, включая сервис параметризации лигандов, совместимый с CHARMM, CHARMM-GUI и среду статистических вычислений, не указаны.

версии сборок пакетов и подверсии инструментов расчета энергии связывания методом молекулярной механики/континуального растворителя не были полностью восстановлены из архивной записи проекта и должны быть указаны как доступные в отдельной Таблице материалов. Значения MM-PBSA являются относительными оценками, не содержат явного члена конфигурационной энтропии и не должны интерпретироваться как экспериментальные показатели сродства. В исследовании отсутствуют метаболомные данные, вследствие чего невозможно определить, изменяется ли доступность лигандов при IBS-C, а также являются ли наблюдаемые изменения экспрессии причинами, следствиями, компенсаторными ответами или несущественными корреляциями.

Будущие применения данного метода должны включать независимую транскриптомную репликацию, количественную полимеразную цепную реакцию (qPCR) и валидацию на уровне белков, определение локализации по типам клеток с помощью одноклеточной или пространственной транскриптомики, метаболомное профилирование соответствующих классов метаболитов, а также функциональные анализы лигандного ответа в колоноидах, полученных от пациентов, эксплантатах слизистой оболочки или сопоставимых моделях. Сравнение с когортами пациентов с СРК с преобладанием диареи, смешанным типом СРК, воспалительными заболеваниями кишечника и запорами без СРК^1,2 поможет установить специфичность заболевания. Что касается структурного компонента, воспроизведение траекторий молекулярной динамики (MD), проведение анализа чувствительности с альтернативными начальными позами, а также полное документирование депозиции топологии, траекторий и входных и выходных файлов MM-PBSA способствовали бы дальнейшему повышению воспроизводимости. Экспериментальные анализы лигандного ответа остаются необходимыми для определения функциональной значимости приоритетных комплексов; представленные результаты не являются основанием для клинических или терапевтических утверждений.

Материалы

Список материалов, использованных в этой статье

Имя	Компания	Каталожный номер	Комментарии
AutoDock Vina	Scripps Research / открытый исходный код	v1.2.7; https://vina.scripps.edu/	Молекулярный докинг метаболитов-лигандов с белками-мишенями.
CGenFF/ParamChem	SilcsBio / Мэрилендский университет	v4.6; https://cgenff.com/	Параметризация силового поля лиганда для молекулярной динамики.
силовое поле CHARMM36m	разработчики CHARMM / открытое исходное приложение	CHARMM36m; https://www.charmm.org/charmm/resources/charmm-force-fields/	Силовое поле белка, используемое для моделирования молекулярной динамики.
Конструктор мембран CHARMM-GUI	CHARMM-GUI / Университет Лихай	Веб-сервер; точная версия не определена; https://www.charmm-gui.org/?doc=input/membrane	Создание и настройка уравновешивания системы явной мембраны POPC.
Discovery Studio Visualizer	BIOVIA (Dassault Systèmes)èмесяцы)	2021; https://discover.3ds.com/discovery-studio-visualizer-download	Двумерный анализ взаимодействий лиганд-остаток.
GEO2R	Архив экспрессии генов NCBI (Gene Expression Omnibus)	Веб-инструмент; доступ в январе-мае 2026 г.; https://www.ncbi.nlm.nih.gov/geo/geo2r/	Анализ дифференциальной экспрессии GSE36701.
GeneCards	Институт науки Вейцмана	Веб-база данных; доступ в январе–мае 2026 г.; https://www.genecards.org/	Верификация символов и информации о генах при стандартизации мишеней.
gmx_MMPBSA	Открытый исходный код (Vald)éс-Tresanco и др.)	1.5.x; https://valdes-tresanco-ms.github.io/gmx_MMPBSA/	Оценка свободной энергии связывания методом MM-PBSA и поостаточное разложение.
GROMACS	Команда разработчиков GROMACS / открытый исходный код	2024.2; https://www.gromacs.org/	Программный комплекс для моделирования молекулярной динамики.
транскриптомный набор данных GSE36701	NCBI Gene Expression Omnibus	GSE36701; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE36701	Общедоступный набор данных экспрессии слизистой оболочки прямой кишки при СРК с запорами.
Open Babel	Открытый исходный код	3.2.0; https://openbabel.org/	Конвертация форматов химических файлов, генерация трехмерных структур лигандов и подготовка лигандов.
База данных OPM	Мичиганский университет	Веб-база данных; доступ в январе-мае 2026 г.; https://opm.phar.umich.edu/	Ориентация белков в мембранах, координаты, использованные для выравнивания HTR2A.
ParmEd	разработчики ParmEd / открытое программное обеспечение	4.x; https://parmed.github.io/ParmEd/html/index.html	Перераспределение массы водорода и обработка топологии молекулярного моделирования.
PyMOL	Шрвödinger / открытое программное обеспечение	2.x; https://www.pymol.org/	Трехмерная структурная визуализация и схемы взаимодействия рецептора с лигандом.
Банк данных белков RCSB (RCSB Protein Data Bank)	RCSB PDB	Веб-база данных; доступ в январе-мае 2026 г.; https://www.rcsb.org/	Источник экспериментальных структур белков и координат PDB.
STITCH	консорциум STITCH (EMBL)	v5.0; https://stitch.embl.de/	Предсказание мишеней химического белкового взаимодействия.
Пожалуйста, предоставьте текст (STRING), который необходимо перевести.	Консорциум STRING / ELIXIR	v12.0; https://version-12-0.string-db.org/	Построение сети белок-белковых взаимодействий и анализ обогащения путей.
SwissADME	SIB Швейцарский институт биоинформатики / Лозаннский университет	Веб-инструмент; доступ в январе-мае 2026 г.; https://www.swissadme.ch/	Хемоинформационные дескрипторы, фармакокинетические прогнозы и оценка PAINS.
SwissTargetPrediction	SIB Швейцарский институт биоинформатики / Лозаннский университет	Веб-инструмент; доступ в январе–мае 2026 г.; https://www.swisstargetprediction.ch/	Прогнозирование белков-мишеней человека на основе лигандов.
Сопоставление идентификаторов UniProt	Консорциум UniProt	Веб-сервис; доступ в январе–мае 2026 г.; https://www.uniprot.org/id-mapping	Сопоставление идентификаторов белков со стандартизированными символами генов, утвержденными HGNC.
NVIDIA RTX 3080	NVIDIA Corporation	RTX 3080; ≥8 ГБ видеопамяти (VRAM); версия CUDA/драйвера в рукописи не указана	Графический процессор с поддержкой CUDA, используемый для моделирования молекулярной динамики.
GPU с поддержкой CUDA	NVIDIA Corporation	версия CUDA toolkit не указана в рукописи; ≥8 ГБ видеопамяти (VRAM)	Графический процессор (GPU) с поддержкой CUDA с ≥8 ГБ видеопамяти (VRAM); также потребовалась рабочая станция ≥32 ГБ оперативной памяти и 6-ядерный процессор.
Ubuntu Linux	Canonical Ltd. / открытое программное обеспечение	22.04 LTS	64-разрядная операционная система Linux.
Python 3.9	Python Software Foundation	3.9	Среда программирования общего назначения, используемая для написания сценариев рабочих процессов и анализа данных.
Архив экспрессии генов (Gene Expression Omnibus, GEO)	NCBI / Национальная медицинская библиотека США	Общедоступный веб-репозиторий; версия программного обеспечения в рукописи не указана	Общедоступный репозиторий данных функциональной геномики.
AutoDockTools/MGLTools	Лаборатория молекулярной графики, Scripps Research	1.5.7	Инструментарий для подготовки молекулярных структур и входных данных для докинга.
Инструменты анализа GROMACS	команда разработчиков GROMACS / открытое программное обеспечение	2024.2	Инструменты анализа траекторий молекулярной динамики.
Шестиэтапный протокол CHARMM-GUI	CHARMM-GUI / Университет Лихай	Веб-протокол; точная версия выпуска не определена	Веб-ориентированный многоэтапный рабочий процесс подготовки и уравновешивания молекулярных систем.
cgenff_charmm2gmx_py3.py	Скрипт конвертации с открытым исходным кодом; источник в рукописи не указан	Версия в рукописи не указана	Скрипт для конвертации топологии силового поля.
Python	Python Software Foundation	3.9	Среда программирования общего назначения.
SciPy	Сообщество SciPy / открытое программное обеспечение	Версия в рукописи не указана	Библиотека для научных вычислений.
scipy.stats.fisher_exact	Сообщество SciPy / открытое программное обеспечение	Версия SciPy не указана в рукописи	Fisher’реализация точного теста.
R	R Foundation for Statistical Computing	4.3.x	Среда статистических вычислений.
Bioconductor	Проект Bioconductor / открытое программное обеспечение	3.18	Программный комплекс для биоинформатики.
limma	Проект Bioconductor / открытое программное обеспечение	Версия не указана в рукописи	Пакет для анализа дифференциальной экспрессии генов.
Беньямини–Процедура Хохберга	Статистический метод	Неприменимо (статистическая процедура)	Метод корректировки ложного уровня значимости.
NVIDIA RTX 3080	NVIDIA Corporation	RTX 3080; ≥8 ГБ видеопамяти (VRAM); версия CUDA/драйвера в рукописи не указана	Графический процессор с объемом видеопамяти не менее 8 ГБ.
GPU с поддержкой CUDA	NVIDIA Corporation	версия инструментария CUDA не указана в рукописи; ≥8 ГБ видеопамяти (VRAM)	Графический процессор с поддержкой параллельных вычислений общего назначения.
Ubuntu Linux 22.04 LTS	Canonical Ltd. / открытое программное обеспечение	22.04 LTS	64-битная операционная система Linux.
TIP3P	Разработчики силового поля CHARMM / открытое программное обеспечение	TIP3P; версия программного обеспечения не применима	Трехточечная модель явной воды.
MM/PBSA	разработчики gmx_MMPBSA / открытое программное обеспечение	gmx_MMPBSA 1.5.x	Молекулярная механика/Пуассон–Метод расчета энергии связывания на основе площади поверхности Больцмана.