Скачать данные
Данные экспрессии генов
Данные секвенирования одноклеточных РНК (scRNA-seq), использованные в настоящем исследовании, были получены из репозитория Gene Expression Omnibus (GEO), поддерживаемого Национальным центром биотехнологической информации (NCBI) (https://www.ncbi.nlm.nih.gov/geo/), а именно из набора данных с номером9 GSE161470 (человеческая сердечная ткань, состоящая из четырёх контрольных образцов и одного патологического образца). Этот набор данных был первоначально опубликован Чжаном и др. в 2022году 10. Основной целью первоначального исследования было изучение клеточной гетерогенности и молекулярных регуляторных механизмов в сердечной ткани человека при сердечной недостаточности. Для текущего анализа были отобраны пять образцов из этого набора данных, каждый из которых содержит комплексные профили экспрессии отдельных клеток, полученные из сердечной ткани человека. Другой набор данных, использованный в этом исследовании, также был получен из публичного репозитория NCBI GEO, в частности матричный файл серии, соответствующий номеру доступа GSE161472, сопровождаемый аннотационным файлом GPL11154. Профиль экспрессии включает в общей сложность 84 образца, включая 37 контрольных образцов и 47 образцов заболеваний. Это исследование включает интегративный мультиомический анализ, при этом все исследования проводятся с использованием общедоступных данных.
Данные eQTL
Данные eQTL, полученные от консорциума eQTLGen, сосредоточены на разъяснении генетической структуры экспрессии генов в крови и генетических факторов, влияющих на сложныепризнаки 3. В настоящее время консорциум находится на втором этапе своего масштабного проекта, проводя мета-анализ общегеномных данных, касающихся экспрессии генов крови.
Данные экспозиции — mQTL
Данные mQTL были получены из опубликованного метаанализа европейской (EUR) когорты, который изучает метилирование ДНК цельной крови в рамках генетической структуры 3 701 образца из популяций европейскогопроисхождения 11. Набор данных содержал информацию о 426 636 признаках mQTL.
Данные экспозиции — pQTL
Плазменные pQTL-данные были получены из базы данных deCODE (https://www.decode.com/summarydata/)4. В этом исследовании использовались данные 2021 года по набору данных deCODE pQTL, включавшему общегеномное ассоциативное исследование (GWAS) уровней белков в плазме, измеряемые с использованием 4 907 аптамеров в когорте из 35 559 особей европейского происхождения.
Данные о результатах
Сводные статистические данные по сердечной недостаточности были получены в результате масштабного геномного ассоциативного исследования (GWAS), в котором в основном участвовали участники европейского происхождения, доступного через базу данных Европейского института биоинформатики (EBI) (GCST90162626). Набор данных по сердечной недостаточности включал 115 150 случаев и 1 550 331 контрольную группу. Каталог GWAS, включающий публикации, ведущие ассоциации и подробную сводную статистику, в настоящее время предлагает данные, сопоставленные с Genome Assembly и dbSNP Build.
Менделевский рандомизационный анализ mQTL, eQTL и pQTL
Для систематического изучения возможных причинно-следственных связей между экспрессией генов, содержанием белка, уровнем метилирования ДНК и риском сердечной недостаточности были проведены менделевские рандомизации (MR) с использованием количественных признаков экспрессии (eQTL), количественных белковых признаков (pQTL) и количественных признаков метилирования (mQTL). Во время предварительной обработки данных о воздействии однонуклеотидные полиморфизмы (SNP), связанные с каждой переменной воздействия (ген, белок или сайт метилирования), были извлечены из соответствующих баз данных при общегеномном пороге значимости P < 1 × 10⁻⁵ для использования в качестве начальных кандидатов инструментальных переменных (IV). Впоследствии для IVs каждого фактора воздействия было проведено скопление дисбаланса связи (LD) с использованием окна размера 10 000 килобаз (kb) и порога LD R² 0,001 для обеспечения независимости между приборами. Эти отобранные IV затем были согласованы с сводными статистическими данными, полученными в исследовании ассоциации по всему сердечной недостаточности (GWAS; ID: GCST90162626) путём использования функции read_outcome_data, сохраняя только те SNP с ассоциативным значением P ниже 5×10⁻⁵ в наборе результатов. Для снижения слабого смещения прибора статистика F для каждого IV была рассчитана как F = (β_exposure/SE_exposure)², и в последующие анализы были включены только инструменты с F > 10. Для оценки причинных эффектов выравнивание аллелей между наборами данных по воздействию и результатам выполнялось с помощью функции harmonize_data из пакета TwoSampleMR. Затем были проведены MR-анализы с использованием четырёх дополнительных статистических подходов: (1) метод обратного взвешенного дисперсии (IVW), который предоставляет мета-анализ оценок коэффициента Вальда между SNP; (2) регрессия MR-Эггера, которая учитывает направленную плейотропию путём включения перехватного члена в предположение Сила инструмента, независимая от прямого эффекта (InSIDE); (3) взвешенный медианный метод, который даёт последовательные причинные оценки даже если до 50% инструментов недействительны; и (4) метод взвешенного режима, который определяет наиболее частый кластер оценки причинно-следственных эффектов, обеспечивая повышенную статистическую мощность и сниженную ошибку типа I относительно MR-Egger. В случаях, когда была доступна только одна инструментальная переменная, применялся исключительно метод коэффициента Вальда. Для оценки надёжности результатов были проведены комплексные анализы чувствительности, включая тестирование на гетерогенность через mr_heterogeneity-функцию, оценку плейотропии с использованием mr_pleiotropy_test и анализы с исключением, реализованные через функцию mr_leaveoneout, которая итеративно исключает каждый SNP для определения влияния отдельных вариантов на общие результаты. Значимые ассоциации визуализировались с помощью графических инструментов, таких как mr_scatter_plot и mr_forest_plot. Этот аналитический конвейер был равномерно применён во всех наборах данных eQTL, pQTL и mQTL для поддержания методологической согласованности на протяжении всего исследования.
Анализ колокализации
Был проведён колокализационный анализ с помощью метода Колока, сводных данных eQTL и GWAS сердечнойнедостаточности 5. Индексный однонуклеотидный полиморфизм (SNP) использовался для вычисления апостериорной вероятности в окне скопления размером 100 кб. В колокализационном (колокальном) анализе гипотеза H3 обозначает постериорную вероятность того, что два признака — экспрессия генов и сердечная недостаточность — коррелированы, но имеют различные причинные варианты. Напротив, гипотеза H4 указывает на постериорную вероятность того, что связь между двумя признаками обусловлена одним общим причинным вариантом. Порог SNP. PP. Для определения колокализации использовалась H4 больше 0,90.
Иммунная инфильтрация
Метод CIBERSORT широко применяется для оценки типов иммунных клеток вмикросреде. Используя принципы регрессии опорных векторов, можно провести анализ деконволюции экспрессионной матрицы подтипов иммунных клеток. Включая 547 биомаркеров, CIBERSORT может дифференцировать 22 фенотипа иммунных клеток человека, включая Т-клетки, В-клетки, плазматические клетки и различные субпопуляции миелоидных клеток. Используя набор данных GSE161472, был проведён анализ с использованием алгоритма CIBERSORT в сочетании с интегрированной матрицей сигнатуры LM22, которая характеризует профили экспрессии генов 22 различных типов человеческих иммунных клеток. Уровни инфильтрации этих 22 популяций иммунных клеток были количественно оценены для каждого отдельного образца. Впоследствии функция cor.test была применена для оценки корреляций между экспрессией ключевых генов и соответствующими уровнями инфильтрации иммунных клеток.
Обработка данных и контроль качества секвенирования одиночных клеток РНК
Данные профиля выражения отдельных ячейок были обработаны с помощью пакета Seurat (V4.3.0) в среде R (V4.3.0) 6. В этом исследовании использовался традиционный рабочий процесс для анализа данных секвенирования одноклеточных РНК. Изначально профили выражений импортировались с использованием пакета Seurat. Клетки фильтровались на основе нескольких показателей качества, включая общее количество UMI в каждой клетке, количество экспрессированных генов, долю митохондриальных показаний и долю рибосомных показаний. Выбросы были выявлены как значения, отклоняющиеся от медианы более чем на три медианные абсолютные отклонения (MAD). Конкретные пороги фильтрации были следующими: nFeature_RNA ≥ 200, percent.mt ≤ 2.15718, nFeature_RNA ≤ 2840.941 и nCount_RNA ≤ 5194.27. Обычно клетки с чрезмерно высоким общим количеством UMI и количеством экспрессированных генов классифицировались как дублеты, тогда как клетки с повышенным процентом митохондриальных или рибосомных показателей считались низкого качества, потенциально подвергающимися апоптозу или фрагментации. После этих этапов фильтрации был использован DoubletFinder (версия 2.0.4) для индивидуального выявления и удаления дублетов из каждого образца, завершив процесс контроля качества ячеек. Изначально нормализация данных проводилась с помощью функции normalizeData. Состояние клеточного цикла впоследствии оценивалось с помощью функции CellCycleScoring, и с помощью метода FindVariableFeatures были выявлены высоковариабельные гены. Затем набор данных был масштабирован с помощью ScaleData для стандартизации данных и снижения влияния митохондриальных генов, рибосомных генов и клеточных циклов на последующие анализы. Линейное уменьшение размерности проводилось с помощью анализа главных компонентов (PCA) через функцию RunPCA, при этом значительные основные компоненты были выбраны для дальнейшего анализа. Для решения пакетных эффектов между различными выборками был использован алгоритм Harmony (версия 1.1.0). Этот подход итеративно кластеризует похожие клетки из отдельных партий в пространстве PCA, сохраняя при этом разнообразие партий внутри кластеров. Учитывая относительно слабые пакетные эффекты, наблюдаемые в наборе данных, были применены стандартные параметры (θ = 2). Впоследствии было проведено нелинейное уменьшение размерности с помощью RunUMAP, затем построено граф окрестностей ячеек с помощью FindNeighbors и кластеризация ячеек через FindClusters. Для аннотаций типов клеток была реализована иерархическая структура аннотирования: первичная ручная аннотация основывалась на характерных паттернах экспрессии генов, основанных на базе данных CellMarker и соответствующей литературе; это дополнялось автоматическими результатами аннотаций, полученными из программного обеспечения SingleR в качестве справочника. Для дальнейшего повышения точности и полноты идентификации типов клеток были проконсультированы несколько авторитетных баз данных, включая Human Primary Cell Atlas (HPCA), BlueprintEncode, MonacoImmune, DatabaseImmuneCell и NovershternHaematopoetic. Аннотация ячеек проводилась путем запроса к базе данных CellMarker (http://117.50.127.228/CellMarker/CellMarkerBrowse.jsp) и изучения литературы, при поддержке автоматической поддержки аннотации, предоставляемой программным обеспечением SingleR (V2.4.0) 13. Её цель — определить типы клеток, присутствующие в соответствующих тканях, и связанные с ними маркерныегены 14.
Анализ взаимодействия лиганда и рецепторов
В этом исследовании CellCall (версия 1.0.7) использовалась для проведения комплексного анализа межклеточных коммуникационныхсетей 15. Используя аннотации по типу клеток, полученные из Сёра, вместе с исходной матрицей счёта, был создан нормализованный аналитический объект с параметрами, настроенными для человеческого генома. Функция TransCommuProfile применялась для количественной оценки силы взаимодействий между клетками с помощью взвешенного алгоритма, реализуя порог значимости значения p < 0,05 для выявления надёжных пар лиганд-рецептор. Значимые пары взаимодействия впоследствии подвергались анализу обогащения путей KEGG с помощью функции getHyperPathway, а взаимосвязи между типами клеток и путями были проиллюстрированы с помощью пузырьковых графиков. Общая сеть связи в конечном итоге визуализировалась с помощью кругового графика, на котором восемь различных цветов обозначали разные типы клеток. Сила взаимодействия и направленность были представлены стрелочными особенностями, что дало детальную характеристику динамики межклеточной сигнализации.
Анализ псевдовремени
Для изучения динамической транскрипционной регуляции макрофагов на протяжении всего прогрессирования сердечной недостаточности в этом исследовании использовался алгоритм Monocle для проведения псевдовременного анализа макрофагов в субпопуляциях. Матрица экспрессии генов, соответствующая целевой субпопуляции клеток, была извлечена для построения объектов анализа траектории отдельных клеток, при этом в качестве упорядоченных признаков выбирались высоковариабельные гены. Используя метод уменьшения размерности DDRTree, ячейки отображались в двумерном пространстве для реконструкции траектории дифференцировки. Были проведены визуализационные анализы для определения распределения клеток вдоль псевдооси времени и выявления генов, экспрессия которых значительно изменилась за псевдовремя. Последующие анализы были сосредоточены на ключевом гене DBNL и охарактеризовали динамику его экспрессии вдоль клеточной траектории, выяснив транскрипционные механизмы перепрограммирования макрофагов при прогрессированиисердечной недостаточности 16.
Анализ обогащения множества генов (GSEA)
В этом исследовании был использован подход GSEA для выяснения регуляторных механизмов, связанных с ключевыми генами, участвующими в сердечной недостаточности. Используя ранее выявленные ключевые гены, образцы были стратифицированы на когорты с высокой и низкой экспрессией на основе медианного значения экспрессии. Дифференциальный анализ экспрессии проводился с использованием лимма-пакета, формируя ранжированный список генов по логическим изменениям складки (logFC). Последующий анализ обогащения путей KEGG проводился с помощью инструмента clusterProfiler, при этом наборы генов, полученные из базы данных MsigDB, служили эталонным фоном. Затем алгоритм GSEA применялся для выявления сигнальных путей, значительно обогащённых между двумя группами выражения, а для определения статистической значимости использовался скорректированный порог p-значения менее 0,05. Для иллюстрации регуляторных функций основных генов в критически важных путях использовались различные методы визуализации, включая многопутевые графики GSEA и круговые сетевые диаграммы.
Анализ вариации генных множеств (GSVA)
GSVA — это непараметрический, неконтролируемый подход, используемый для оценки обогащения генных наборов в транскриптомных данных. Этот метод преобразует вариации на уровне генов в вариации на уровне пути, рассчитывая составные оценки для конкретных наборов генов, что облегчает оценку биологических функциональных изменений в различных образцах. В настоящем исследовании наборы генов были получены из базы данных молекулярных подписей. Алгоритм GSVA был использован для вычисления составных оценок для каждого набора генов, что позволило оценивать потенциальные биологические функциональные изменения на различных образцах. Результаты анализа обогащения GSVA приведены в дополнительном материале (Дополнительная таблица 1).
Прогнозирование препаратов CTD
Целевой ген (DBNL) был введён в поле поиска Сравнительной токсикогеномной базы данных (CTD), выбрана категория заболевания «сердечно-сосудистое заболевание», и был выполнен запрос для получения данных прогнозирования лекарств, связанных с состоянием «сердечная недостаточность». Полученные результаты прогнозирования впоследствии импортировались в программное обеспечение Cytoscape для облегчения визуализации данных и создания карты сети взаимодействия ген-химических.
Методы молекулярного стыковки
Из-за неразрешённой трёхмерной кристаллической структуры человеческого белка DBNL (UniProt ID: Q9UJU6) это исследование предсказало трёхмерную структуру DBNL на основе AlphaFold317. Пириниксиновая кислота (WY-14643) доступна для скачивания из базы данных PubChem (PubChem CID: 5694). Затем структура белка была предварительно обработана с помощью программного обеспечения MGLTools (версия 1.5.7)18, включая такие этапы, как добавление атомов водорода. В то же время белки и малые молекулы были преобразованы в формат PDBQT, необходимый для стыковки. Программное обеспечение AutoDock Vina (версия 1.1.2)19 использовалось для глобального молекулярного стыковки (экзогенность=16, num_modes=30) с целью изучения возможных режимов связывания. После завершения расчётов стыковки для последующих молекулярно-динамических моделей должна быть выбрана комплексная конформация с наибольшей аффинностью, обозначенная минимальной свободной энергией связывания.
Метод моделирования молекулярной динамики
Для систематического изучения стабильности связывания и механизмов взаимодействия между кандидатами и белками были проведены традиционные молекулярнодинамические симуляции (MD) с использованием программного пакета GROMACS (версия 2024.03)20. Параметры белка были сгенерированы с помощью силового поля Amber14SB21, модель молекулы воды — с помощью моделиTIP3P 22, а параметры топологии лиганда — с помощью инструмента Antechamber Python Parser Interface (ACPYPE), основанного на поле общего янтаря (GAFF). Система комплекса лиганд-белок впоследствии размещалась в периодической граничной октаэдрической коробке, заполненной молекулами воды TIP3P. Были введены ионы натрия (Na⁺) и хлорида (Cl⁻) для достижения концентрации 0,15 моль/л и нейтрализации общего заряда системы. После завершения строительства системы первым шагом было минимизация энергии с использованием самого крутого спуска на 50 000 шагов с целью устранения потенциально необоснованных форм внутри конструкции. Впоследствии были проведены два этапа уравновесения системы: симуляция NVT с постоянной частотой частиц (постоянное количество частиц, объём и температура), а затем симуляция NPT с постоянной частотой частиц, давления и температуры на 100 пс. Во время этих моделей на тяжёлые атомы белкового основного основания применялись позиционные ограничения для сохранения структурной целостности белка. Температура поддерживалась на уровне 300 К с помощью термостата с V-образным масштабированием, а давление регулировалось на уровне 1 бар с помощью баростата Парринелло-Рахмана для сцепления давления. После завершения фазы уравновесия была проведена симуляция производственной фазы продолжительностью 100 наносекунд, в ходе которой все позиционные ограничения устраняются. Траектория была интегрирована с использованием двухфемтосекундного временного шага, а для точного управления дальнодействующими электростатическими взаимодействиями применялся метод Эвальда с частицами-сеткой (PME). Траектория сохранялась каждые 10 ps, а для последующего анализа выводилось всего 10000 кадров. Кроме того, из симуляции были извлечены стабильные траектории в интервале 90–100 нс, а свободная энергия связывания белковых комплексов лигандов рассчитана с помощью инструмента gmx MMPBSA23.
Статистический анализ
Достоверность этого менделевского рандомизационного анализа (МР) зависит от трёх фундаментальных предположений. (1) Релевантность: Инструментальные переменные (IV) должны иметь сильную связь с воздействием. (2) Независимость: IV должны быть независимы от любых факторов, влияющих как на воздействие, так и на результат. (3) Ограничение на исключение: IV должны влиять на результат исключительно через их влияние на воздействие. Нарушение этого предположения, при котором внутривенное введение влияет на результат через пути, не связанные с воздействием, называется горизонтальной плейотропией. Все статистические анализы проводились с использованием версии R 4.3.0 с двусторонними тестами, и значение p менее 0.05 обычно считалось статистическим показателем значимости.