In silico-приоритизация потенциальных белковых взаимодействий для глутатион-чувствительных пролекарств, блокирующих апуриниковые/апиримидиновые сайты, при немелкоклеточном раке легкого

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

В данном исследовании были интегрированы методы сетевой фармакологии, молекулярного докинга, молекулярной динамики и MM-PBSA для определения приоритетных потенциальных белковых взаимодействий пролекарств, улавливающих апуриновые/апиримидиновые сайты и реагирующих на глутатион, при НМРЛ. Соединение 5 было выбрано в качестве приоритетного для последующего изучения взаимодействий с MMP9 и SRC; тем не менее, они остаются расчетными кандидатами. Полученные результаты носят гипотетический характер и требуют биохимической и клеточной валидации.

Аннотация

Немелкоклеточный рак легкого (НМРЛ) остается одной из основных причин смертности от онкологических заболеваний, а эффективность лечения часто ограничена развитием приобретенной резистентности и системной токсичностью. В предыдущих экспериментальных исследованиях пролекарства с ловушками для апурининовых/апиримидиновых сайтов, реагирующие на глутатион, продемонстрировали селективную противоопухолевую активность, однако вопрос о том, проявляют ли их высвобождаемые метаболиты значимое взаимодействие с белками, имеющими отношение к раку, остается открытым. В данной работе представлен интегрированный in silico был использован рабочий процесс, объединяющий сетевую фармакологию, молекулярный докинг, молекулярную динамику (МД) и анализ площади поверхности Пуассона — Больцмана методом молекулярной механики (MM-PBSA), для приоритизации проверяемых гипотез белкового взаимодействия для двух глютатион-чувствительных пролекарств (соединение 1 и соединение 2), их аминоокси-содержащих продуктов (соединение 4 и соединение 5), а также соответствующей пары незахватывающих контрольных соединений (соединение 3 и соединение 6). Были выявлены двадцать один пересекающихся мишеней, связанных с соединениями и заболеванием; на основании топологии белок-белковых взаимодействий приоритет был отдан серин-треониновой киназе 1 AKT (AKT1), рецептору эпидермального фактора роста (EGFR), фактору некроза опухоли (TNF), матриксной металлопротеиназе 9 (MMP9) и нерецепторной тирозинкиназе протоонкогена SRC (SRC). Соединение 5 показало наиболее благоприятный единичный балл AutoDock Vina с MMP9 (-8,418 ккал·моль⁻¹) и продемонстрировала сравнительно стабильные позы, полученные в результате докинга, на протяжении траекторий молекулярной динамики (МД) для MMP9 и SRC. Однако наиболее высокорейтинговая поза для MMP9 не показала прямой координации каталитического Zn²⁺ взаимодействие или прямое связывание с His401, Glu402, His405 или His411 не наблюдается, поэтому данному соединению не может быть приписан канонический способ ингибирующего связывания с MMP9. Анализы методом молекулярной динамики (MD) и MM-PBSA характеризуют только поведение траекторий и относительный энергетический рейтинг этих комплексов; они не демонстрируют внутриклеточного связывания с мишенью, ингибирования фермента или регуляции сигнального пути. Таким образом, установленная активность по захвату абазических сайтов и новые предсказанные белковые взаимодействия рассматриваются как отдельные, потенциально параллельные гипотезы, а не как доказанная механистическая цепочка. В целом, полученные результаты позволяют выделить приоритетные пары «соединение-мишень» для дальнейшего тестирования, но не устанавливают многомишеньный механизм борьбы с НМРЛ.

Введение

Немелкоклеточный рак легкого (НМРЛ) является наиболее распространенным гистологическим подтипом рака легкого и остается основной причиной смертности от онкологических заболеваний во всем мире¹. Несмотря на значительный прогресс в области таргетной терапии и прецизионной онкологии, долгосрочная эффективность лечения все еще часто снижается из-за развития приобретенной резистентности, ограниченной длительности ответа и связанной с лечением токсичности². Препараты, нацеленные на рецептор эпидермального фактора роста (EGFR), улучшили исходы у пациентов с определенным молекулярным профилем, однако в ходе лечения почти неизбежно развивается резистентность, что создает острую необходимость в терапевтических стратегиях, действующих через альтернативные или дополняющие механизмы^3,4. Традиционная химиотерапия на основе препаратов платины остается важным компонентом лечения, но ее клиническая польза ограничена кумулятивной токсичностью и развитием резистентности при длительном применении^5,6. В совокупности эти ограничения подчеркивают необходимость поиска противоопухолевых агентов, которые отличались бы по механизму действия и селективно активировались в опухолевом микроокружении.

Среди эндогенных повреждений ДНК абазические или апиримидиновые сайты встречаются очень часто и имеют большое биологическое значение: ежедневно в каждой клетке образуются тысячи таких повреждений⁷. При отсутствии эффективного восстановления эти повреждения могут превратиться в разрывы цепей, что способствует геномной нестабильности и гибели клеток^8,9. Апуриновая/апиримидиновая эндонуклеаза 1 является центральным ферментом пути эксцизионной репарации оснований, поскольку она расщепляет абазические сайты и обеспечивает последующую репарацию^10,11. Эта зависимость от репарации сделала повреждения, связанные с абазическими сайтами, привлекательной мишенью для разработки противоопухолевых препаратов¹². Исходя из этого обоснования, ранее были разработаны реагирующие на глутатион пролекарства для захвата абазических сайтов, чтобы использовать повышенную концентрацию глутатиона в клетках опухоли. Продукты их расщепления под действием глутатиона содержат аминоокси-функциональную группу, способную захватывать альдегидные абазические сайты; предыдущие экспериментальные работы продемонстрировали селективную цитотоксичность, остановку клеточного цикла и апоптоз в клетках H1299¹³. Эти данные подтверждают эффективность компонента захвата повреждений ДНК в дизайне соединения. Однако они не показывают, что MMP9, SRC, EGFR, AKT1 или TNF регулируются ниже по каскаду после захвата абазических сайтов. Следовательно, любое взаимодействие с белком-мишенью, выявленное в данном вычислительном анализе, должно рассматриваться как отдельная гипотеза, если только оба процесса не будут продемонстрированы в одной и той же биологической системе.

Настоящее исследование было спланировано с учетом этого различия. Первый уровень биологического обоснования основывается на ранее установленных механизмах глутатион-зависимого высвобождения и химического захвата в апиримидиновых/апуриновых сайтах¹⁴. Второй уровень, рассматриваемый здесь, представляет собой исследовательский вопрос: совместимы ли исходные пролекарства или продукты их высвобождения с выбранными рако-ассоциированными белками с точки зрения компьютерного моделирования. Для приоритизации белков-кандидатов использовался метод сетевой фармакологии, после чего были проведены докинг, молекулярная динамика (MD) и анализ MM-PBSA для изучения выбранных белок-лигандных комплексов^15,16. Данные расчеты не ставили целью доказать, что предсказанные белки опосредуют известный фенотип повреждения ДНК, или установить причинно-следственную связь между захватом в апиримидиновых/апуриновых сайтах и онкогенной сигнализацией. Вместо этого данный рабочий процесс был использован для создания ранжированного набора экспериментально проверяемых гипотез, которые в дальнейшем могут быть оценены с помощью анализов прямого связывания, ферментативной активности, сигнальных путей, повреждения ДНК и анализа фенотипа.

Протокол

Одобрение институционального наблюдательного совета не требовалось, поскольку данное исследование было полностью вычислительным и не включало участия людей, позвоночных животных, биообразцов, полученных от пациентов, или идентифицируемых персональных данных. Таким образом, информированное согласие не применялось. Все базы данных, программные пакеты, силовые поля и вычислительные ресурсы, использованные в данном протоколе, перечислены в Таблице материалов.

Изучаемые соединения и аналитический рабочий процесс

