Работа над ошибками: неудача в создании дизайна нового калькулятора MDLP World K23

К сожалению, после полгода работы над новой версией генетико-этно-популяционного калькулятора MDLP World K23, должен признать, что из-за неправильных методологических решений, от релиза калькулятора придется отказаться.

В чем же кроются oшибки?

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

Полученные путем слияния данные хранятся в стандартном формате Plink. Но, как говорится, не все так просто.  Хотя в некоторых из вышеприведенных источниках данные поставляются уже в готовом формате Plink, большинство все же имеют другой формат. Ситуацию усложняет и тот факт, что данные (генотипы) в  разных файлах могут даваться в разных ДНК-цепях; геномные координаты соотвествуют разным  геномным билдам человеческого геномного референса (b36 или b37). Наконец, вместо классических генотипов, могут приводится данные по частоте минорных аллель.
Но самая главная проблема — это существенное расхождения числа генотипированных снипов в разных базах данных. В самых больших проектах, вроде 1000genomes и HapMap, число прошедших QC-фильтры снипов приближается к миллионам, в то время как в научных данных из публикаций, число снипов в панели генотипирования выборки часто не превышает 100-300 тысяч. Вторая важная проблема — использование разных платформ генотипирования, например Illumina или Affymetrix; и по причине этого только часть снипов выборки, типированной с помощью машин Illumina перекрывается со снипами выборки, типированной на другой платформе. Поскольку большинство полногеномных этно-популяционных анализов требует наличие однородной выборки, то это означает что при собирании референсной популяций неизбежно приходится придерживаться четкого правила —  удаление снипов имеющих низкую частоту (менее 0.01) генотипирования в одной из соединяемых выборок. Совершенно очевидно, что при поэтапном слиянии выборок, большинство снипов теряется; так, вместо первоначальных 2 миллионов снипов, на конечном этапе мы можем получить всего 20 000 снипов.

Импутирование отсутствующих генотипов —  потенциальное решение проблемы?

Одним из возможных вариантов решения проблемы сохранения баланса между количеством и качеством снипов является использование алгоритмов импутирования (определения) отсутствующих генотипов. Обычно для этих целей я пользуюсь программой Beagle, однако в этот раз я решил использовать более мощную комбинацию программ ShapeIT и IMPUTE и  метод импутирования снипов за счет использования большой референсной панели аутосомных гаплотипов (из 1000 genomes),  гарантирующей более аккуратное определение снипов.

Я не буду загружать читателя избытком технической информации касательно аспектов фазирования. Скажу лишь, что после нескольких месяцев использования эвристического метода проб и ошибок, я пришел к оптимальному сочетанию снипов в файле полученном путем соединения вышеуказанных выборок.  В целях сохранения совместимости с данными генотипирования от коммерческих компаний, я снизил  до примерно 130 000 (именно такое число снипов является общим при слиянии данных 1000genomes, файлов 23andme и Family Finder от FTDNA). Именно эти снипы и были задействованы при параллельных процессах импутирования в IMPUTE.

В результате я рассчитал аллельные частоты примерно 130000  снипов по референсным популяциям, взятым из академических источников. С целью сохранения совместимости с коммерческими данными, были выбраны только те снипы, которые присутствуют либо в последнем чипсете 23andme, и/либо в последнем чипсете FAmilyFinder от FTDNA. Я высчитал средние значения каждого из 23 компонентов для всех референсных популяций. Также имеются данные по кластерному определению каждого из снипов. Судя по разбивке популяций, они выглядят весьма правдоподобными. Но некоторые кластеры вызвают вопросу. Поэтому перед тем как делать релиз очередного этно-популяционного ДНК-калькулятора, я решил воспользоваться услугами бета-тестеров. В первую очередь, я надеялся на помощь Александр Бурнашев и Srkz (Сергея Козлова) поскольку они лучше остальных в русскоязычном секторе любителей ДНК-генеалогии понимают принципы работы калькуляторов на основе Додекад DIY.

Затем я произвел запуск тестовой версии ДНК-калькулятора с целью проверки надежности предковых компонент. Названия пока носят условный характер.
West-African 0,44

Caucausian 28,87

Near-East-Mediterranean 0,07

Samoyedic 1,66

Polinesian 1,27

North-European 11,36

Amerindian 1,12

East-European 36,89

South-African 0,34

Papuan-Melanesian 0,85

East-African 0,12

Central-Asian-Caucasian 2,37

Arabic 3,51

Nigerian 0

Central-African 0,04

Indian 9,98

Sino-Tibetan 0,13

Indo-Chineese 0

Sub-Saharian 0

West-Mediterranean 0

