ИСПОЛЬЗОВАНИЕ МЕТОДА ПОЛНОГЕНОМНОГО АССОЦИАТИВНОГО ИССЛЕДОВАНИЯ В СЕЛЕКЦИИ СВИНЕЙ

Ускорение генетического улучшения основных хозяйственно-полезных признаков, имеющих значение в секторе экономики, путем применения селекции с помощью маркеров повышает эффективность отбора и разведения свиней. Однако анализ генетической архитектуры этих признаков остается сложной и трудной задачей, и для любого признака у свиней было выявлено несколько определенных причинных вариантов – это миссенс-мутация гена рецептора меланокортина 4 (MC4R), которая участвует в энергетическом гомеостазе и соматическом росте

Геномная селекция, основанная на прогнозировании племенной ценности (геномной племенной ценности или GEBV) животных с учетом геномной информации, меняет стратегии и подходы к разведению молочного КРС и ряда других видов животных. Отправными точками для применения геномной селекции у свиней стали разработка первой коммерческой панели однонуклеотидного полиморфизма для высокопроизводительного генотипирования, секвенирование генома свиньи и применение статистических и методологических подходов, впервые разработанных для молочного скота

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

Материалом исследования являлись данные геномных оценок племенной ценности показателя среднесуточного прироста 408 хряках породы дюрок 2017-2020 годов рождения, проходивших откорм на автоматических кормовых станциях. В тексте встречается следующее сокращение признака среднесуточного прироста – ADG, г.

Геномный BLUP (GenomicBLUP, GBLUP) также идентичен по применяемой методологии решения уравнения обычному BLUP, однако матрица родства A заменяется матрицей геномного сходства G, учитывающей как генотип по каждому SNP, так и частоту этого генотипа:

где

С – «центрированная матрица» SNP-маркеров:

C = M - 2P, где

M – матрица содержания SNP-файла (закодированная зиготность SNP-маркеров, где 0 – гомозиготность AA, 1 – гетерозиготность AB, 2 – гомозиготность BB), P – матрица частот генотипов, где Pi=2(pi-0,5); pi – частота встречаемости минорного аллеля (MAF) для i-го SNP.

C’ – транспонированная матрица C.

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

Исходя из этого, решение уравнения GBLUP будет иметь вид:

 где

G-1 – обратная матрица геномного сходства.

Оценки племенной ценности, получаемые в результате применения GBLUP, называют GEBV (Genomic EBV).

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

где y – наблюдаемое значение признака;

µ – популяционная константа X – фиксированный эффект, учитывающий группу содержания, год и сезон;

w – живая масса индивида, рассматриваемая, как коварианта;

c – совокупный эффект SNP;

a – вектор случайных аддитивных генетических эффектов при σ_α^2 ~ N(0, Gσ_α^2), где G – матрица геномного сходства особей, вычисленная на основании данных генотипирования с использованием ДНК-чипов средней или высокой плотности, а σ_α^2 – дисперсия аддитивных полигенных эффектов;

K – коэффициент регрессии оцениваемого признака на живую массу индивидов;

e – остаток модели при σ_e^2 ~ N(0, Iσ_e^2), где I – единичная диагональная матрица, σe2 – остаточная дисперсия.

В данном исследовании GWA-анализ проводили с использованием ДНК-чипа Porcine GGP HD (платформа GeneSeek Genomic Profiler, «Neogene», США), содержащим ~70 тыс. SNP. Контроль качества и фильтрацию данных генотипирования для каждого SNP и каждого образца выполняли с использованием программного пакета PLINK 1.9, применяя следующие фильтры:

- call-rate по всем исследуемым SNP для индивидуального образца не ниже 90% (--mind 0,10);

- call-rate для каждого из исследованных SNP по всем генотипированным образцам не ниже 90% (--geno 0,10);

- частота встречаемости минорных аллелей (MAF) ≥0,05 (--maf 0,05);

- отклонение генотипов по SNP от распределения по Харди-Вайнбергу в совокупности протестированных образцов (--hwe 1e-3).