В качестве исследуемых молекул были использованы ранее описанные соединения, способные захватывать апуринированные сайты в присутствии глутатиона¹⁷. Соединение 1 и соединение 2 были выбраны в качестве исходных пролекарств, поскольку расщепление, инициируемое глутатионом, приводит к образованию продуктов, содержащих аминоокси-группу, — соединения 4 и соединения 5 соответственно. Соединения 4 и 5 были выбраны для структурного анализа, так как они представляют собой высвобождаемые формы, сохраняющие аминоокси-группу, реакционноспособную по отношению к апуринированным сайтам. Соединение 3 было включено в качестве соответствующего структурного контроля, чувствительного к глутатиону; его продукт расщепления, соединение 6, не обладает аминоокси-функциональностью, необходимой для ковалентного захвата альдегидов апуринированных сайтов. Соответственно, соединения 1–3 были включены в обратный поиск мишеней для сравнения исходных каркасов, соединения 4 и 5 оценивались в отношении приоритетных белков, а соединение 6 использовалось в качестве лиганда отрицательного контроля при сравнении методом молекулярной динамики (MD) с SRC. Траектория SRC–соединение 6 была включена для обеспечения сопоставимого структурного сравнения с системой соединения 5, содержащего аминоокси-группу, а не в качестве доказательства того, что сама аминоокси-группа определяет связывание с SRC. Химические структуры и взаимосвязи активации соединений 1–6 показаны на Рисунке 1. Такой дизайн намеренно отделяет установленную химию захвата повреждений ДНК от настоящего анализа по поиску гипотез о возможных белковых взаимодействиях.

Предсказание мишеней, связанных с соединением

Двумерные структуры соединения 1 (Compound 1), соединения 2 (Compound 2) и соединения 3 (Compound 3) были сохранены в формате MDL MOL и преобразованы в канонические строки упрощенного представления системы ввода молекул (SMILES) с использованием Open Babel версии 3.1.1 с форматом вывода канонических SMILES¹⁸. Каждая экспортированная строка была повторно импортирована, и регенерированная структура подверглась визуальной перекрестной проверке с соответствующей двумерной структурой перед отправкой в SwissTargetPrediction с ограничением по виду Homo sapiens¹⁹. Были сохранены предсказанные мишени с ненулевыми значениями вероятности. Списки мишеней, полученные для трех соединений, были объединены, дубликаты удалены, а оставшиеся мишени приведены к официальным символам генов человека перед дальнейшим анализом. Стандартизированные пары «соединение-мишень» были импортированы в программное обеспечение для визуализации и анализа сетей в виде сетевой таблицы, где соединения и предсказанные мишени представлены в виде узлов, а связи между соединением и мишенью — в виде ребер, для визуализации предсказанных взаимосвязей мишеней²⁰.

Поиск мишеней, связанных с заболеванием, и идентификация общих мишеней

Мишени, связанные с заболеванием, были извлечены из базы данных GeneCards с использованием поискового запроса «lung cancer H1299»²¹. Для отбора генов с более высокой степенью соответствия данному запросу был применен порог показателя релевантности более 0,27. Символы генов были стандартизированы, а дубликаты удалены вручную. Пересечение между предсказанными мишенями соединения и мишенями, связанными с заболеванием, было определено с помощью онлайн-инструмента для построения диаграмм Венна и анализа пересечений²². Для последующего анализа белок-белковых взаимодействий (PPI), обогащения и приоритизации мишеней были отобраны только пересекающиеся мишени.

Анализ белок-белковых взаимодействий и скрининг основных мишеней

Пересекающиеся мишени были введены в инструмент поиска взаимодействующих генов/белков (STRING) версии 11.5 с ограничением по виду до Homo sapiens и минимальный требуемый балл взаимодействия установлен на уровне 0,400²³Полученные данные о белок-белковых взаимодействиях (PPI) были импортированы в Cytoscape версии 3.10.0 для визуализации и топологического анализа. В рамках анализа топологии сети для каждого узла в исходной сети белок-белковых взаимодействий (PPI) из 21 узла были рассчитаны такие показатели, как plugin, степень, центральность по посредничеству и центральность по близости; полученные значения затем использовались для последовательной фильтрации на основе медианы, описанной ниже²⁴Показатели центральности были рассчитаны для исходной сети белок-белковых взаимодействий (PPI) из 21 узла, а затем использованы для последовательной фильтрации. Медианная степень исходной сети составила 12; были отобраны узлы со степенью ≥ 12, что позволило выявить 13 кандидатов. Среди этих 13 кандидатов медианная центральность по посредничеству составила 0,031293, а медианная центральность по близости — 0,769231. При втором этапе фильтрации были отобраны узлы с центральностью по посредничеству ≥ 0,031293 и центральностью по близости > 0,769231, что позволило выделить пять потенциальных хабов: AKT1, EGFR, TNF, MMP9 и SRC. Сообщенные значения центральности были перенесены из исходной сети из 21 узла, а не пересчитаны после выделения подмножества. Эти пять итоговых потенциальных хабов были использованы для последующего структурного анализа. Исходные показатели топологии соответствуют ненаправленной сети белок-белковых взаимодействий (PPI), содержащей 21 узел и 116 ребер.

Анализ обогащения с использованием Gene Ontology и Киотской энциклопедии генов и геномов

Пересекающиеся мишени были подвергнуты анализу обогащения по базе данных Gene Ontology (GO) и Киотской энциклопедии генов и геномов (KEGG) с использованием Database for Annotation, Visualization and Integrated Discovery (DAVID), при этом поиск был ограничен видом Homo sapiens²⁵. Обогащение GO оценивали по биологическим процессам (BP), клеточным компонентам (CC) и молекулярным функциям (MF), а также по сигнальным путям KEGG. Для данного разведочного анализа пороговым значением для включения записей об обогащении было номинальное p < 0.10; скорректированные по Бенджамини p-значения, значения Бонферрони, показатели частоты ложноположительных результатов и значения точного теста Фишера были приведены в дополнительных таблицах, но не использовались для определения итогового набора. Отобранные записи были ранжированы по номинальному p-значению. Для последующей визуализации и интерпретации были отобраны 20 основных путей KEGG и 10 основных терминов из каждой категории GO. Столбчатые диаграммы и пузырьковые диаграммы были построены с помощью онлайн-инструмента для биоинформатической визуализации.

Подготовка рецепторов и лигандов для молекулярного докинга

Кристаллические структуры приоритетных белков были получены из Protein Data Bank (PDB) RCSB: AKT1, PDB ID 3O96; EGFR, PDB ID 5UWD; TNF, PDB ID 2AZ5; MMP9, PDB ID 1GKC; и SRC, PDB ID 2H8H. Структуры белков, использованные для докинга, были подготовлены в программе для молекулярной визуализации путем удаления кокристаллизованных лигандов и молекул воды с последующей обработкой в программном обеспечении для подготовки к докингу. Для MMP9 была сохранена отдельная немодифицированная копия PDB 1GKC в качестве кристаллографического эталона для окружения каталитического Zn²⁺ и способа связывания N²-[(2R)-2-{[formyl(hydroxy)amino]methyl}-4-methylpentanoyl]-N,3-dimethyl-L-valinamide (NFH). В цепи A структуры 1GKC каталитический Zn²⁺ координируется His401, His405 и His411 на расстояниях 2.21, 2.23 и 2.22 Å соответственно, в то время как два атома кислорода NFH координируют Zn²⁺ на расстояниях 2.07 и 2.38 Å; Glu402 является каталитическим кислотно-основным остатком. Эти кристаллографические контакты использовались в качестве положительного структурного эталона для оценки карты взаимодействий при докинге Соединения 5. Для расчетов молекулярной динамики (MD) SRC использовалась структура PDB 2H8H. Эти структурные сравнения использовались только для интерпретации и не рассматривались как доказательство ингибирования фермента или взаимодействия с внутриклеточной мишенью.

Трехмерные структуры соединения 4 и соединения 5 были получены из базы данных PubChem и оптимизированы путем минимизации энергии в программе для молекулярного моделирования²⁶. Лиганды были протонированы в физиологических условиях и минимизированы с использованием молекулярного силового поля Merck 94 (MMFF94) до тех пор, пока градиент энергии не стал ниже 0,01 kcal·mol⁻¹·Å⁻¹. В программе для подготовки к докингу были добавлены атомы водорода, назначены заряды Гастейгера и определены вращающиеся связи. Соединение 6, продукт расщепления глутатионом контрольного несекуляризующего соединения 3, было подготовлено по той же схеме и подвергнуто докингу к SRC только для создания исходной позы для траектории МД соответствующего контроля. Таким образом, соединение 6 не вводилось в качестве дополнительного прогнозируемого терапевтического лиганда и не использовалось для подтверждения мультитаргетного механизма.

Процедура молекулярного докинга