Для сравнения -усредненные результаты референсных беларусов

East-European 41,6645%

Caucausian 26,3905%

Indian 12,1348%

North-European 11,0739%

Near-East-Mediterranean 2,0315%

Indo-Chineese 1,4123%

Austronesian 1,1291%

Samoyedic 1,1257%

West-African 1,0845%

Polinesian 0,6104%

Nigerian 0,5530%

Arabic 0,4754%

South-African 0,3044%

North-African 0,0010%

North-European-Caucausian 0,0010%

Amerindian 0,0010%

Papuan-Melanesian 0,0010%

East-African 0,0010%

Central-Asian-Caucasian 0,0010%

Central-African 0,0010%

Sino-Tibetan 0,0010%

Sub-Saharian 0,0010%

West-Mediterranean 0,0010%

Критика и очевидные ошибки

Бета-тестеры указали на очевидные несуразности в результатах калькулятора:

У меня есть RAW Data своих родственников и еще пары человек, я прогнал их через калькулятор, добавил результаты ваши и вашей мамы и попробовал сравнить с таблицей. Пока на первый взгляд европейские популяции делятся на две группы — «выглядящие странно» и «выглядящие нормально». Russian_Center смотрится странно, Center-Russian — нормально. Russian_North странно, North Russian нормально. В чем между ними разница, какая популяция откуда?

Странность в том, что прослеживаются две отдельные системы географического распределения компонентов в Европе, при этом с данными реальных протестированых (7 человек, из них неродственных образцов 4) пока что коррелирует только одна из них. 4 образца достаточно, чтобы строить предположения, но для уверенных выводов лучше бы несколько десятков.

Вспоминается «эффект калькулятора» имени Поляко, когда в популяциях, которые использовались для выделения компонентов, их баланс оказывается смещен. При этом сами компоненты реальны и полезны. На примере World-22, компонент North-European Mesolithic у финнов эталонных популяций оказался Finnish 23,6; Finnish-South 30,1; Finnish-North 37,5 . У финнов же с форумов и моих совпаденцев на Gedmatch компонент стабильно между 15 (ближе к западу) и 20 (ближе к востоку), у большинства 16-17. То есть он действительно хорош для Оракула, но только с этой поправкой. Далее, когда я смотрел результаты Оракула World-22 для восточных славян, их почти всех стабильно смещало в сторону запада. Компонент Athlantic_Mediterranean у них оказывался выше на 2-3 пункта, чем в русских эталонах, а North-East European слегка ниже.

Поэтому по World-23 я сразу захотел сравнить распределение с результатами людей, не участвовавших в расчете. У нас есть три реперные точки — Белоруссия в вашем лице, северо-восточные русские в лице моей семьи, есть образец центрального русского и дополнительно лицо смешанного славяно-ашкеназского происхождения. У всех распределение компонентов укладывается в одну систему, и близко к распределению во многих европейских популяциях в таблице. Основные компоненты East-European, Caucausian, Indian, North-European, Samoyedic.

А вот компоненты, по которым географическое распределение удивляет: Пара West-Mediterranean/North-European (North-European, когда его больше 15%). Есть ряд популяций, которые состоят только из этих двух компонентов, в треугольнике Сардиния-Оркнейские острова-Россия. Сформирован клин Сардиния-Прибалтика. У русских показано Russian North-Euro 72,59%, West-Med 27,05%. При этом у всех семерых протестированых West-Mediterranean ноль либо практически ноль, North-European в диапазоне 7-11. French_Basque получили всего 8,92% West-Mediterranean, вместо этого компонента Средиземноморье у них представлено 38,6% Caucausian

North-European-Caucausian встречается в нескольких популяциях, зато по всей планете, очень непонятный компонент. У протестированых в пределах 0-2 процента.

Пара мыслей по другим компонентам: По Indian очень интересное распределение, на первый взгляд это давно известный Gedrosia/Baloch. Однако если Gedrosia распространяется в Европу больше по маршруту Турция-Балканы-Германия-Британия, то Indian World-23 преобладает у финно-угорских народов (маршрут Каспий-Поволжье-Финляндия), а в направлении Западной Европы спадает. А ведь давно замечено в других калькуляторах неожиданное проявление небольших количеств Индии у финно-угров. Пока что результаты протестированых распределение Indian не опровергают и не подтверждают.

Samoyedic пик у юкагиров, скорее его можно назвать Siberian или North-Eurasian