Достоверность результатов GWAS оценивается по P-значению, которое подвергается коррекции на множественное тестирование. Пороговое значение может быть разным в разных исследованиях, но чаще всего для анализа, в котором рассматривается десятки и сотни тысяч SNP, оно принимается равным 5×10-8. Для корректировки результирующего значения P используются методы проверки нескольких гипотез, такие как Бонферрони, Бенджамина-Хохберга или расчет ожидаемой доли ложных отклонений (FDR)

В результате проведенного анализа, для проведения анализа ассоциаций были отобраны 40 069 SNP, прошедшие фильтрацию по всем параметрам качества.

Для выявления ассоциаций SNP маркеров с изучаемыми признаками применяли регрессионный анализ, реализованный в PLINK 1.90 (--assoc --adjust --qt-means). Для подтверждения достоверного влияния SNP и определения значимых регионов в геноме свиней использовали тест для проверки нулевых гипотез по Бонферрони при пороговом значении p <1,25×10-5, 0,05/40069. Данные визуализировали в пакете qqman с помощью языка программирования R.

Поскольку стратификация популяции сильно влияет на надежность GWAS, анализ участка квантиль-квантильности (Q-Q) считается эффективным способом определения надежности результатов GWAS. В графике Q-Q горизонтальная ось представляет ожидаемое значение -log10P, а вертикальная ось представляет наблюдаемое значение -log10P. Диагональная линия представляет зависимость вида y=x, показывает 95% доверительный интервал, основанный на бета-распределении. Общее отклонение над диагональной линией идентичности обычно указывает на сильную стратификацию популяции. Отклонения от диагональной линии указывают на то, что предполагаемое распределение неверно или что образец содержит значения, возникающие каким-либо другим образом, аналогичные тем, которые порождаются истинным объединением. График Q-Q строится с использованием средств визуализации языка программирования R.

Для поиска генов-кандидатов, локализованных в области идентифицированных SNP, использовали геномный ресурс Sscrofa11.1. Функциональные аннотации генов выполняли с привлечением веб-программы PANTHER 19.0, базы данных DAVID и Pig QTL data base. Чтобы получить представление о функциональном обогащении генов, были впоследствии проведены анализы путей обогащения термина Gene Ontology (GO)

Впервые в России апробация GWAS анализа по хозяйственно-полезным признакам свиней породы дюрок произошла в 2019 г.

По результатам геномных оценок племенной ценности массив животных составил 408 голов, что составляет 47,2% от общего количества изучаемого поголовья. На рисунке 1 представлена диаграмма распределения GEBV по ADG. Наибольшее значение составило +124,5 г, наименьшее – +1,2 г.

Рисунок 1 - Диаграмма распределения лучших геномных оценок племенной ценности по показателю среднесуточного прироста хряков породы дюрок

DOI:10.60797/JAE.2024.50.2.1

По данным значениям GEBV, описанным выше, на 408 головах свиней породы дюрок, проведено полногеномное ассоциативное исследования, с учетом фильтрации данных, представленных в разделе «Материал и методы исследования», по результатам которых осталось 40 069 SNP и 400 особей. Визуализация Манхэттен плота представлена на рисунке 2. Общее количество достоверно значимых SNP составила 28, из них 3 выходящие за порог полногеномных значений (p<1,25×10-5) и 25 относящихся к диапазону достоверности (от p<1,25×10-4 до p=1,25×10-5).