Молекулярный докинг проводили с помощью программного обеспечения для молекулярного докинга по протоколу полугибкого докинга, при котором рецепторы оставались жесткими, а лигандам позволялось сохранять гибкость^27,28. Для каждого рецептора центр окна докинга совмещали с положением кокристаллизованного лиганда, чтобы область поиска соответствовала экспериментально определенному сайту связывания. Размеры окна докинга были установлены на уровне 24 Å × 24 Å × 24 Å для AKT1, EGFR и SRC, 26 Å × 26 Å × 26 Å для MMP9 и 28 Å × 28 Å × 28 Å для TNF. Параметр exhaustiveness (полнота поиска) был установлен на 32, количество выходных поз — на 20, а энергетический диапазон — на 4 kcal·mol⁻¹. Конформация с наивысшим рейтингом из каждого цикла докинга была отобрана для анализа взаимодействий.

Для оценки внутренней надежности настроек докинга каждый кокристаллизованный лиганд был заново задокирован в соответствующий карман рецептора с использованием тех же параметров; критерием приемлемости служило среднеквадратичное отклонение (RMSD) тяжелых атомов ниже 2.0 Å. Финальные позы анализировали в программном обеспечении для молекулярной визуализации. Для комплекса MMP9–Compound 5 сохраненная карта взаимодействий была специально проверена на наличие аннотированной прямой координации лиганда с каталитическим ионом Zn²⁺ и контактов с His401, Glu402, His405 и His411, а положение лиганда сравнивали с кристаллографической позой NFH в PDB 1GKC. Количественное расстояние металлической координации указывалось только в тех случаях, когда прямая координация лиганда-Zn²⁺ была очевидна в сохраненной позе; в противном случае она обозначалась как «неприменимо» вместо того, чтобы быть выведенной косвенно. Таким образом, поза докинга MMP9, в которой отсутствовали эти канонические особенности каталитического центра, классифицировалась как неканоническая и не интерпретировалась как доказательство ферментативного ингибирования MMP9. В более общем смысле показатели AutoDock Vina и позы докинга использовались для относительной приоритизации и генерации гипотез, а не в качестве доказательства аффинности связывания или взаимодействия с мишенью внутри клетки.

Протокол молекулярной динамики

Комплексы белок-лиганд, выбранные для МД-анализа, были сконструированы на основе поз связывания, полученных в результате докинга. МД-расчеты проводились с помощью GROMACS²⁹. Для параметризации белков использовалось силовое поле CHARMM36, в то время как типы атомов и параметры лигандов были назначены с помощью общего янтарного силового поля второго поколения (GAFF2)^30,31. Частичные заряды по модели Остина 1-bond charge correction (AM1-BCC) были созданы с помощью инструментов параметризации малых молекул и генерации топологии лигандов; затем программное обеспечение для генерации топологии лигандов создало файлы топологии лигандов, совместимые с GROMACS^32,33,34. Каждый комплекс был помещен в бокс с водой TIP3P с периодическими граничными условиями при минимальном расстоянии от растворенного вещества до стенки бокса 1,0 nm. Ионы натрия и хлорида добавлялись для нейтрализации суммарного заряда в каждой системе, а дополнительно добавляли NaCl для достижения конечной ионной силы 0,15 M.

Минимизация энергии проводилась с использованием алгоритма наискорейшего спуска до тех пор, пока максимальная сила не упала ниже 1000 kJ·mol⁻¹·nm⁻¹. Затем минимизированные системы уравновешивали в условиях постоянного числа частиц, давления и температуры (NPT) при 310 K и 1 bar с применением позиционных ограничений к остову белка. Контроль температуры осуществлялся с помощью термостата V-rescale, а давления — с помощью баростата Parrinello-Rahman. Дальнодействующие электростатические взаимодействия рассчитывали методом Particle Mesh Ewald. Радиус отсечки для короткодействующих электростатических взаимодействий и взаимодействий Ван-дер-Ваальса был установлен на уровне 1.0 nm; все связи с участием атомов водорода были ограничены с помощью алгоритма линейного решения ограничений (LINCS). Продукционные траектории генерировали в течение 150 ns с шагом интегрирования 2 fs, при этом координаты сохраняли каждые 10 ps для последующего анализа.

Анализ траекторий

Анализ траекторий проводился для уравновешенных участков продукционных траекторий. RMSD остова белка и RMSD лиганда рассчитывали после применения метода наименьших квадратов для совмещения с начальной референсной конформацией. Значения среднеквадратичных флуктуаций (RMSF) вычисляли для каждого остатка по атомам Cα. Анализ водородных связей между каждым лигандом и его рецептором проводили с использованием порогового значения расстояния между донором и акцептором 3.5 Å и порогового значения угла донор-водород-акцептор 30°. Степень заполнения водородной связи определяли как доля проанализированных кадров, в которых присутствовала данная водородная связь. Эти показатели использовали для характеристики структурной стабильности, гибкости на уровне остатков и устойчивости межмолекулярных контактов в трех траекториях MD.

Расчет свободной энергии связывания методом молекулярной механики и поверхности Пуассона-Больцмана (MM-PBSA)

Свободная энергия связывания была рассчитана с помощью метода MM-PBSA, реализованного в программном обеспечении для MM-PBSA анализа, на основе уравновешенных сегментов траекторий МД³⁵. Для каждого комплекса белок-лиганд для анализа свободной энергии использовались последние 50 ns из 150 ns продуктивной траектории. Всего было отобрано 500 равномерно распределенных кадров с интервалом 100 ps в диапазоне от 100.0 ns до 149.9 ns; конечная точка 150.0 ns была исключена из выборки. Общая свободная энергия связывания рассчитывалась как сумма энергии ван-дер-Ваальса, электростатической энергии, энергии полярной сольватации и энергии неполярной сольватации:

ΔG_bind = ΔE_vdW + ΔE_ele + ΔG_polar + ΔG_nonpolar (1)

Поправка на энтропию не применялась. Средние значения свободной энергии связывания и стандартные отклонения были рассчитаны по всем выбранным кадрам.

Контроль воспроизводимости и компьютерная валидация

Все названия генов, использованные на этапах прогнозирования мишеней, поиска связей «болезнь-мишень» и анализа пересечений, были приведены к официальным символам генов человека перед проведением последующего анализа. Для всех пар рецептор-лиганд поддерживались идентичные параметры докинга. Техническая проверка настроек докинга проводилась путем повторного докинга соответствующего кокристаллизованного лиганда в карман каждого рецептора с использованием тех же параметров поиска и заранее определенного порога допустимого RMSD по тяжелым атомам <2.0 Å. Кристаллографический комплекс MMP9–NFH был сохранен в качестве положительного структурного эталона для каталитического окружения Zn²⁺, в то время как SRC–Compound 6 служил соответствующим неловушечным компаратором для анализа траекторий SRC. Перед включением в последующие анализы траектории MD проверялись на стабильность температуры и давления, а также на отсутствие аномального дрейфа объема ячейки. Все входные файлы, структуры рецепторов, файлы топологии лигандов, файлы конфигурации докинга, файлы параметров MD и записи выбора кадров MM-PBSA были заархивированы для обеспечения воспроизводимости вычислений.

Результаты

Предсказание мишеней, связанных с соединением, и идентификация пересекающихся мишеней

Обратное прогнозирование мишеней для соединений 1–3 позволило выявить 212 потенциальных мишеней человека. Параллельно с этим поиск мишеней заболеваний в базе данных GeneCards по ключевому слову «lung cancer H1299» с пороговым значением показателя релевантности выше 0,27 позволил определить 188 мишеней, связанных с НМРЛ. Анализ пересечения между прогнозируемыми мишенями соединений и пулом мишеней, ассоциированных с заболеванием, выявил 21 общую мишень, которые были отобраны для всех последующих анализов. Пересечение двух наборов мишеней показано на рисунке 2. Эти результаты указывают на то, что исследуемые соединения воздействуют на ограниченный поднабор значимых для заболевания мишеней, а не на диффузное неспецифическое пространство мишеней.

Анализ белок-белковых взаимодействий и скрининг основных мишеней