Результаты (названия компонентов вписал в калькулятор сам, может слегка не совпадать) Русская из Сибири (близка к северным) 0.00% North-African 0.88% North-European-Cauca 0.10% Austronesian 1.74% West-African 22.08% Caucasian 1.09% Near_East-Med 8.18% Samoedic 2.02% Polinesian 7.96% North-European 0.99% Amerindian 35.79% East-European 0.91% South-African 0.76% Papuan-Melanesian 0.00% East-African 4.32% Central-Asian-Caucas 2.45% Arabic 0.00% Nigerian 0.00% African-1 10.72% Indian 0.00% Sino-Tibetan 0.00% Indo-Chinese 0.00% African-2 0.01% West-Mediterranean

Русский из центральных областей: 0.00% North-African 1.98% North-European-Cauca 0.38% Austronesian 0.75% West-African 25.45% Caucasian 1.57% Near_East-Med 2.48% Samoedic 1.12% Polinesian 9.70% North-European 1.03% Amerindian 35.59% East-European 0.01% South-African 2.32% Papuan-Melanesian 0.00% East-African 2.81% Central-Asian-Caucas 3.44% Arabic 0.00% Nigerian 0.00% African-1 11.37% Indian 0.00% Sino-Tibetan 0.00% Indo-Chinese 0.00% African-2 0.00% West-Mediterranean.

я прогнал через ваш калькулятор около 20 человек из разных мест европейской части России/СНГ. Ну явно что-то не так с частью исходных популяций. ИМХО, до тех пор, пока непонятны причины этого, результаты так и будут странными. Ощущение, что просто вместе объединены несколько разных калькуляторов с разными популяциями. Почему, к примеру, Russian_Center, Russian_cossack, Russian_North на выходе имеют 46-50% таджикско-лакского компонента Central-Asia? Или это так и задумано? В чем тогда цель? Коротко пробегусь по дереву компонентов, что получилось на данных людей, не участвовавших в расчете: Sub-Saharian (+ кое-где East-African) — выскочило практически у всех в диапазоне 1-3%. Уральско-сибирская веточка смотрится адекватно, за исключением компонента Sino-Tibetan, который проявился в диапазоне 0-2,2 без какой-либо видимой закономерности. Дальше у нас самая странная ветка, где объединились Центральная Азия, Европа и Африка. Как раз к ней-то и относятся странные популяции. У реальных же людей из этой ветки наиболее нормально проявилась Европа — Балтика+Северо-Западная Европа растут при приближении к Литве, Central-Asia уже более разбросана, но общее распределение правдоподобно, однако у всех присутствует и Африка, причем у большинства сумма трех компонентов North African, West African, Central African в районе 4-5 процентов. Indian проявляется в очень небольших количествах, в основном у тех же людей, что и в других калькуляторах. Выглядит адекватно. Ближневосточная ветка явных несообразностей не содержит, но распределение не очень четкое, ощущение зашумленности. Jewish у женщины 1/4 ашкеназкого происхождения не проявился (3,13), вместо него 8,41 Near-East. Почему-то на втором месте по Near-East вы, видимо, шум. Европейская ветка — North-Eastern-Europe и West-Southern-Europe выглядит отменно. Правда, вызывает сомнение, что в этом калькуляторе эталоном северо-восточных европейцев стали финны, а литовцы, наоборот, оказались аутсайдерами со своим отдельным компонентом. Как-то привычнее наоборот )) . Но так даже интереснее. Таблицу исходных популяций почти не прорабатывал по вышеизложенным причинам.

Надеюсь, вы не в обиде на мою критику )) ИМХО, единственное, чем я могу помочь — это именно посмотреть взглядом человека, который много возится с Оракулами. В чем все-таки может быть причина задвоения кластеров? Я не понимаю специфику работы с Admixture, поэтому сам предположить не могу. Разве что другой набор снипов в лишних популяциях? Или расчет кластеров идет несколькими отдельными кусками?

 Как я и предполагал при первоначальных тестовых запусках выявились странности с африканским и центрально-азиатским компонентом.  Впрочем, то что будут странности было видно уже на дереве Fst между компонентами. Я считаю, что корень проблемы  -это фазированные и импутированные данные. Думаю причина именно в этом. Если это так то придется убрать практически все странные или неправдоподобные популяции из коллекции. Если и после изменений с новыми данными ничего не получится, придется убирать провальные популяции.

Насчет африканских компонентов есть еще одно объяснение. Когда я брал выборку с неандертальцами и денисовцами, они практически всегда оказывались в пигмейском кластере (он видимо самый архаичный и древний). Так что скорее в этих 4-5% содежится и реальный африканский след, и «прокси» сигнала древнего адмикса с неандерами и денисовцами. Как показывают неандертальские тесты на 23andme и NG, диапазон такого адмикса с евразийцами составялет 2-4 % В целях поиска первоисточника  ошибки,  я переделал еще раз данные под калькулятор. К сожалению, в первых исправлениях калькулятора пришлось отказаться от популяций Panasian dataset.

