Mintagyűjtés és teljes genomszekvenálás
Az összesen 105 házi baktriai teve Ázsiában, 19 vadon élő baktriai tevét a Góbi-Altai régióból MG-ből, valamint 4 dromedárt IRÁNból gyűjtöttünk össze a tanulmányhoz (Kiegészítő ábra. 1. és 1. kiegészítő táblázat). A házi baktriai tevéket úgy választottuk ki, hogy a lehető legtöbb nagy földrajzi régiót lefedjék, köztük 55-öt a kínai Belső-MG-ből (IMG), Hszincsiangból (XJ) és Qinghaiból, 28-at MG-ből, 6-ot KAZA-ból, 10-et Oroszországból (RUS) és 6-ot IRÁNból. Mivel Kínában és MG-ben a háziasított baktriai tevék széles körű felhasználása miatt sokféle helyi fajta alakult ki, a régiókból nyolc különböző reprezentatív fajtát választottak ki. A többi közép-ázsiai házi baktriai teve a Kaszpi-tenger környékén élt.
A DNS extrakció után az egyes genomokat átlagosan 13×-os lefedettséggel szekvenálták (2. kiegészítő ábra és 2. kiegészítő táblázat). A szekvenciaolvasásokat a variánsok megnevezéséhez a bakteriális teve korábbi genom-összeállításához29 igazítottuk. Szigorú szűrés után (3. kiegészítő ábra) összesen 13,83 millió egynukleotid-polimorfizmust (SNP) és 1,41 millió kis indelt azonosítottunk. Figyelemre méltó, hogy a nyers SNP-k átmenet-transzverzió aránya (2,29) alacsonyabb volt, mint a dromedárok esetében (2,31-2,34)31 , de a szűrési eljárások után 2,44-re nőtt, ami az azonosított SNP-k minőségének javulására utal. A variánsok funkcionális annotációja azt mutatta, hogy körülbelül 63,10%-uk intergenikus, 33,62%-uk intronikus és 0,94%-uk exonikus volt (3. kiegészítő táblázat). A házi baktriai tevéknél 13,73 millió, a vadon élő baktriai tevéknél 6,39 millió, illetve a dromedároknál 10,55 millió variánst azonosítottak. Bár a dromedárok filogenetikailag jobban eltértek mindkét baktriai tevefajtól, a házi baktriai tevék több változatot osztottak meg a dromedárokkal (66,73%), mint a vad baktriai tevékkel (39,31%) (4. kiegészítő ábra), ami a genetikai változatoknak a ma élő vad baktriai tevéknél megfigyelt óriási csökkenésének, valamint a dromedárok és a házi baktriai tevék közötti génáramlásnak köszönhető. A házi baktriai tevék körében 12,68 millió, illetve 11,61 millió variánst azonosítottak a kelet-ázsiai és a közép-ázsiai populációkban (4. kiegészítő ábra). Bár a kelet-ázsiai házi tevékből több mintát vettünk, mint a közép-ázsiaiakból, az egyes populációkhoz tartozó saját variánsok száma nem mutatott jelentős eltérést a két terület között (P-érték = 0,77, kétfarkú t-próba; 4. kiegészítő táblázat).
Genetikai diverzitás és differenciálódás
A különböző populációk közötti genetikai diverzitás részletesebb összehasonlításához először eltávolítottunk 14 egyedet, amelyek szoros genetikai kapcsolatot mutattak a többi egyeddel (5. kiegészítő táblázat). A dromedárok páros nukleotiddiverzitása π (1a. ábra) (1,54 × 10-3) jelentősen magasabb volt, mint az összes földrajzi régióból származó baktriai tevéké (0,88 × 10-3-1,11 × 10-3; 6. kiegészítő táblázat), ami ellentétben állt a korábbi, egyedi genomokon alapuló heterozigozitási becslésekkel10. Ennek egyik fontos oka lehet a dromedárok és a baktriai tevék közötti hibridizációs gyakorlat Közép-Ázsiában30. A baktriai tevék közül a vadon élő populáció mutatta a legalacsonyabb π-t (0,88 × 10-3) az összes hazai populációhoz képest (1a. ábra és 6. kiegészítő táblázat). Bár ez az eredmény sok esetben sérti azt, hogy a vadon élő állatok általában nagyobb genetikai diverzitással rendelkeznek, mint házi társaik, mint például a kutyák24, sertések27 és nyulak32 , ez a veszélyeztetett, rendkívül kis populációmérettel rendelkező vadon élő állatok, például a lovak33 esetében fordulhat elő. Ezenkívül a Közép-Ázsiában élő háziállat-populációk általában nagyobb diverzitást mutattak (1,03 × 10-3-1,11 × 10-3), mint a Kelet-Ázsiában élők (0,95 × 10-3-1,02 × 10-3; 1a. ábra). Ezt a tendenciát a Watterson θ is megerősítette (5. kiegészítő ábra). Ismét a közép-ázsiai dromedárokkal való hibridizáció magyarázhatja a régió magasabb diverzitását.
a Nukleotiddiverzitás π. A boxplot a π-t mutatja 2,0 × 105 10 kb-os csúszóablakokra a genomban. Az egyes populációk földrajzi eredete és mintanagysága a bal oldalon, a π átlaga pedig a jobb oldalon látható. Az MG, IMG és XJ esetében több helyi fajtából vettek mintát. A szoros genetikai kapcsolatokkal rendelkező egyedeket eltávolítottuk. A boxplot elemei a következőképpen vannak meghatározva: középvonal, medián; dobozhatárok, a harmadik és az első kvartilis; whiskerek, 1,5 × interkvartilis tartomány. b Páros populációs differenciálódás Fst. A hőtérkép 2,0 × 105 10 kb-os csúszóablakok átlagos Fst értékét mutatja. A dendrogram a populációk hierarchikus klaszterezését mutatja az Fst alapján.
A továbbiakban a tevepopulációk közötti páros genetikai távolságot mértük a Weir-féle Fst segítségével (1b. ábra). Az eredmény jól megfelelt az ismert filogenetikának, amely szerint a dromedároknak volt a legnagyobb Fst-jük a baktriai tevékkel (0,54-0,64), a vadon élő baktriai tevéknek pedig a második legnagyobb Fst-jük a házi tevékkel (0,27-0,31). A házi baktriai tevék közötti differenciálódás sokkal kisebb volt, ami összhangban van a közelmúltbeli egységes eredetükkel. Érdekes módon a házi baktriai tevék közül az IRÁN-ból származó tevék mutatták a legnagyobb eltérést az összes többi tevével (0,05-0,06). A populációs differenciálódás validálása érdekében az összes egyedre vonatkozóan szomszéd-összekötő (NJ) fát készítettünk a páronkénti azonosság-az-állapot (IBS) mátrix alapján (6. kiegészítő ábra). Az NJ-fa is támogatta az összes házi baktriai teve monofiletikus kládját, amelyen belül az IRAN képezte a legmélyebb ágat.
A populáció genomikai becsléseivel kapcsolatos lehetséges probléma a referencia torzítás volt, ahol egyetlen referencia genom használata alacsony hatékonyságot eredményezne a variánsok megnevezésében az attól nagymértékben eltérő egyedek esetében34. A torzítás vizsgálatához összehasonlítottuk a variánsok hiányzó számát a három faj között, a szekvenálási mélységet kovariátumnak tekintve (7. kiegészítő táblázat). A varianciaanalízis (ANOVA) azt mutatta, hogy bár a házi baktriai tevéknél nem volt szignifikáns különbség a vadon élő egyedekhez képest (P-érték = 0,50), a hiányzó variánsok száma valóban alacsonyabb volt, mint a dromedároknál (P-érték = 4,38 × 10-3). Hogy értékeljük a torzítás hatását a becsléseinkre, újra kiszámítottuk a genetikai diverzitást és az Fst-t csak szinonim SNP-kkel (7. kiegészítő ábra), mivel a kódoló szekvenciák nagyobb valószínűséggel változatlanok a fajok között. Ennek eredményeként a szinonim SNP-ken alapuló becslések minden faj esetében jó összhangban voltak a teljes genommal, ami arra utal, hogy a referencia torzításnak csak kisebb hatása volt a populációgenomikai becsléseinkre.
Populációszerkezet keveredéssel
A potenciális keveredéssel rendelkező teljes populációszerkezet feltárása érdekében az SNP-ket a magas kötési egyenlőtlenséggel és potenciális funkcionális hatással rendelkező SNP-k eltávolításával metszettük. A többdimenziós skálázási (MDS) elemzés a metszett részhalmazon alapulva hasonló eredményt adott, mint a teljes halmaz (2a. ábra és 8. kiegészítő ábra). A várakozásoknak megfelelően a dromedárokat és a vadon élő baktriai tevéket az első, illetve a második koordináták alapján lehetett elkülöníteni. Amikor a harmadik koordinátát is bevonták az MDS-be, az IRAN elkülönült az összes többi házi baktriai tevétől (2a. ábra).
a Többdimenziós skálázási (MDS) diagram az 1-4. koordinátával (C1-C4). b K különböző számú ősiséget feltételező keveredési elemzés. Az egyén genomjának az egyes ősökhöz rendelt arányát különböző színekkel ábrázoljuk. c TreeMix-elemzés különböző migrációs eseményeket feltételezve m. A migrációs súly a donorpopulációból kapott ősök aránya.
A különböző ősarányok becsléséhez populációszerkezeti elemzést végeztünk Admixture35 segítségével, K őspopulációt feltételezve (2b. ábra). A keresztvalidálási eljárás alátámasztotta, hogy a K = 3 optimális (9. kiegészítő ábra), ami egyértelmű felosztást mutat a dromedárok, a vadon élő baktriai tevék és a házi baktriai tevék között. A házi baktriai tevéknek az iráni dromedárokba való nyilvánvaló introgressziója volt megfigyelhető, legalábbis egy dromedár esetében. Ennek megfelelően a dromedár származás a közép-ázsiai baktriai tevepopulációkban, beleértve IRÁNt, KAZA-t és RUS-t is, 1-10%-ra becsült arányban volt elterjedt. Ezenfelül több vadon élő egyedben is megfigyeltük a házi baktriai tevék származását, 7-15%-os arányban. Ez eredhet ősi polimorfizmusokból, de okozhatja introgresszív hibridizáció is, amelyet az mtDNS36 és az Y kromoszómák15 esetében figyeltek meg, és amely a vadon élő fajok genetikai különállóságát veszélyeztetheti. Meglepő módon a vadon élő tevék szinte semmit sem járultak hozzá a hazai populációk őshazájához, még a mongol populációkhoz sem, amelyek közeli élőhelyeken osztoznak a vadon élő tevékkel. Bár a legtöbb házi baktriai teve nem differenciálódott, amikor K nőtt, az IRAN volt az első populáció, amely egyedi ősiséggel különült el (K = 5; 2b. ábra).
A keveredéssel járó populációszerkezet vizsgálatának másik módszereként a TreeMix37 segítségével következtettünk a tevék populációfájára (2c. ábra). Amikor nem adtunk hozzá migrációs nyomot, a fa topológiája ismét azt jelezte, hogy az IRAN volt az első elkülönült populáció az összes hazai baktriai teve között (10. kiegészítő ábra). A migrációs nyomvonalak növelésével (m = 1-3) jelentősen javítani lehetett a modell illeszkedését (11. kiegészítő ábra), amely a dromedárokból a házi baktriai tevékbe történő génáramlást azonosította Közép-Ázsiában, beleértve KAZA-t, RUS-t és IRÁN-t, a migrációs súlyok 4% és 9% között mozogtak (8. kiegészítő táblázat). Érdemes volt megemlíteni, hogy bár a migrációs nyomvonal a dromedár ág végén IRÁN felé mutatott, a dromedár ág közepén KAZA és RUS felé mutatott (2c. ábra). Ez egy, az iráni dromedárral rokon szellempopulációra utalhat, amely hozzájárult KAZA és RUS felmenőihez. További migrációs nyomvonal (m = 4) tovább javíthatta a modell illeszkedését, amely a dromedár XJ fajtába való vándorlását jelezte (2c. ábra). Bár a TreeMix nem észlelt erős migrációs jelet a vadon élő és a házi baktriai tevék között, a maradékok mérsékelt keveredést mutattak a vadon élő és a kelet-ázsiai fajták között (11. kiegészítő ábra). Ezután a kevésbé paraméterezett három- és négypopulációs (F3/F4) tesztet38 alkalmaztuk ezen keveredési események statisztikai szignifikanciájának értékelésére. Az F3 teszt ismét erősen alátámasztotta a dromedárok és a baktriai tevék keveredését KAZA-ban, RUS-ban és IRÁNban (9. kiegészítő táblázat). Az érzékenyebb F4 teszt megerősítette a dromedárok és a baktriai tevék közötti keveredés szignifikánsan nagyobb mértékét Közép-Ázsiában, mint Kelet-Ázsiában (10. kiegészítő táblázat). Az utóbbiak közül a dromedárokkal való keveredés nagyobb mértékű volt az XJ-ben, mint az MG/IMG-ben.
A közép-ázsiai eredet bizonyítéka az introgresszió eltávolításával
Kelet- és Közép-Ázsia volt a két alternatív háziasítási régió a baktriai tevék számára a régészeti bizonyítékok alapján1,12,17, de a legvalószínűbb megoldatlan maradt. Bár a legnagyobb genetikai differenciálódást az iráni populáció és az összes többi háziasított populáció között figyeltük meg, a dromedárok és a baktriai tevék közötti keveredés megléte Közép-Ázsiában gyengítené a származási következtetés alátámasztását. E hatás csökkentése érdekében a “BABA/ABBA” teszt39 segítségével megpróbáltuk eltávolítani a dromedárok introgresszált szegmenseit a baktriai tevék genomjából. Csoportosítottuk a kelet- és közép-ázsiai populációkat, és összehasonlítottuk a két csoport közötti allélmegosztást a dromedárokkal (3a. ábra). Mivel a baktriai tevék származása egy dromedárban, valamint a házi tevék származása három vad tevében (2b. ábra) zavaró tényező lenne, ezt a négy egyedet eltávolítottuk az elemzésből. Az fd statisztikát, a Patterson-féle D robusztus változatát használtuk az introgresszált szegmensek40 felkutatására, és a Jackknife-eljárás alkalmazásával szigorú Z-score = 2 szignifikancia-szintet alkalmaztunk (12. kiegészítő ábra). A genomban összesen 21 153 nem átfedő 100 kb-os szegmensben a közép-ázsiai populációkban sokkal több olyan szegmens volt, amely introgresszió feltételezett jeleit mutatta (11 711, Z-score > 2), mint a kelet-ázsiai populációkban (3891, Z-score < -2), ahogy az várható volt. Ezután elvégeztük az admixtúraelemzést a fennmaradó szegmensek alapján, és megerősítettük, hogy a dromedárok introgressziója hatékonyan csökkent (13. kiegészítő ábra). A páronkénti Fst újraszámítása az introgresszió eltávolítása után még mindig azt mutatta, hogy az IRAN volt a leginkább differenciált (0,04-0,06) az összes hazai populáció között (3b. ábra). Hogy további betekintést nyerjünk a populáció filogeniájába, rekonstruáltuk az NJ-fát a páros Fst alapján, és elvégeztük a bootstrap-tesztet (3b. ábra). Ez megerősítette, hogy az összes hazai baktriai populáció közül elsőként az IRAN különült el, majd a KAZA és a RUS következett. A filogenián alapuló Bayes-féle bináris Markov-lánc Monte Carlo (MCMC) elemzés határozottan alátámasztotta Közép-Ázsiát mint a házi baktriai tevék őshazáját (valószínűség = 99,78%) és a Közép-Ázsiából Kelet-Ázsiába történő későbbi elterjedési útvonalat (14. kiegészítő ábra).
a BABA/ABBA elemzés a dromedároknak a házi baktriai tevékbe való introgressziójára. Hogy erre az introgresszióra összpontosítsunk, egy dromedárt a baktriai tevék ősével és három vad tevét a házi tevék ősével távolítottunk el. A szignifikáns fd-vel (|Z-score| > 2) rendelkező 100 kb szegmensek száma az egyes fa konfigurációk esetében látható. b A populációk szomszéd-összekötő (NJ) fája az introgressziós szegmensek eltávolítása után. A hőtérkép az 5,1 × 104 10 kb-os csúszóablakra vonatkozó átlagos páros Fst értékeket mutatja. Az NJ-fa Bootstrap-értékeit ötezer 10 kb-os ablak 100-szoros véletlenszerű mintavételezésével számoltuk. c A teljes hosszúságú mtDNS-ek maximális valószínűségű fája. A populációkat különböző színekkel, a Genbankból származó szekvenciákat pedig pontokkal jelöljük. A főágak Bootstrap-értékei fel vannak tüntetve.
Független bizonyítékként a teljes hosszúságú mtDNS-ek maximális valószínűségű fáját is rekonstruáltuk a jelen tanulmányban szekvenált 128 minta, valamint a Genbankból elérhető további 39 minta alapján (3c. ábra és 11. kiegészítő táblázat). Az mtDNS-ek introgressziója könnyen azonosítható és kizárható volt a fából. Például két újonnan szekvenált mtDNS a KAZA-ból és a RUS-ból a dromedárokkal klasztereződött. A házi baktriai tevék kládján belül, bár a különböző földrajzi régiókból származó tevék többsége keveredett, két IRÁN-ból származó mtDNS képezte a hazai populációk legbazálisabb ágát (3c. ábra). A Bayes-féle bináris MCMC-elemzés ismét alátámasztotta a házi baktriai tevék közép-ázsiai eredetét (valószínűség = 76,43%).
A baktriai tevék demográfiai története
A tevék történeti demográfiai dinamikájára való következtetés érdekében több parametrikus modellezési elemzést végeztünk. A korábbi tanulmányokkal10 összhangban a baktriai tevék hosszú távú pályája a páronkénti szekvenciális Markov-koaleszcencia (PSMC)41 modell alapján az ős tevék effektív populációméretének óriási csökkenését mutatta ki korábban, mint egymillió évvel ezelőtt (15. kiegészítő ábra). Bár a vadon élő és a házi baktriai tevék hosszú távú pályája általában hasonló volt, nyilvánvaló volt, hogy már 0,4 millió évvel ezelőtt eltértek egymástól, ami kizárja, hogy az előbbiek az utóbbiak közvetlen őseit képezték volna, ahogyan azt a korábbi mtDNS-elemzések9,14 is tették.
A tevepopulációk közötti divergenciaidő részletesebb feltárásához az általánosított filogenetikai koaleszcens mintavételezőt (G-PhoCS)42 használtuk. A tevepopulációk filogeniáját figyelembe véve a G-PhoCS képes volt megbecsülni a mutációval skálázott populációméretet és a populációk divergenciaidejét az egyes populációk egyes genomjaiban található, nem kapcsolt semleges lókuszok alapján. A modell komplexitásának csökkentése érdekében csak a dromedárokat, a vadon élő baktriai tevéket és a házi baktriai tevék három reprezentatív populációját (IRAN, KAZA és MG) vettük figyelembe (16. kiegészítő ábra és 12. kiegészítő táblázat). A 3b. ábra szerint az IRAN és a KAZA volt az első két közép-ázsiai populáció, amely elkülönült, az MG kettéválása pedig a Közép-Ázsiából Kelet-Ázsiába történő szétszóródásra utalhat. A kort a Timetree adatbázis43 szerinti 5,73 millió éves baktriai teve és dromedár divergenciát feltételezve kalibráltuk. Ha nem építettünk be migrációs sávot, az összes paraméterbecslés konvergenciája könnyen elérhető volt (17. kiegészítő ábra és 13. kiegészítő táblázat). A PSMC-eredményekhez hasonlóan az effektív populációméret általában csökkent az ősi populációktól a modern populációkig (4. ábra). A vadon élő és a házi baktriai tevék közötti divergencia idejét 0,43 millió évvel ezelőttre becsültük (95%-os konfidenciaintervallum : 0,13-0,73 Mya; 4. ábra), ami valamivel későbbi, mint az mtDNS alapján becsült idő (0,714 vagy 1,1 Mya9). A hazai populációk közül az IRAN körülbelül 4,45 ezer évvel ezelőtt (95%-os CI: 0,07-17,6 Kya), majd a közép- és kelet-ázsiai populációk körülbelül 2,40 ezer évvel ezelőtt (95%-os CI: 0,01-7,84 Kya; 4. ábra) különültek el a többiektől.
A mutációval skálázott hatékony populációméret θ változását a hőszínek mutatják. Az években kifejezett időt a dromedárok és a baktriai tevék közötti divergenciaidővel kalibráltuk. A 95%-os konfidenciaintervallumokat az időtengelyen elhelyezett sávok mutatják. A piros és a kék sáv az IRAN-MG, illetve a KAZA-MG divergenciát jelzi. Ezek a becslések a migráció nélküli modellen alapulnak.
A génáramlás figyelembevétele érdekében megpróbáltunk migrációs sávokat is bevezetni a dromedárokból a baktriai tevepopulációkba (16. kiegészítő ábra és 12. kiegészítő táblázat). A becslések csak akkor tudtak konvergálni, ha egy iráni dromedárokból IRÁN-ba irányuló migrációs sávot és egy szellempopulációból KAZA-ba irányuló migrációs sávot vezettünk be (18. kiegészítő ábra és 13. kiegészítő táblázat). Bár a vadon élő és a háziasított baktriai tevék közötti divergenciaidő nem változott a migrációs modellel (0,46 Mya, 95% CI: 0,24-0,71 Mya), a háziasított populációk első divergenciaideje (0,19 Mya, 95% CI: 0,08-0,31 Mya) irreálissá vált, mert messze meghaladta a haszonállatok háziasításának ismert történetét (11,5 Kya44). Emellett a teljes migrációs arányt csak 0,27%-ra, illetve 0,16%-ra becsülték az IRÁN és a KAZA felé irányuló migrációs sáv esetében, ami jóval alacsonyabb, mint a keveredéssel becsült érték (1-10%). A gyenge becslés egyik lehetséges oka az lehet, hogy a G-PhoCS által feltételezett állandó arányú folyamatos migrációs modellnél bonyolultabb keveredési előzményekről van szó.