21 пересекающихся мишени были импортированы в STRING для построения сети белок-белковых взаимодействий (PPI). Полученная сеть содержала 21 узел и 116 ребер. Сумма значений степени узлов в Дополнительной таблице 1 составляет 232, что соответствует 116 ненаправленным ребрам. Для этой исходной сети были рассчитаны показатели центральности. Медианная степень составила 12; применение критерия степени ≥ 12 позволило оставить 13 кандидатов. Среди этих 13 кандидатов медианная посредническая центральность (betweenness centrality) составила 0,031293, а медианная центральность по близости (closeness centrality) — 0,769231. Применение критериев посреднической центральности ≥ 0,031293 и центральности по близости > 0,769231 позволило выделить пять кандидатов в хабы: AKT1, EGFR, TNF, MMP9 и SRC. Значения, указанные для отобранных подмножеств, представляют собой метрики центральности исходной сети из 21 узла, перенесенные в процессе фильтрации; они не пересчитывались для подсетей из 13 или 5 узлов. Их сетевая центральность использовалась только для ранжирования кандидатов для последующей оценки на основе структуры и не должна интерпретироваться как доказательство того, что они являются биологическими мишенями соединений. В Дополнительной таблице 1, Дополнительной таблице 2 и Дополнительной таблице 3 приведены метрики для исходной сети из 21 узла, 13 кандидатов, оставшихся после первого этапа скрининга, и пяти финальных кандидатов в хабы соответственно. Последовательный скрининг сети PPI и пять финальных кандидатов в хабы показаны на Рисунке 3A, Рисунке 3B, Рисунке 3C, и Рисунке 3D.

Анализ обогащения по Gene Ontology и Kyoto Encyclopedia of Genes and Genomes

Функциональный анализ обогащения 21 общего мишени выявил 121 путь KEGG, соответствующих номинальному уровню значимости p < 0.10 критерии включения. 20 наиболее значимых путей представлены на Рисунок 4AСреди них особенно выделялись эндокринная резистентность, раковые пути, протеогликаны при раке, резистентность к ингибиторам тирозинкиназы EGFR и сигнальный путь ErbB. Эти пути тесно связаны с пролиферацией опухоли, выживаемостью, инвазией и резистентностью к терапии при НМРЛ. Полная статистика обогащения для всех 121 сохраненного пути KEGG, соответствующих номинальному p < 0,10 представлены в Дополнительная таблица 4.

Анализ обогащения GO дополнительно выявил 177 терминов BP, 29 терминов CC и 61 термин MF, соответствующих тому же номинальному критерию включения p < 0,10. 10 наиболее значимых терминов из каждой категории показаны на Рисунке 4B, Рисунке 4C, и Рисунке 4D. Доминирующие биологические процессы включали положительную регуляцию пролиферации гладкомышечных клеток сосудов, переход G2/M митотического клеточного цикла, сигнализацию рецептора инсулиноподобного фактора роста, фосфорилирование белков, отрицательную регуляцию апоптоза и передачу сигнала. Основными терминами клеточного компонента были ядро, мембранный рафт, фокальная адгезия, плазматическая мембрана и теломерная область хромосомы. Напротив, основные термины молекулярной функции, представленные на Рисунке 4D, включали протеинкиназную активность, серинпротеинкиназную активность, связывание ATP, серин/треонинопротеинкиназную активность, тирозинпротеинкиназную активность, киназную активность гептапептидных повторов CTD РНК-полимеразы II, киназную активность, киназную активность по H2AXY142 гистона, киназную активность по H3Y41 гистона и связывание идентичных белков. Полные данные по обогащению для BP, CC и MF приведены в Дополнительной таблице 5, Дополнительной таблице 6 и Дополнительной таблице 7. В совокупности эти результаты обогащения указывают на то, что пересекающийся набор мишеней сосредоточен в процессах сигналинга, регуляции выживаемости и онкогенного ответа, имеющих значение для прогрессирования НМРЛ.

Приоритизация комплексов метаболит-мишень на основе молекулярного докинга

Молекулярный докинг был проведен для Соединения 4 и Соединения 5 с пятью белками, приоритетность которых была определена на основе топологии сети. Значения, приведенные в Таблице 1, представляют собой показатели докинга AutoDock Vina и не являются экспериментально измеренной свободной энергией связывания. Соединение 5 показало наиболее благоприятный единичный результат с MMP9 (-8.418 kcal·mol⁻¹), за ним следовало Соединение 4 с MMP9 (-7.840 kcal·mol⁻¹). Соединение 5 также имело более благоприятные показатели, чем Соединение 4, для AKT1 и EGFR, тогда как Соединение 4 показало несколько более благоприятные результаты для SRC (-6.549 против -6.204 kcal·mol⁻¹) и TNF (-5.436 против -5.299 kcal·mol⁻¹). Таким образом, Соединение 5 не продемонстрировало единообразного преимущества по показателям для всех пяти белков. Результаты докинга использовались только для отбора репрезентативных комплексов для дальнейшего структурного анализа.

Репрезентативные конформации докинга показаны на рисунке 5. На карте взаимодействий MMP9–Соединение 5, занявшей первое место в рейтинге, отображаемые контакты располагались вблизи Ala417 и Pro421 с расстояниями примерно 3,0 Å и 2,4 Å соответственно. На сохраненной карте взаимодействий для лучшей позы не было отмечено прямой координации Соединения 5 с ZN²⁺, и не было показано прямых контактов с His401, Glu402, His405 или His411. Это контрастирует с кристаллографическим эталоном 1GKC, в котором His401, His405 и His411 координируют каталитический ZN²⁺ на расстояниях 2,21, 2,23 и 2,22 Å, а реверсивный гидроксаматный ингибитор NFH координирует тот же ZN²⁺ через два атома кислорода на расстояниях 2,07 и 2,38 Å. Поскольку прямая координация Соединения 5-ZN²⁺ не была отмечена на сохраненной карте взаимодействий, расстояние координации Соединения 5-ZN²⁺ не было назначено; это интерпретируется как отсутствие продемонстрированной прямой координации на сохраненной карте, а не как измеренное значение разделения металлов. Трехмерное сравнение с кристаллографическим эталоном MMP9–NFH показано на дополнительном рисунке 1. Следовательно, геометрия отличается от канонического цинк-зависимого ингибирующего режима связывания, и настоящий результат докинга не подтверждает классификацию Соединения 5 как ингибитора MMP9. MMP9 была сохранена для МД-анализа только для того, чтобы определить, сохранится ли эта специфическая неканоническая геометрия докинга в ходе траектории. Для других комплексов Соединение 5 показало прогнозируемые контакты с AKT1 и SRC, в то время как Соединение 4 также образовало определенные взаимодействия при докинге с MMP9 и SRC. В соответствии с рисунком 5F, данные наблюдения описывают прогнозируемые взаимодействия и относительные баллы докинга, а не экспериментально подтвержденное сродство.

Для сравнительного, а не подтверждающего анализа методом молекулярной динамики (МД) были выбраны три комплекса. Комплекс MMP9–Compound 5 был выбран из-за наиболее благоприятного индивидуального значения скоринга при докинге, однако он имел неканоническую позу MMP9, что требовало тщательного структурного анализа. В качестве второго комплекса-кандидата был выбран SRC–Compound 5, а SRC–Compound 6 был включен в качестве соответствующей контрольной траектории с неловушечным соединением. Соответствующие начальные конформации показаны на Рисунок 6A, Рисунок 6B, и Рисунок 6CТакой дизайн позволил сравнить устойчивость выбранных геометрий докинга, не рассматривая стабильность молекулярной динамики как доказательство связывания с мишенью или функциональной регуляции.

Анализ молекулярной динамики

Для сравнения стабильности выбранных геометрий, полученных в результате докинга, в динамических водных условиях были созданы траектории молекулярной динамики (МД) продолжительностью 150 ns для комплексов MMP9–соединение 5 и SRC (PDB 2H8H)–соединение 5, при этом траектория SRC (PDB 2H8H)–соединение 6 была включена в качестве соответствующего отрицательного контроля. Начальные конформации представлены на рисунке 6A, рисунке 6B, и рисунке 6C. В ходе этих траекторий соединение 5 продемонстрировало более низкое значение RMSD лиганда (рисунок 6D>) в системах MMP9 и SRC по сравнению с соединением 6 в SRC. Траектория MMP9–соединение 5 достигла режима сравнительно низких флуктуаций, траектория SRC–соединение 5 стабилизировалась после начального периода адаптации, а траектория SRC–соединение 6 характеризовалась более выраженными флуктуациями. Эти различия указывают на большую стабильность выбранных поз докинга соединения 5 во время МД. Они не доказывают, что соединение 5 связывается с MMP9 или SRC в клетках, а траектория MMP9 не компенсирует отсутствие канонического взаимодействия с каталитическим-ZN²⁺ в исходной позе.

RMSD остова белка продемонстрировал схожую сравнительную динамику. Траектория MMP9–Compound 5 entered в относительно стабильный режим остова примерно через 30 ns, в то время как траектория SRC–Compound 5 вышла на плато позже, а траектория SRC–Compound 6 характеризовалась более выраженными колебаниями. Эти наблюдения описывают только поведение траекторий. Стабильность остова белка или траектории лиганда не может служить доказательством внутриклеточного связывания с мишенью, ингибирования фермента или модуляции передачи сигнала. Профили RMSD остова белка для всех трех систем представлены на Дополнительном рисунке 2.

