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

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

Данный протокол объединяет прогнозирование мишеней микробных метаболитов, транскриптомику слизистой оболочки прямой кишки, анализ белок-белковых взаимодействий и обогащения путей, молекулярный докинг, моделирование молекулярной динамики и оценку свободной энергии связывания методом молекулярной механики/поверхности Пуассона — Больцмана (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₂ кратности изменения (fold change) более 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, что отражает почти повсеместное подавление транскрипции при сравнении слизистой оболочки прямой кишки при СРК с запорами. Все 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 прогнозируемых мишеней были последовательно и значимо подавлены в ткани слизистой оболочки прямой кишки при СРК с запорами, а не как доказательство статистического обогащения или независимой валидации на геномном фоне. Формальное тестирование на обогащение потребовало бы репликации в транскриптомных наборах данных с более селективными профилями дифференциальной экспрессии, в которых значимость достигают существенно меньше половины всех генов. Следует подчеркнуть, что единообразное снижение экспрессии всех 17 перекрывающихся генов является описательным наблюдением, а не отдельно валидированным статистическим результатом, поскольку дифференциально экспрессируемый фон данного набора данных сам по себе преимущественно подавлен; общая направленность снижения экспрессии перекрывающихся генов была ожидаемой и не подвергалась формальному тесту на направленность. Следовательно, эта единообразная направленность не должна интерпретироваться как независимое статистическое доказательство скоординированной, метаболит-специфичной регуляции.

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

Пропионат имел наибольшее количество перекрывающихся генов, включая CASR, FFAR2, GCG, GNAQ, GPHN, GPR68, MLN и TBXA2R, что указывает на возможную вовлеченность сигналинга, чувствительного к короткоцепочечным жирным кислотам и связанного с Gq. Бутират перекрывался с HDAC3, что согласуется с биологией гистондеацетилаз, связанной с бутиратом, хотя снижение экспрессии мРНК само по себе не подтверждает изменение чувствительности к бутирату. Перекрывания, связанные с желчными кислотами, включали ядерные рецепторы VDR и NR1H4, которые являются признанными эффекторами сигналинга желчных кислот в кишечнике^38,39. Триптамин перекрывался с HTR1B, HTR2A, HTR2B и HTR6, что указывает на серотонинергический сигналинг как на потенциальный модуль — систему с хорошо изученной ролью в моторике и секреции желудочно-кишечного тракта⁴⁰. Индол-3-лактат и индол-3-пропионат перекрывались с KYAT1, а уролитин А перекрывался с 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, чем комплексы уролитин А-CASP3 и триптамин-HTR2A. Наилучший показатель был получен для LCA-VDR (−10.0 kcal/mol), за которым следовали LCA-NR1H4/FXR (−9.9 kcal/mol) и UDCA-NR1H4/FXR (−9.4 kcal/mol). Уролитин А-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).

В комплексе urolithin 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 по нумерации Баллестероса-Вайнштейна), который фиксирует протонированный амин аминергических лигандов в серотониновых и родственных рецепторах^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: Вычислительный рабочий процесс приоритизации генов хозяина, связанных с метаболитами, при СРК-З. Схематическое представление восьмиэтапного рабочего процесса, объединяющего отбор метаболитов, предсказание мишеней, транскриптомный анализ дифференциальной экспрессии, анализ перекрытий, анализ обогащения сетей и путей, молекулярный докинг, симуляцию молекулярной динамики и расчет свободной энергии связывания методом MM-PBSA. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

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

Рисунок 3: Статистический анализ 17 предсказанных генов-мишеней метаболитов по сравнению с GSE36701. (A) log2 кратность изменения для каждого из 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) Трехмерное представление поверхности и структуры в виде лент, где литохолевая кислота показана в виде сфер. (B) Двумерная карта взаимодействий, демонстрирующая водородную связь с Ser278, а также окружающие гидрофобные контакты и контакты Ван-дер-Ваальса. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

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