Как я не пытался, популяции в этой выборке упорно не кластеризуются с другими юго-восточными азиатскими популяциями )).

В итоге, они практически всегда образовывали в результатах Admixture свой собственный кластер.

Maya Auca Quechua Karitiana Pima Ami Atayal Melanesians Zhuang Han_Cantonese Hmong Jiamao Jinuo Han_Shanghai Uyghur Wa Alorese Dayak Javanese Batak_Karo Lamaholot Lembata Malay Mentawai Manggarai Kambera Sunda Batak_Toba Toraja Andhra_Pradesh Karnataka Bengali-Assamese Rajasthan Uttaranchal Uttar Pradesh Haryana Spiti Bhili Marathi Japanese Ryukyuan Korean Bidayuh Jehai Kelantan Kensiu Temuan Ayta Agta Ati Iraya Minanubu Mamanwa Filipino Singapore_Chinese Singapore_Indian Singapore_Malay Hmong (Miao) Karen Lawa Mlabri Mon Paluang Plang Tai_Khuen Tai_Lue H’tin Tai_Yuan Tai_Yong Yao Hakka Minnan

Итак в переделанной версии нового калькулятора, этно-генографическое происхождение моего генома выглядит следущим образом (в процентах): North-East-European 37,03 Caucasian 25,36 Balto-Scandinavian 10,73 East-Mediterranean 8,67 Near-Eastern 5,46 South-Indian 4,07 Amerindian 2,25 Central-Asian 2,09 Siberian 1,62 East-African 0,65 Dagestan 0,63 Austronesian 0,54 Tibeto-Burmese 0,41 Pygmy 0,39 North-Indian 0,07 West-African 0,03 Paleo-Siberian 0,01 North-African 0 South-Asian 0 Central-African 0 Equatorial-African 0 Sub-Saharan 0 Atlantic-Mediterranean 0 Иными словами, в своей генетической основе я на 33% северо-восточно-европейцев (наследие кордедов-шнуриков и уральцев ямочно-гребенчатой культуры), на 25,36 — носитель анатолийско-кавказских аллелей (в Европе львиная доля этих аллелей — наследие потомков неолитических поселенцев-земледельцев) . East-Mediterranean 8,67 + Near-Eastern 5,46 = 14% ближневосточных аллелей (семито-арабско-еврейская кровь). 11% балто-скандинавских аллелей — наследие мезолитических охотников-собирателей северной Европы. 4% дравидских генов — это явный цыганский след. Ну и наконец, Amerindian 2,25 + Central-Asian 2,09 + Siberian 1,6% — примерно 5% тюркско-сибирско-монгольских генов.

В ходе последующих изменений и уменьшении числа снипов всего до 900 (!)  получилось впихнуть и паназиатскую выборку, богатую популяциями восточной и юго-восточной Азии, и данные из статьи Xing et al.2010. Для коррекции данных, добавил еще в выборку генотипы неандертальцев и денисовского человека. На этот раз уменьшил и число априорных компонентов до 15. На выходе я получил матрицу дивергенции (Fst) между прогнозируемыми предковыми компонентами. Поскольку это есть ни что иное, как матрица дистанции, то я преобразовал ее в дерево с помощью метода Wheeler (2009) NINJA large scale Neighbour-Joining clustering method.

1069319_10201480061153220_287033667_n

В результате получилось дерево, которое весьма неплохо коррелирует с топологиями Y-хромосомного и митохондриального древ человечества. За значение 1.0 принят возраст корня, все остальные цифры в дробях на ребрах дерева (ветвях) означают возраста относительно корня. Например, цифра 0.27 на узле обозначающем диверегенцию австронезийского и восточноазиатского компонента, означает что компоненты дивергировали относительно друг друга 0.27*(возраст дерева) лет тому назад. Теоретически, если мы сможем откалибрировать возраст дерева (или принять априорную оценку возраста корня), можно просчитать возраст диверегенции. Было бы любопытно услышать современные оценки возраста расхождения расовых стволов и возраста субрас в современной антрополгии, так как аутосомные маркеры достоверно лучше коррелируют с антропологическими или расовыми признаками, чем однородительские маркеры.  Предположим, что  возраст корня составит 100000 лет. огда время выделение западных африканцев -80 000 лет
центральных афрриканцев — 80 000 лет, выход из Африки 55 000 лет; разделение на западный и восточный ствол — 50 000 лет; выделение ближневосточно-афроазиатской ветви 45 000 лет; европейской ветви — 48 000 лет. Oколо этого же времени выделился восточный ствол, от которого практически сразу 50 000 лет тому назад ответвились дравиды и чуть позже австралоиды, и так далее.