Анализ траекторий

Для описания устойчивости контактов в траекториях МД использовались анализы заполненности водородных связей и флуктуаций остатков (Рисунок 7A). Соединение 5 продемонстрировало высокую заполненность водородной связи с Arg95 в траектории MMP9 (>85%) и повторяющееся взаимодействие с Leu325 в траектории SRC (>70%), тогда как репрезентативные контакты в контроле SRC–Соединение 6 имели более низкую заполненность. Эти остатки представлены не как доказательство функциональной модуляции мишени; значения заполненности лишь указывают на частоту возникновения указанных контактов в ходе анализируемых траекторий.

Анализ RMSF остатков в кармане связывания выявил специфические для каждой системы различия в локальной гибкости (Рисунок 7B). Траектория SRC–Соединение 6 продемонстрировала несколько более выраженных локальных флуктуаций по сравнению с траекторией SRC–Соединение 5, в то время как траектория MMP9–Соединение 5 показала сравнительно ограниченный профиль флуктуаций в наборе остатков собственного кармана связывания. Поскольку MMP9 и SRC являются разными белками, значения RMSF на уровне остатков не интерпретировались как прямое поостаточное сравнение. Вместе с данными о занятости водородных связей эти результаты характеризуют стабильность контактов и локальную гибкость, что помогает приоритизировать комплексы для экспериментального тестирования. Эти данные не определяют MMP9 или SRC в качестве внутриклеточных мишеней и не доказывают, что какой-либо из этих белков опосредует противораковый фенотип соединений.

Расчет свободной энергии связывания методом молекулярной механики и поверхности Пуассона-Больцмана (MM-PBSA)

Оценки MM-PBSA, рассчитанные на основе уравновешенных сегментов траекторий, представлены на Рисунке 8. Для комплекса MMP9–Compound 5 расчетное значение ΔG_bind составило -19.65 ± 6.43 kcal·mol⁻¹, для комплекса SRC–Compound 5 — -17.72 ± 6.84 kcal·mol⁻¹, а для контрольного комплекса SRC–Compound 6 — -10.37 ± 5.61 kcal·mol⁻¹. Таким образом, в рамках данного вычислительного протокола относительный энергетический рейтинг был следующим: MMP9–Compound 5, затем SRC–Compound 5 и SRC–Compound 6. Эти значения являются зависимыми от метода оценками, полученными из конечного сегмента траектории, и поправка на энтропию не применялась. Следовательно, они использовались только для внутриисследовательского сравнения и не должны интерпретироваться как экспериментально измеренные аффинности связывания или как доказательство функциональной модуляции белка.

На основании анализа докинга, МД, занимаемости контактов, флуктуации остатков и MM-PBSA, соединение 5 было выбрано в качестве приоритетного для дальнейшего изучения в комплексах с MMP9 и SRC. Согласованность этих расчетов обосновывает выбор данных пар для последующих экспериментов, однако она не подтверждает, что MMP9 или SRC являются прямыми внутриклеточными мишенями. В частности, неканоническая поза MMP9 и отсутствие продемонстрированной координации с каталитическим ZN²⁺ не позволяют сделать вывод о каноническом механизме ингибирования MMP9 на основе представленных структурных данных.

Выводы, сделанные на основании результатов

В ходе вычислительного анализа были приоритезированы 21 перекрывающихся мишени, связанных с заболеванием, и определены AKT1, EGFR, TNF, MMP9 и SRC как кандидаты в топологические хабы. Анализы на основе структуры позволили дополнительно выделить соединение 5 (Compound 5) для последующих экспериментов с комплексами MMP9 и SRC. Данные результаты не демонстрируют прямого связывания с мишенью, ингибирования MMP9 или SRC, регуляции сигнальных путей или причинно-следственной механистической связи между белковыми взаимодействиями и ранее установленным эффектом захвата апуринированных/апиримидинированных сайтов. Таким образом, исследование подтверждает совокупность проверяемых вычислительных гипотез, а не экспериментально установленный мультимишеневой механизм борьбы с НМРЛ (NSCLC).

ДОСТУПНОСТЬ ДАННЫХ:

Набор данных, подтверждающий результаты данного исследования, находится в открытом доступе: Wang X, Peng Z, Xing Y, Xue L. In silico prioritization of potential protein interactions for glutathione-responsive abasic site-trapping prodrugs in non-small cell lung cancer [dataset]. Figshare; 2026. doi:10.6084/m9.figshare.33313620.v1.