Рисунок 2 - Распределение однонуклеотидных мутаций по хромосомам свиней породы дюрок в связи с уровнем достоверности (-log10 (p) по вероятностному суггестивному значению (синяя линия, p<1,25×10-4) и критерию Бонферрони (красная линия, p<1,25×10-5) для GEBV среднесуточного прироста

DOI:10.60797/JAE.2024.50.2.2

Рисунок 3 - Квартиль вероятностного распределения ожидаемого и наблюдаемого отклонений от нормального распределения для значений достоверности

DOI:10.60797/JAE.2024.50.2.3

Аннотация выявленных SNP по сборке Sscrofa11.1 идентифицировало 79 генов, расположенных на хромосомах 1, 3, 4, 5, 9 и 13 (всего 6 из 18 хромосом). Самая высокая достоверность, превышающая критерий Бонферрони у SNP WU_10.2_1_180284104 - P=6,15 х 10-7, наименьшая – у SNP DRGA0009891 – P=10,6 х 10-4. (табл.1).

Наибольшее количество генов расположены на 4 хромосоме (n=60).

Далее все выявленные гены были соотнесены с биологическими функциями в веб-программе PANTHER 19.0 (рис. 3). Идентифицировано 11 групп генов GO, включая метаболический процесс, который из всех представляет наибольший интерес с сопряженностью со среднесуточным приростом.

При этом, распределение по группам было следующее: GO:0009987 – 22%, GO:0051179 – 4,9%, GO:0044419 – 1,1%, GO:0065007 – 18,1%, GO:0050896 – 9,3%, GO:0042592 – 2,2%, GO:0032502 – 7,7%, GO:0048511 – 0,5%, GO:0032501 – 8,2%, GO:0008152 – 7,1%, GO:0002376 – 1,1% и не была присвоена категория – 17,6%.

Рисунок 4 - Распределение идентифицированных генов по биологическим функциям (GO): 

1 - категория не присвоена; 2 - биологический процесс, участвующий в межвидовом взаимодействии организмов (GO:0044419); 3 - биологическая регуляция (GO:0065007); 4 - клеточный процесс (GO:0009987); 5 - процесс развития (GO:0032502); 6 - гомеостатический процесс (GO:0042592); 7 - процесс иммунной системы (GO:0002376); 8 - локализация (GO:0051179); 9 - метаболический процесс (GO:0008152); 10 - многоклеточный организменный процесс (GO:0032501); 11 - ответ на стимул (GO:0050896); 12 - ритмический процесс (GO:0048511)

Примечание: выделенная часть в диаграмме – метаболический процесс (GO:0008152)

DOI:10.60797/JAE.2024.50.2.5

Цель GO – отразить современное состояние знаний в области биологии, поэтому он постоянно пересматривается и расширяется по мере накопления биологических знаний. По итогу, аннотация генов в с привлечением баз данных DAVID и Pig QTL идентифицировала 30 и 6 генов, соответственно, по которым есть биологический путь и которые зарегистрированы ранее по изучаемому виду животных (табл.2). Интерес представляют гены, внутри которых локализованы SNP, а именно: NEDD4L, KIDINS220, NAXE, S100A1, PLA2G4A и GRM8. На них в дальнейшем можно разработать и апробировать тест-системы, которые смогут детектировать хряков с высоким показателем среднесуточного прироста.

Анализ обогащения показал наибольшую локализацию генов по кластеру 1 (коэффициент обогащения = 6,88), принадлежащего к функциям: связывание ионов кальция – 12 генов, Р=4,1х10-5; кальций-зависимое связывание с белками – 10 генов, Р=2,1х10-12 и перинуклеарная область цитоплазмы – 8 генов, Р=9,4х10-4.

В группу связывания ионов кальция (GO:0005509) входят гены: S100A2, S100A16, S100A1, S100A6, S100A13, PLA2G4A, S100A5, S100A4, S100A14, S100A3, S100A7, THBS3. В кальций-зависимое связывание с белками (GO:0048306) – гены S100A2, S100A16, S100A1, S100A6, S100A13, S100A5, S100A4, S100A14, S100A3, S100A7. И в группу, отвечающую за перинуклеарную область цитоплазмы (GO:0048471) – гены SNAPIN, S100A16, S100A1, S100A6, S100A13, S100A4, S100A14, THBS3.

Данное научно-практическое исследование можно использовать как методические аспекты структуризации расчета геномной оценки в области свиноводства, применяя актуальные программы для получения достоверного результата.

По выявленным генам, внутри которых локализованы SNP, а именно, NEDD4L, KIDINS220, NAXE, S100A1, PLA2G4A и GRM8, в дальнейшем будут разработаны тест-системы и проведен массовый скрининг на остальных коммерческих породах свиней.

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