Комментарий  Евгения Марчукова:

«А. А. Зубов (2004) писал, что начальный этап дифференциации Homo sapiens не был единомоментным, он охватывал период от раннего верхнего палеолита до неолита (включительно), а, если учесть вклад периферийных форм – процесс уходит корнями в средний палеолит. Важно, что вклад «восточного центра» в формирование современного населения отвергать нельзя, хотя линия его собственной эволюции угасла в конце среднего плейстоцена. Она была восстановлена за счет мигрантов с запада представителей «позднего архаичного сапиенса». Позже новая миграция на восток сформировавшегося сапиенса (Мальта и алтайские находки) привели к метисации с местными группами и приняла участие в формировании монголоидного ствола (Чжоукоудянь сходен с Сунгирем). «Сетевидная» эволюция ответственна за появление разных «евро-монголоидных» и «монголо-австралоидных» вариантов. Маркина Гора, Гримальди (на западе) и Лючан (на востоке) указывают на наличие периода, когда Homo был разделен на два надрасовых конгломерата – «западный» и «восточный» без четкой границы.

Доказательством этого является следующее: реликтовые одонтоглифические признаки объединяют европейцев и негроидов, австралийцев с монголоидами. По данным генетики первое деление внутри вида имеет древность 92 тыс.лет, деление внутри больших стволов – 39 тыс.лет. Бушмены и австралийцы выделились первыми- 60-90 тыс. лет, а сохранились благодаря изоляции. В направление к нашим дням происходит изживание архаики в морфологии людей. Чаще всего, идентификация субтаксонов возможно в мезолите — (неолите) – верхнем палеолите. В первичном «восточном стволе» сформировались монголоидный и австралоидный компоненты, а «западный» ствол дал европеоидный и негроидный компоненты сапиенса. Процесс дифференциа-ции имеет непрерывный характер.

Дивергенция расовых стволов происходила, по данным генетики, в пределах последних 100 тыс.лет. Ранее всех выделилась ветвь монголоидов (92 тыс.лет), протонегроиды и протоевропейцы разделились 33 тыс.лет назад. В Африке формирование экваториальных комплексов происходило еще раньше. Монголоидный и негроидный комплексы оформились окончательно позже. Европеодная и монгоидная расы наиболее метисированы.»

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

26,09% West-Southern-Europe
24,27% North-Eastern-Europe
15,62% Caucasian-Asian
13,28% Baltic
8,14% Near-Eastern
4,44% Uralic
2,70% Central-Asia
2,32% Amerindian
1,24% Sino-Tibetan
1,16% Central-African
0,41% Sub-Saharian
0,26% North-African
0,03% Indian
0,03% West-African
0,02% Altaic-Siberian
0,01% South-Asian
0,00% East-African
0,00% West-Northern-Europe
0,00% Jewish

1000875_10201544272598466_1144492335_n

 

За кулисами: как создавался этно-популяционный калькулятор World-22

Летом 2011 года я создал целый рядсобственных модификаций получившего широкую известность калькулятора DIY Dodecad гениального грека Диенека Понтикоса. К моему приятному удивлению, за прошедшее время калькулятором успело воспользоваться несколько тысяч людей, некоторые из которых даже выложили свои результаты в Интернете.  Разумеется, многие также разместили и свои собственные интерпретации полученных результатов. Некоторые из приведенных в комментариях интерпретации выделялись (в хорошем смысле этого слова) высоким академическим уровнем, но мне попадались и такие комментарии, при чтении которых становилось понятно, что авторы не только не понимают принципов и сути парадигмы анализа, предложенного Понтикосом, но и — что гораздо хуже — выдавали свои фантазии за действительности. Особенно часто мне попадались подобные фантастические рассуждения в русскоязычном секторе Интернета.Пример такого невежества можно найти в рассуждениях само-провозглашенного академика ДНК-генеалогии Анатолия  Клесова:

Но и в этом случае различия все равно будут между русскими и монголами. Качественно и как-то полуколичественно его можно рассматривать, но не в виде профанации, как это делает Понтикос. Более того, это рассмотрение – если правильно – надо проводить не на выбранных маленьких фрагментах, а действительно по всему геному. На маленьких фрагментах будут вылезать отдельные особенности – то присущие в основном, например, гаплогруппам Y-I2 и мтДНК-Н, то кому-то еще. И это еще будет зависеть от разрешения, которые и обозначают индексами К=4, К=8 и другими. То есть берут маленький фрагмент генома, да еще с малым (или бóльшим) разрешением, стягивают в точку, и все равно получают в целом ерунду. Но для коммерции годится. Годятся для коммерции и вот такие, в частности, «открытия» того же Понтикоса: Перевод: Интересно то, что европейская популяция показывает присутствие американских индейцев, что показывает и f-статистика, и она же показывает присутствие компонента с Сардинией. Как видим, Понтикос уже забыл, что названия им придуманы как попало, и уже придает им абсолютные значения. Про Сардинию Понтикос уже вошел в состояние экзальтации. Он придает Сардинии некую пра-европейскую значимость, на основании, конечно, этой ерунды с «геномом», который анализирует как хочет. Пример – он трубил по всему свету, что Отци, «ледовый человек», имел геном «Сардинии». Однако только что опубликована статья о том, что Отци – никакая не Сардиния, а типичная Центральная Европа. Ну, и что делать будем? Понтикос, с его страстным желанием сенсаций, каждый раз наступает на одни и те же грабли. Впрочем, фарс продолжается. Теперь тем же занялся некто российский Веренич, а именно тоже насчитывает «польскую компоненту», пользуясь подходом своего гуру-Понтикоса.

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

Предварительные замечания

Как вы возможно знаете, MDLP блог не обновлялся с февраля 2012 года.  Полгода тому назад я пообещал себе, что я не буду писать новые сообщения на MDLP блоге до те пор пока я не напишу краткую научный отчет о проделенной работе. Так как приоритеты завершения научной работы были важнее рутиного обновления блога,  то  в связи с нехваткой времени, я был не в состоянии продолжать обновление блога на регулярной основе, в связи с нехваткой времени, я должен был внести изменения в свой исследовательский график. Поэтому я решил воздерживался от размещения новых данных на блоге в течение нескольких месяцев, фокусируясь на более важных вопросах. Несмотря на все ограничения, я продолжал втайне работать  на проектом MDLP, сбором необходимых данных и выполением различных ‘геномных’ экспериментов в целях достижения своей конечной цели. Однако с течением времени, некоторые результаты секретных экспериментов с новыми полногеномными популяционными выборками и инструментами в конечном итоге просочились в Интернет,  порождая огромный интерес к моему проекту. После выпуска новой версии моей собственной модификации DIYDodecad калькулятор на сайте Gedmatch.com, я был буквально завален письмами пользователями сервиса Gedmatch.com.
Тогда я осознал свою основную стратегическую ошибку, которая заключалась в  отсутствии подробной документации к выпущенными мной данными и результатам анализа, и почувствовал себя обязанным разместить более подробные разъяснения. Очевидно, я начну новую серию публикацию в своем блоге,  которая будет тесным образом связанна с теми аспектами моей работы, которая наиболее интересует общественность, то есть с калькулятором MDLP World22.

Основы отбора референсных популяций калькулятора MDLP World22.

Референсный набор  популяций в этом калькуляторе был собран в программе PLINK   методом «intersection&thinning» ( дословно «пересечением и истончением») образцов из различных источников данных: HapMap 3 (отфильтрованный набор данных КЕС, YRI, JPT, CHB), 1000genomes,   Rasmussen et al. (2010),   HGDP (кураторская база данных Стэнфордского университета), Metspalu et al. (2011),  Yunusbayev et al (2011), Chaubey et al. (2010) и т.д.
Кроме того, я отобрал произвольным образом по 10 сэмплов (или максимальное количество доступных сэмплов в тех случаях, когда общее число сэмплов в популяции было меньше 10) от каждой европейской страны, представленной в панеле базе данных POPRES. Наконец, для того чтобы оценить степень корреляции между современным и древним генетическим разнообразием населения Европы, я также включил в выборку образцы древней ДНК Эци (Keller et al. (2012)) ,  образцы житлей шведского неолита Gök4, Ajv52, Ajv70, Ire8, STE7 ( Skoglund et al. (2012)) и 2 образца La Braña  — останков мезолитических жителей Пиренейского полуострова (Sánchez-Quinto et al.(2012)).
Затем я добавил 90 образцов — анонимизированных данных — участников моего проекта. После слияния вышеупомянутых наборов данных и истончения набора SNP с  помощью особой команды PLINK, я исключил SNP-ы с  более чем 1% минорных аллелей. После чего я отфильтровал дубликаты, лиц с высоким уровнем общих по происхождению идентичных сегментов (IBD). В качестве критерия фильтрации  были использованы расчеты IBD в Plink, где IBD представлена как средняя доля аллелей общих между двумя людьми по всем анализируемым локусам.  Затем я удалил из выборки лиц с  высоким коэффициентом предпологаемого родства (коэффициенты родства были вычислены в программном обеспечении King). Для получения более стабильных результатов, я также отфильтровал сэмплы с более чем 3 стандартными отклонениями от средних данных  по популяции. Поскольку коэффициент родства может быть надежно определен с помощью оценки HWE (ожидания, вытекающего из закона Харди-Вайнберга) между SNP-ами с той же базовой частотой аллелей, то SNP-ы с существенным отклонением (p < 5.5 x10−8) от  ожидания Харди-Вайнберга были удалены из объединенного набора данных. После этого я выделил те SNP-ы, которые присутствовали в чипах Illumina / Affymetrix, и затем произвел фильтрацию снипов на основе расчетов степени неравновесного сцепления  (в этой я использовал хромосомное ‘окно’ размером в 50 базовых пар, с шагом 5 базовых пар и пороговым значением уровня сцепления R ^ 2, равным 0,3).
По окончанию этой сложной последовательности операций, я получил окончательноый набора данных, который включал в себя 80 751 снипов,  2516 человек и 225  референсных популяций.