Рисунок 1: Химические структуры и взаимосвязи превращения исследуемых соединений, запускаемого глутатионом. Соединение 1 и соединение 2 представляют собой чувствительные к глутатиону пролекарства, которые высвобождают содержащие аминокси-группу метаболиты — соединение 4 и соединение 5 соответственно. Соединение 3 является соответствующим чувствительным к глутатиону структурным контролем, который образует соединение 6, лишенное функциональности по захвату абезитных сайтов с помощью аминокси-группы. Соединения 1–3 использовались для обратного предсказания мишеней, соединения 4 и 5 — для докинга с узловыми мишенями, а соединение 6 — в качестве лиганда отрицательного контроля при сравнении молекулярной динамики SRC. Сокращения: SRC, протоонкоген SRC, нерецепторная тирозинкиназа. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 2: Пересечение между предсказанными мишенями соединений и мишенями, связанными с немелкоклеточным раком легкого. (A) Сеть «соединение-мишень», созданная на основе результатов обратного предсказания мишеней для соединений 1–3. (B) Диаграмма Венна, показывающая перекрытие между предсказанными мишенями соединений и мишенями, связанными с заболеванием, которые были получены с помощью запроса по раку легкого H1299. 21 пересекающаяся мишень была отобрана для анализа белок-белковых взаимодействий, анализа обогащения и последующей приоритезации на основе структуры. Нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 3: Сеть белок-белковых взаимодействий и скрининг основных мишеней. (A) Сеть белок-белковых взаимодействий 21 пересекающейся мишени (116 ребер). (B) Первый этап скрининга с использованием критерия Degree ≥ 12, в результате которого было отобрано 13 кандидатов. (C) Второй этап скрининга 13 отобранных кандидатов с использованием показателей betweenness centrality ≥ 0,031293 и Closeness centrality > 0,769231, что позволило выделить пять кандидатов в хабы. (D) Итоговые пять кандидатов в хабы: AKT1, EGFR, TNF, MMP9 и SRC. Значения центральности, использованные для последовательных фильтров, были рассчитаны для исходной сети из 21 узла и 116 ребер и переносились без пересчета после формирования каждого подмножества. Сокращения: AKT1, серин/треониновая киназа AKT 1; EGFR, рецептор эпидермального фактора роста; TNF, фактор некроза опухоли; MMP9, матриксная металлопротеиназа 9; SRC, протоонкоген SRC, нерецепторная тирозинкиназа. Нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 4: Функциональный анализ обогащения пересекающихся мишеней. (A) Пузырьковая диаграмма 20 наиболее обогащенных путей Kyoto Encyclopedia of Genes and Genomes. (B) Столбчатая диаграмма 10 наиболее обогащенных терминов биологических процессов Gene Ontology. (C) Столбчатая диаграмма 10 наиболее обогащенных терминов клеточных компонентов Gene Ontology. (D) Столбчатая диаграмма 10 наиболее обогащенных терминов молекулярных функций Gene Ontology. Степень обогащения (fold enrichment) показана на оси x; размер пузырька на панели A отражает количество генов. Все лежащие в основе записи KEGG и GO, включенные в Дополнительные таблицы 4–7, соответствовали номинальному критерию включения p < 0.10; на графиках представлены пути/термины с наивысшим рангом по номинальному значению p. Значения, скорректированные с учетом множественного тестирования, приведены в дополнительных таблицах, но не использовались для включения. Нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 5: Позы молекулярного докинга и показатели AutoDock Vina для Соединения 4 и Соединения 5 с приоритетными белками. (A) Прогнозируемая поза докинга Соединения 5 с AKT1. (B) Прогнозируемая поза докинга Соединения 4 с MMP9. (C) Лучшая по рангу прогнозируемая поза Соединения 5 с MMP9; отображенные контакты находятся рядом с Ala417 и Pro421, при этом прямая координация каталитического ZN²⁺ или прямой контакт с His401, Glu402, His405 или His411 не отмечены. Таким образом, данная поза не представлена как канонический ингибиторный режим связывания с MMP9. Сравнение с контрольной картой взаимодействий кристаллической структуры MMP9, связанной с NFH (PDB 1GKC), включая кристаллографические расстояния координации ZN²⁺, приведено в Дополнительном рисунке 3. (D) Прогнозируемая поза докинга Соединения 4 с SRC. (E) Прогнозируемая поза докинга Соединения 5 с SRC. (F) Тепловая карта показателей докинга AutoDock Vina (kcal·mol⁻¹) для Соединения 4 и Соединения 5 в отношении пяти приоритетных белков. Более отрицательные значения указывают на более благоприятные показатели Vina в рамках данного протокола докинга; они не являются экспериментально измеренными значениями аффинности связывания. Сокращения: AKT1, серин/треониновая киназа 1 AKT; MMP9, матриксная металлопротеиназа 9; NFH, N²-[(2R)-2-{[формил(гидрокси)амино]метил}-4-метилпентаноил]-N,3-диметил-L-валинамид; SRC, протоонкоген SRC, нерецепторная тирозинкиназа; PDB, Protein Data Bank; Ala, аланин; Pro, пролин; His, гистидин; Glu, глутамат. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 6: Структурный обзор и анализ стабильности лигандов в комплексах МД. (A) Исходная конформация докинга соединения 5 с MMP9, использованная в качестве стартовой структуры для МД. (B) Исходная конформация докинга соединения 5 с SRC (PDB 2H8H). (C) Исходная конформация докинга соединения 6 с SRC (PDB 2H8H); соединение 6 является продуктом расщепления глутатионом контрольного соединения 3 и не обладает аминоокси-функциональностью для захвата абазических сайтов. (D) Среднеквадратичное отклонение (RMSD) лиганда относительно исходной позы докинга на траекториях длительностью 150 ns для комплексов MMP9–соединение 5, SRC–соединение 5 и SRC–соединение 6. На панели сравнивается стабильность позы во время МД; данные не демонстрируют связывание с мишенью внутри клеток. Сокращения: MD — молекулярная динамика; MMP9 — матриксная металлопротеиназа 9; SRC — протоонкоген SRC, нерецепторная тирозинкиназа; PDB — Protein Data Bank; RMSD — среднеквадратичное отклонение. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 7: Характеристики динамического взаимодействия в ходе молекулярно-динамического анализа. (A) Заполненность репрезентативных водородных связей между лигандом и белком в течение траекторий длительностью 150 ns для MMP9–Compound 5, SRC–Compound 5 и SRC–Compound 6. (B) Среднеквадратичное отклонение (RMSF) остатков связывающего кармана. Профиль MMP9 интерпретируется в рамках системы MMP9, тогда как профили SRC–Compound 5 и SRC–Compound 6 обеспечивают прямое сопоставление внутри SRC. Данные анализы описывают устойчивость контактов и локальную гибкость в ходе MD и не демонстрируют внутриклеточное взаимодействие с мишенью или функциональную модуляцию MMP9 или SRC. Сокращения: MD, молекулярная динамика; MMP9, матриксная металлопротеиназа 9; SRC, протоонкоген SRC, нерецепторная тирозинкиназа; RMSF, среднеквадратичное отклонение. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Рисунок 8: Оценки энергии методом молекулярной механики и поверхности Пуассона-Больцмана (MM-PBSA) для анализируемых комплексов. Рассчитанные значения ΔG_bind, полученные методом молекулярной механики и поверхности Пуассона-Больцмана на основе сегментов уравновешенных траекторий комплексов MMP9–Соединение 5, SRC (PDB 2H8H)–Соединение 5 и SRC (PDB 2H8H)–Соединение 6. Значения представлены как среднее ± стандартное отклонение и используются для относительного сравнения в рамках исследования, а не в качестве экспериментально измеренных показателей сродства к связыванию. На Рисунке 8 используется обозначение оси y ΔG_bind (kcal·mol⁻¹), что соответствует уравнению и терминологии, использованной в разделах «Методы» и «Результаты». Сокращения: MMP9 — матриксная металлопротеиназа 9; SRC — протоонкоген SRC, нерецепторная тирозинкиназа; PDB — Protein Data Bank; ΔG_bind — свободная энергия связывания; MM-PBSA — молекулярная механика и поверхность Пуассона-Больцмана. Пожалуйста, нажмите здесь, чтобы просмотреть увеличенную версию этого рисунка.

Соединение	AKT1 (ккал·моль⁻¹)	EGFR (ккал·моль⁻¹)	MMP9 (ккал·моль⁻¹)	SRC (ккал·моль⁻¹)	ФНО (ккал·моль⁻¹)
Соединение 4	-5.658	-4.913	-7.840	-6.549	-5.436
Соединение 5	-5.960	-5.188	-8.418	-6.204	-5.299

Таблица 1: Показатели докинга соединений 4 и 5 с пятью приоритетными белками в AutoDock Vina. Показатели докинга (kcal·mol⁻¹) соединений 4 и 5 с AKT1, EGFR, MMP9, SRC и TNF в AutoDock Vina. Более отрицательные значения указывают на более благоприятные показатели в рамках заданного протокола докинга. Данные значения являются расчетными и не должны описываться как экспериментально измеренные свободные энергии связывания или аффинности. Сокращения: AKT1 — серин/треонинкиназа 1 AKT; EGFR — рецептор эпидермального фактора роста; MMP9 — матриксная металлопротеиназа 9; SRC — протоонкоген SRC, нерецепторная тирозинкиназа; TNF — фактор некроза опухоли.

Дополнительный рисунок 1: Структурное сравнение кристаллографического эталона MMP9–NFH и наиболее высокоранжированной позы докинга MMP9–Соединение 5. (A) Окружение каталитического ZN²⁺ в эталонном комплексе MMP9–NFH (PDB 1GKC) с указанием His401, His405, His411, Glu402 и приведенных расстояний координации NFH. (B) Наиболее высокоранжированная поза докинга Соединения 5 с указанием контактов с Ala417 и Pro421. (C) Альтернативный трехмерный вид той же позы докинга Соединения 5. Данное сравнение приведено в качестве структурного эталона и не доказывает ингибирование MMP9. Сокращения: MMP9, матриксная металлопротеиназа 9; NFH, N²-[(2R)-2-{[формил(гидрокси)амино]метил}-4-метилпентаноил]-N,3-диметил-L-валинамид; PDB, Protein Data Bank; Ala, аланин; Pro, пролин; His, гистидин; Glu, глутамат.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

Дополнительный рисунок 2: Среднеквадратичное отклонение остова белка в ходе анализа молекулярной динамики. Профили среднеквадратичного отклонения остова белка для систем MMP9–Соединение 5, SRC (PDB 2H8H)–Соединение 5 и SRC (PDB 2H8H)–Соединение 6 на протяжении всех траекторий молекулярной динамики. На финальном графике используются стандартизированные обозначения MMP9–Соединение 5, SRC–Соединение 5 и SRC–Соединение 6; оси обозначены как RMSD (nm) и Time (ns). Данные профили описывают временную конформационную зависимость в ходе MD и не должны интерпретироваться как доказательство клеточного связывания или регуляции белка. Сокращения: MMP9, матриксная металлопротеиназа 9; SRC, протоонкоген SRC, нерецепторная тирозинкиназа; PDB, Protein Data Bank; RMSD, среднеквадратичное отклонение; MD, молекулярная динамика.Пожалуйста, нажмите здесь для загрузки файла.

Дополнительный рисунок 3: Сравнение окружения каталитического ZN²⁺ в эталонном комплексе 1GKC-NFH и в наиболее ранжированной позе докинга Соединения 5. В кристаллографическом эталонном комплексе MMP9–NFH (PDB 1GKC) His401, His405 и His411 координируют каталитический ZN²⁺ на расстояниях 2.21, 2.23 и 2.22 Å соответственно, а два атома кислорода NFH координируют ZN²⁺ на расстояниях 2.07 и 2.38 Å; Glu402 является каталитическим кислотно-основным остатком. Напротив, карта взаимодействий для наиболее ранжированной позы Соединения 5 показывает контакты с Ala417 (3.0 Å) и Pro421 (2.4 Å), но не содержит аннотированной прямой координации ZN²⁺ или прямых контактов с His401, Glu402, His405 или His411. Соответственно, расстояние координации Соединение 5-ZN²⁺ не было определено. Данное сравнение подтверждает классификацию позы Соединения 5 как неканонической предсказанной ассоциации, а не как канонического цинк-зависимого ингибирующего режима связывания. Сокращения: MMP9, матриксная металлопротеиназа 9; PDB, Protein Data Bank; Ala, аланин; Pro, пролин; His, гистидин; Glu, глутамат.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

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