Рисунок 7: Трехмерное и двумерное структурное представление урсодезоксихолевой кислоты в комплексе с NR1H4/FXR (PDB ID: 3DCT). (A) Трехмерное представление поверхности и упрощенная модель (cartoon). (B) Двумерная карта взаимодействий, демонстрирующая водородные связи с His447 и Gly322, π-анионное взаимодействие с 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: Анализ траектории молекулярной динамики комплекса urolithin 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. Рецептор представлен в виде ленточной модели, пересекающей бислой; липиды 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). Все перечисленные перекрывающиеся гены были подавлены. Таблица 1 представлена отдельно в виде электронной таблицы и содержит для каждой мишени: соответствующий метаболит(ы), функциональную категорию, источник прогнозирования мишени (прогноз взаимодействия химическое вещество-белок, программа молекулярного докинга или и то, и другое), комбинированный показатель взаимодействия прогноза взаимодействия химическое вещество-белок и вероятность программы молекулярного докинга, где применимо, уровень прогнозирования, log2 fold change и FDR с уровнем значимости экспрессии. Источник: значения экспрессии генов были получены из таблицы дифференциальной экспрессии GSE36701 с объединением по генам (зонд с наименьшим FDR на ген). Источники прогнозирования мишени и значения достоверности были скомпилированы на основе результатов прогноза взаимодействия химическое вещество-белок и программы молекулярного докинга с использованием пороговых значений комбинированного показателя взаимодействия прогноза взаимодействия химическое вещество-белок ≥ 0.700 и вероятности программы молекулярного докинга ≥ 0.70. Показатели прогноза взаимодействия химическое вещество-белок представляют собой комбинированные баллы по шкале 0–1; STP обозначает вероятность программы молекулярного докинга. Уровень 1 (Tier 1) = строгая поддержка прогноза взаимодействия химическое вещество-белок; Уровень 1+ (Tier 1+) = строгая поддержка прогноза взаимодействия химическое вещество-белок, подтвержденная программой молекулярного докинга.

Комплекс	Белок (PDB ID)	Лиганд	Балл Vina (ккал/моль)	Размер полости (A^3)	Центр сетки X,Y,Z (A)	Область поиска (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
UDCA-NR1H4/FXR	NR1H4/FXR (3DCT)	Урсодезоксихолевая кислота	−9.4	3395	137, 31, 78	25 x 25 x 25
Urolithin A-CASP3	CASP3 (2DKO)	Уролитин A	−7.1	233	37, 34, 32	25 x 25 x 25
Tryptamine-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).

Тип взаимодействия	Остаток(и)	Расстояние (A)	Примечания
Классическая водородная связь	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
Пи-анионное взаимодействие	Val325	4.96
Углерод-водородная связь	Trp469	4.51
Неблагоприятный донор-донорный контакт	Arg395	3.89
Неблагоприятный донор-донорный контакт	Gln396	3.40

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

Тип взаимодействия	Остаток(и)	Расстояние (A)	Примечания
Классическая водородная связь	Gln161	3.78
Классическая водородная связь	Gln161	4.19	вторичный контакт
Классическая водородная связь	Ser120	3.95
Обычная водородная связь	Arg207	3.05
Обычная водородная связь	Arg207	3.77	второй контакт
$\pi$-катионное взаимодействие	Arg207	-
Водородная связь с π-донором	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 (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 (nm), среднее + / – SD (диапазон)	Rg (nm), среднее + / – SD (диапазон)	SASA (nm^2), среднее + / – SD (диапазон)	H-связи, среднее + / – SD (диапазон)	RMSF (nm), среднее (макс.)
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)
Tryptamine-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 (мембрана) — количественный анализ вклада отдельных остатков
Остаток	Суммарный вклад ddG (kcal/mol), среднее ± 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	Смешанная/распределенная сеть
Уролитин A-CASP3/2DKO	Arg64	Сильно благоприятное (полярное/электростатическое)
Уролитин A-CASP3/2DKO	Arg207	Сильно благоприятное (полярное/электростатическое)

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

Источник: данные предыдущего моделирования молекулярной динамики Previous_MD 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	Молекулярная механика/Пуассон–Метод расчета энергии связывания на основе площади поверхности Больцмана.