Анализ этно-популяционного адмикс

 В ходе следующенго этапа, окончательный набор данных по референсным популяциям (которые я храню в linkage-формате PLINK) был обработан в программе Admixture.  Во время выбора подходящей модели проведения теста на этно-популяционный адмикс, я столкнулся с крайне трудной задачей: как было показано в профильных научных исследованиях (Patterson et al.2006) количество маркеров, необходимых для надежной стратификации популяций в анализе обратно пропорциональна генетическому расстоянию (фСТ) между популяциями. Согласно рекомендациям пользователей программы Admixture, считается что примерно 10 000 генетических SNP-маркеров достаточно для выполнения интер-континентальной GWAS -коррекции обособленных популяций (например, уровень дивергенции между африканскими, азиатскими и европейскими популяциями  FST > 0.05), в то время как для аналогичной коррекции между внутриконтинентальными популяциями требуется более чем 100000 маркеров (в Европе, например, ФСТ < 0.01). Для повышения точности результатов Admixture я решил использовать метод, предложенный Dienekes. Этот метод позволяет преобразовать частот аллелей в «синтетические» индивиды (см. также пример Зака Аджмала из проекта HarappaDNA). Идея метода довольно проста: сначала необходимо запустить unsupervised анализ Admixture с целью вычисления частот аллелей в так называемых предковых компонентов, а затем на основании аллельных частот сгенерировать «фиктивные популяции».  Именно эти фиктивные популяции и индивиды будут использоваться в ходе чистых референсов в ходе последующего анализа этно-популяционного анализа.
Впрочем, как и любой другой исследователь, работающий над четким решением проблемы этно-популяционного адмикса, я вынужден считаться с ограничениями этого подхода. Хотя я и отдаю себе отчет в существовании явных методологических подвохов в использовании смоделированных искусственных индивидов для определения адмикса в реальной популяции, я все же скорее склонен согласиться с Понтикосом, которые считаeт полученных в ходе  аллельно-частотного моделирования «фиктивных индивидов» лучшей аппроксимацией древних генетических компонентов мирового народонаселения.Как бы то не было, моделирующий подход, предложенный Диенеком и Заком, сослужил свою хорошую службу, поскольку были мной были получены  значимые результаты в ходе создания нового калькулятора. Сначала я произвел unsupervised Admixture (при значении К = 22, т.е 22 кластера частот аллель или предковых компонентов). По выполнению анализа нами были получены оценки коэффициентов адмикса в каждой из этих 22 аллельных кластеров, а также частоты аллелей для всех SNP-ов в каждой из 22 родовых популяций.
Затем я использовал мнемонические обозначения для каждого компонента (имена для каждого из компонентов выведены в  порядке их появления). Нужно помнить, что обозначения этих компонентов носят скорее мнемонический условный характер:
Pygmy
West-Asian
North-European-Mesolithic
Tibetan
Mesomerican
Arctic-Amerind
South-America_Amerind
Indian
North-Siberean
Atlantic_Mediterranean_Neolithic
Samoedic
Proto-Indo-Iranian
East-Siberean
North-East-European
South-African
North-Amerind
Sub-Saharian
East-South-Asian
Near_East
Melanesian
Paleo-Siberean
Austronesian
Вышеупомянутые частоты аллель, вычисленные в ходе unsupervised (безнадзорного) анализа (Admixture K = 22) объединенного набора данных, были затем использованы для симуляции синтетических индивидов, по 10 индивидов на каждую из 22 предковых компонент.  Это симуляционное моделирование проводилось с помощью PLINK команды -simulate Когда моделирование было закончено, я сделал визуализацию расстояния между симулированными индивидами с использованием многомерного масштабирования.
На следущем этапе, я включил группу смоделированных индивидов (220 индивидов) в новую эталонную популяцию. После чего я запустил новый анализ А, на этот раз в полном «поднадзорном» режиме для K = 22, причем полученные в ходе симуляционного моделирования фиктивные популяции фиктивных индивидов использовались в качестве новых референсных эталонных групп.  На конвергенцию 22 априорно заданых предковых компонентов было затрачено  31 итераций (3 7773,1 сек) с окончательным loglikelihood: -188032005,430318 (ниже приведена таблица значений Fst  между расчетными ‘предковыми’ популяциями):
Приведенная выше матрица  Fst дистанций  была использована для определения наиболее вероятной топологии NJ-дерева всех 22 предковых компонентов ( примечание: в качестве outgroup-таксона использовался South-African component). Индивидуальные результаты ‘поднадзорного’ анализа этно-популяционных миксов (в формате Excel) для участников проекта были загружены на GoogleDrive.