Дополнительная таблица 2: Топологические метрики исходной сети для 13 кандидатов, отобранных после скрининга по степени узла. Топологические параметры для 13 узлов, оставшихся после применения критерия степени к исходной сети из 21 узла и 116 ребер. Эти значения представляют собой метрики исходной сети (21 узел, 116 ребер), перенесенные на последующий этап фильтрации по центральности по посредничеству и центральности по близости; они не пересчитывались для подсети из 13 узлов.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

Дополнительная таблица 3: Топологические метрики исходной сети для пяти финальных кандидатов в хабы, отобранных после последовательной фильтрации. Топологические параметры исходной сети из 21 узла и 116 ребер для пяти финальных кандидатов в хабы — AKT1, EGFR, TNF, MMP9 и SRC, отобранных после последовательной фильтрации. Эти перенесенные значения служат только для приоритизации на основе сети и не являются метриками, пересчитанными для подсети из пяти узлов, а также не определяют данные белки как экспериментально подтвержденные мишени для лекарственных средств. Сокращения: AKT1 — серин/треониновая киназа 1 AKT; EGFR — рецептор эпидермального фактора роста; TNF — фактор некроза опухоли; MMP9 — матриксная металлопротеиназа 9; SRC — протоонкоген SRC, нерецепторная тирозинкиназа.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

Дополнительная таблица 4: Полные результаты анализа обогащения по базе данных Kyoto Encyclopedia of Genes and Genomes для 121 пути, соответствующих номинальному критерию включения p < 0.10. Полная статистика обогащения по Kyoto Encyclopedia of Genes and Genomes для всех 121 отобранного пути среди 21 пересекающихся мишени (номинальное p < 0.10), включая отношение генов (gene ratio), количество генов (gene counts), общее количество в списке (list totals), количество попаданий в популяции (population hits), общее количество в популяции (population totals), значения p, значения Бенини, кратность обогащения (fold enrichment), значения Бонферрони, частоту ложных открытий (false discovery rates) и значения точного теста Фишера. 20 наиболее значимых путей представлены на рисунке 4A. Критерием включения служило номинальное значение p; значения Бенини, Бонферрони и частота ложных открытий приведены для прозрачности и не использовались для определения отобранного набора.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

Дополнительная таблица 5: Полные результаты анализа обогащения биологических процессов Gene Ontology (177 терминов, соответствующих номинальному значению p < 0.10). Полная статистика обогащения для всех 177 отобранных терминов биологических процессов Gene Ontology (номинальное p < 0.10), включая отношение генов, количество генов, общий размер списка, количество попаданий в популяции, общий размер популяции, p-значение, значение Бенджамини, кратность обогащения, значение Бонферрони, ложноположительный уровень (FDR) и значение точного теста Фишера. 10 наиболее значимых терминов визуализированы на Рисунке 4B. Критерием включения служило номинальное p-значение; скорректированные значения приведены для прозрачности и не использовались для определения отобранного набора.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

Дополнительная таблица 6: Полные результаты обогащения компонентов клеточного раздела Gene Ontology (29 терминов, соответствующих номинальному значению p < 0,10). Полная статистика обогащения для всех 29 отобранных терминов компонентов клеточного раздела Gene Ontology (номинальное p < 0,10), включая отношение генов, количество генов, общее число в списке, попадания в популяции, общее число в популяции, значение p, значение Бенджамини, кратность обогащения, значение Бонферрони, ложноположительную частоту (FDR) и значение точного теста Фишера. 10 наиболее значимых терминов визуализированы на Рисунке 4C. Критерий включения определялся номинальным значением p; скорректированные значения приведены для прозрачности и не использовались для определения отобранного набора.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

Дополнительная таблица 7: Полные результаты анализа обогащения молекулярных функций Gene Ontology (61 термин, соответствующий номинальному значению p < 0,10). Полная статистика обогащения для всех 61 сохраненного термина молекулярной функции Gene Ontology (номинальное p < 0,10), включая соотношение генов, количество генов, общий объем списка, количество попаданий в популяции, общий объем популяции, значение p, значение Бениджамини, кратность обогащения, значение Бонферрони, ложноположительный результат (FDR) и значение точного теста Фишера. 10 наиболее значимых терминов визуализированы на рисунке 4D. Критерий номинального значения p определял включение; скорректированные значения приведены для прозрачности и не использовались для определения сохраненного набора.Пожалуйста, нажмите здесь, чтобы скачать этот файл.

Обсуждение

Результаты настоящего исследования следует интерпретировать в рамках двухуровневой механистической модели. Первый уровень опирается на предыдущие экспериментальные работы: активация в ответ на глутатион высвобождает аминоокси-содержащий агент, способный улавливать альдегидные абезитные сайты, а данный класс соединений продемонстрировал селективную цитотоксичность, влияние на клеточный цикл и индукцию апоптоза в клетках H1299. Второй уровень носит исследовательский характер и является основной целью данной работы: компьютерное прогнозирование мишеней и структурный анализ позволяют предположить, что та же серия химических соединений может взаимодействовать с определенными белками, ассоциированными с раком. Ни один из экспериментов в данном исследовании не доказывает, что предсказанные белковые взаимодействия происходят в клетках или что они являются следствием улавливания абезитных сайтов. Таким образом, эти два уровня намеренно разделены, а не объединены в доказанный мультитаргетный механизм.

На сетевом уровне 21 пересекающаяся мишень была обогащена в путях, связанных с раком и резистентностью, при этом AKT1, EGFR, TNF, MMP9 и SRC занимали центральные позиции в сети PPI. Эти результаты полезны для приоритизации кандидатов, однако центральность сети и обогащение не могут доказать, что соединение физически взаимодействует с белком или изменяет сигнальный путь. Следовательно, хабовые белки следует рассматривать как кандидатов для таргетной валидации. Для установления функциональной роли потребуются эксперименты по прямому воздействию или анализы взаимодействия с мишенью на той же клеточной линии NSCLC, которая использовалась для измерения клеточного фенотипа.

MMP9 иллюстрирует, почему это различие важно. Соединение 5 имело наиболее благоприятный единичный балл докинга с MMP9, однако на сохраненной карте взаимодействий для лучшей позы не была отмечена прямая координация каталитического иона Zn²⁺ или прямое взаимодействие с His401, Glu402, His405 или His411. В кристаллографическом эталоне 1GKC His401, His405 и His411 координируют Zn²⁺ на расстояниях 2.21, 2.23 и 2.22 Å, в то время как два атома кислорода NFH координируют Zn²⁺ на расстояниях 2.07 и 2.38 Å; Glu402 является каталитическим кислотно-основным остатком. Карта Соединения 5 вместо этого показала контакты рядом с Ala417 и Pro421. На Дополнительном рисунке 3 представлено сравнение с эталонной картой взаимодействий. Поскольку прямая координация Соединения 5-Zn²⁺ не была отмечена, расстояние координации Соединения 5-Zn²⁺ не было определено. Соответственно, текущий структурный результат лучше всего описывается как неканоническая предсказанная ассоциация с MMP9, а не как доказательство ингибирования каталитического центра. Последующая траектория МД лишь проверяет, остается ли эта конкретная поза докинга устойчивой с течением времени; она не может превратить неканоническую геометрию докинга в доказательство ферментативного ингибирования^36,37. Следовательно, MMP9 должна оставаться кандидатом для биохимического тестирования, а не основным или валидированным мишенью. MMP9 широко обсуждалась в связи с инвазией рака, ремоделированием микроокружения опухоли и терапевтическими стратегиями, направленными на MMP9^38,39.

Такие же ограничения в доказательной базе применимы и к SRC. Траектория SRC–Соединение 5 характеризовалась более стойкими контактами, чем контрольная траектория SRC–Соединение 6, однако стабильность молекулярной динамики (MD) не эквивалентна внутриклеточной степени заполнения мишени^40,41. Использование Соединения 6 обеспечивает сопоставимое сравнение для продукта расщепления, не обладающего способностью к захвату, и усиливает внутреннее структурное сравнение, однако оно не доказывает, что аминоокси-функциональность ответственна за связывание с SRC или что сигнализация SRC в клетках изменяется. Для подтверждения такого утверждения потребовалась бы прямая оценка общего уровня SRC, соотношения p-SRC/SRC, маркеров нижележащей сигнализации и ортогональные измерения взаимодействия с мишенью.

Метод MM-PBSA аналогичным образом дает относительно энергетическую оценку для исследуемых траекторий, а не является экспериментальным измерением аффинности. Более благоприятные оценки для MMP9–Compound 5 и SRC–Compound 5 по сравнению с SRC–Compound 6 согласуются с результатами сравнительного анализа траекторий, однако расчеты чувствительны к выбранным конформациям и методологическим приближениям, а настоящий анализ не включал поправку на энтропию. Таким образом, согласованность данных докинга, МД и MM-PBSA повышает лишь внутреннюю вычислительную согласованность; она не устанавливает функциональную модуляцию MMP9 или SRC.

Взаимосвязь между захватом апуриниковых/апиримидиновых сайтов и прогнозируемыми белковыми взаимодействиями остается невыясненной. Один из вариантов заключается в том, что аминоокси-каркас, высвобождаемый под действием глутатиона, сохраняет свою установленную активность по захвату повреждений ДНК, имея при этом независимые параллельные взаимодействия с определенными белками. Другой вариант состоит в том, что некоторые прогнозируемые белковые взаимодействия не происходят при биологически значимых концентрациях или не влияют на фенотип. Текущие данные не позволяют разграничить эти возможности. Для демонстрации механической связи потребуется одновременное измерение ответов на повреждение ДНК и изменений в белковых путях после обработки соединением с последующим проведением экспериментов по пертурбации, чтобы показать, что изменение потенциальной мишени изменяет противораковый фенотип.

Последующая экспериментальная валидация должна быть проведена на клетках H1299 — той же клеточной модели, которая использовалась при предыдущей экспериментальной характеристике данной серии химических соединений. Поэтапная стратегия должна начинаться со сравнения Соединения 5 с контролем (растворителем) и соответствующим неловушечным Соединением 6 с использованием анализов жизнеспособности (зависимость «концентрация-ответ»), апоптоза, клеточного цикла, миграции/инвазии и анализа γH2AX для определения фенотипического контекста и степени повреждения ДНК⁴². Затем следует исследовать MMP9 с помощью желатиновой зимографии и анализа экспрессии белков, а сигнализацию SRC — через определение общего уровня SRC и фосфорилированного по Tyr416 SRC, используя соотношение p-SRC/SRC в качестве основного показателя сигнализации. Прямое взаимодействие соединения с белком следует оценить независимо с помощью ортогонального метода анализа связывания с мишенью, такого как поверхностный плазмонный резонанс⁴³. Гипотеза о взаимодействии считается экспериментально подтвержденной только в том случае, если данные о прямом связывании согласуются с соответствующим функциональным клеточным ответом. Генетическое или фармакологическое воздействие на MMP9 или SRC позволит дополнительно проверить, вносит ли какой-либо из этих кандидатов причинный вклад в фенотип H1299. Такая поэтапная схема позволяет сохранить различие между ранее установленной активностью по захвату апуриновых/апиримидиновых сайтов и приоритетными с точки зрения компьютерного моделирования гипотезами о взаимодействии с MMP9/SRC, определяя при этом прямой путь для последующей экспериментальной проверки.

Материалы

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

Имя	Компания	Каталожный номер	Комментарии
ACPYPE	Разработчики ACPYPE / Bio2Byte	Версия 2023.11.14	Генерация и конвертация файлов топологии низкомолекулярных лигандов в формат, совместимый с GROMACS.
AmberTools (Antechamber)	Группа разработки AMBER	AmberTools 23.3	Присвоение частичных зарядов методом AM1-BCC и генерация параметров лиганда с помощью GAFF2.
AutoDock Vina	Лаборатория Форли, Scripps Research	Версия 1.2.5	Молекулярный докинг и генерация ранжированных поз белок-лигандных комплексов.
AutoDockTools (MGLTools)	Центр вычислительной структурной биологии, Scripps Research	Версия 1.5.7	Подготовка рецептора и лиганда, назначение зарядов Гастейгера, определение вращающихся связей и конвертация в формат PDBQT.
BIOVIA Discovery Studio Visualizer	Dassault Systèmesèmes BIOVIA	Версия 2025	Визуализация поз докинга и анализ взаимодействий белок-лиганд.
силовое поле CHARMM36	Лаборатория МакКерелла / разработчики силового поля CHARMM	CHARMM36	Параметризация белков для расчетов молекулярной динамики; на протяжении всего исследования последовательно использовалось силовое поле CHARMM36 (не CHARMM36m).
Chem3D Ultra (ChemOffice Professional)	Программное обеспечение Revvity Signals	Версия 22.2	Минимизация энергии лиганда с использованием силового поля MMFF94.
Соединения 1–6	Синтезировано ранее, как описано в работе Li et al., ACS Chemical Biology (2022)	Нет данных	Глутатион-чувствительные исходные пролекарства, продукты высвобождения и соответствующая пара незахватывающего контроля, использованные в вычислическом рабочем процессе.
cytoHubba	разработчики cytoHubba / Магазин приложений Cytoscape	Версия 0.1	Анализ степени, центральности по посредничеству и центральности по близости для приоритизации пар «хаб — мишень».
Cytoscape	Консорциум Cytoscape	Версия 3.10.0	Визуализация сетей «соединение-мишень» и белок-белковых взаимодействий (PPI), а также анализ их топологии.
Биоинформационные ресурсы DAVID	Лаборатория ретровирологии человека и иммуноинформатики, Национальная лаборатория по исследованию рака имени Фредерика	Веб-ресурс	Анализ обогащения Gene Ontology и KEGG; номинальное значение p &В качестве предварительного критерия включения использовалось значение < 0,10; скорректированные значения были представлены, но не использовались для определения итогового набора данных.
GAFF2	Команда разработчиков AMBER	GAFF2	Параметризация силового поля лиганда.
База данных генов человека GeneCards	Комплект GeneCards / LifeMap Sciences, Inc. / Институт науки Вейцмана	Веб-ресурс	Поиск мишеней, ассоциированных с заболеванием, с использованием запроса «lung cancer H1299».
gmx_MMPBSA	Команда разработчиков gmx_MMPBSA	Версия 1.6.3	Расчеты свободной энергии связывания методом MM-PBSA на основе траекторий молекулярной динамики GROMACS; из 100,0 нс было отобрано 500 кадров.–149,9 нс с интервалом 100 пс, за исключением конечной точки 150,0 нс; коррекция энтропии не применялась.
GROMACS	Команда разработчиков GROMACS	Версия 2024.4	Генерация и анализ траекторий молекулярной динамики.
Онлайн-платформа визуализации Microbioinfo	Shanghai Newcore Biotechnology Co., Ltd. / Microbioinfo	Веб-ресурс	Построение столбчатых диаграмм и пузырьковых графиков обогащения.
Open Babel	Команда разработчиков Open Babel	Версия 3.1.1	Преобразование файлов структур в канонические SMILES для ввода в SwissTargetPrediction.
PubChem	Национальный центр биотехнологической информации, Национальная медицинская библиотека США, NIH	Веб-ресурс	Получение трехмерных структур лигандов.
Система молекулярной графики PyMOL	Шр.ödinger, LLC	Версия 2.5.4	Подготовка белка, визуализация структуры и анализ поз докинга.
RCSB Protein Data Bank	Исследовательская коллаборатория по структурной биоинформатике (RCSB)	ID PDB: 3O96; 5UWD; 2AZ5; 1GKC; 2H8H	Поиск структур белков AKT1, EGFR, TNF, MMP9 и SRC соответственно.
Пожалуйста, предоставьте текст для перевода.	Консорциум STRING	Версия 11.5	Построение сети белок-белковых взаимодействий; Homo sapiens; минимальный требуемый балл взаимодействия 0,400; исходные метрики топологии соответствуют ненаправленной сети белок-белковых взаимодействий, состоящей из 21 узла и 116 ребер.
SwissTargetPrediction	Группа молекулярного моделирования, Лозаннский университет / SIB Швейцарский институт биоинформатики	Интернет-ресурс	Обратное предсказание мишеней для соединений 1–3; виды, ограниченные Homo sapiens.
модель воды TIP3P	Реализовано в GROMACS	TIP3P	Трехточечная модель воды с явным представлением растворителя, используемая для сольватации белок-лигандных комплексов.
Диаграммы Венна	BioinfoGP, Национальный центр биотехнологииíа (CNB-CSIC)	Версия 2.1	Пересечение списков мишеней, предсказанных для соединения, и мишеней, ассоциированных с заболеванием.