MDLP World22 DIYcalculator

Выходные файлы «поднадзорного» анализа  Admixture K=22 (средние значения коэффициентов адмикса в референсных популяциях и значения Fst) были использованы для разработки новой версии DIYcalculator MDLP, который более известен под кодовым названием «World22» (онлайн версия доступна разделе Admixture-утилит на сервисе Gedmatch в рамках проекта MDLP). Как я уже упоминал выше, MDLP DIYcalculator работает на коде Dodecad DIY calculator (c) Dienekes Pontikos.
В свою очередь,  реализованная на сервисе  Gedmatch модификация DIYcalculator ‘World22’ комбинирована с  Oracle ‘World22’ MDLP, который также работает на коде Диенека и Зака Аджмала ​​(Хараппа/DodecadOracle). Программа «Oracle» работает в двух режимах. В режиме single population программа определяет ближайщие (к анализируемому геному) референсные популяции калькулятора Word22. В смешанном режиме, Oracle рассматривает все пары населения, и для каждой из пар вычисляет минимальное Fst-взвешенное расстояние между парой и анализируемым геномом, а также  коэффициенты сходства.
Предковые популяции (т.е. полученные в ходе симуляционное моделирования популяции — см. выше) обозначены в результатах Oracle суффиксом anc, в то время реальные современные и древние популяции обозначены суффиксом der.
Если у Вас возникли проблемы с пониманием/интерпретацией результатов Oracle и DIYcalculcator,  то я настоятельно рекомендую обратится к соответствующим темам в блогах  Dodecad и НаrappaWorld . Я полагаю, что не имеет особого практического смысла заново изобретать велосипед и слово в слово повторять то, что уже было написано более компетентными в этом вопросе людьми.

Что представляют собой компоненты MDLP World-22?

Один из наиболее частых вопросов, которые задают мне пользователи калькулятора, напрямую касается практической интерпретации референсных популяций и предковых компонентов в моих калькуляторах K = 12 и World-22 анализов в виду. Чуть выше по тексту я уже привел часть ответа на этот вопрос , но — как гласит старинная китайская пословица — одна картинка стоит десять тысяч слов. Вот почему я решил визуализировать компоненты на поверхности земного шара путем отображения коэффициентов адмикса. Избегая излишних премудростей, я воспользовался готовым рецептом Франсуа Оливье, который предложал  использовать графическую библиотеку статистического программного обеспечения R для отображения пространственной интерполяции  коэффициентов адмикса (Q матрица) в двух измерениях (где пространственные координаты записываются как географические долгота и широта).  Благодаря этому решению, мне удалось создать по 2 контурные карты на каждый из предковых компонентов.Pygmy (модальный компонент в  популяциях африканских пигмеев Biaka и Mbuti)

West-Asian (бимодальный компонет с пиком на Кавказе и юго-восточной части Ирана, приблизительно идентичен компонентам Caucasian/Gedrosia Диенека Понтикоса)
North-European-Mesolithic (локальный архаичный компонент с пиком в популяции древних европейских жителей Иберийского полуострова La_Brana и современной популяции саамов).
Tibetan (Indo-Burmese) component (Гималаи-Тибет)
Mesomerican (главный генетический компонент  у мезоамериканских америндов)

 

North-Amerind (нативный компонент северо-американских америндов)

South-Amerind (нативный компонент южно-американских индейцев)
  Atlantic-Mediterranean-Neolithic (доминируюший компонент  в западной и юго-западной Европе)

Контурные карты прочих компонентов можно скачать здесь.