Ganzgenomsequenzierung von 128 Kamelen in Asien offenbart Herkunft und Migration der einheimischen Trampeltiere

Probensammlung und Ganzgenomsequenzierung

Insgesamt 105 einheimische Trampeltiere in Asien, 19 wilde Trampeltiere aus der Gobi-Altai-Region in MG sowie 4 Dromedare aus dem Iran wurden für diese Studie gesammelt (Ergänzende Abb. 1 und ergänzende Tabelle 1). Die einheimischen Trampeltiere wurden so ausgewählt, dass sie möglichst viele große geografische Regionen abdecken, darunter 55 aus Inner MG (IMG), Xinjiang (XJ) und Qinghai in China, 28 aus MG, 6 aus KAZA, 10 aus Russland (RUS) und 6 aus dem IRAN. Da sich durch die breite Nutzung von Trampeltieren in China und MG eine Vielzahl lokaler Rassen gebildet hat, wurden acht verschiedene repräsentative Rassen aus den Regionen ausgewählt. Die anderen baktrischen Hauskamele aus Zentralasien lebten rund um das Kaspische Meer.

Nach der DNA-Extraktion wurden die einzelnen Genome mit einer durchschnittlichen 13-fachen Abdeckung sequenziert (siehe ergänzende Abb. 2 und ergänzende Tabelle 2). Die Sequenzlesungen wurden an unsere frühere Genomassemblierung des Trampeltiers29 angeglichen, um Varianten zu erkennen. Nach einer strengen Filterung (s. Abb. 3) konnten wir insgesamt 13,83 Millionen Einzelnukleotid-Polymorphismen (SNPs) und 1,41 Millionen kleine Indels identifizieren. Bemerkenswert ist, dass das Verhältnis von Transition zu Transversion bei den Roh-SNPs (2,29) niedriger war als bei Dromedaren (2,31-2,34)31, aber nach den Filterverfahren auf 2,44 anstieg, was auf eine Qualitätsverbesserung der identifizierten SNPs schließen lässt. Die funktionelle Annotation der Varianten zeigte, dass etwa 63,10 % von ihnen intergenisch, 33,62 % intronisch und 0,94 % exonisch waren (ergänzende Tabelle 3). Es wurden 13,73 Mio., 6,39 Mio. bzw. 10,55 Mio. Varianten bei baktrischen Hauskamelen, wilden baktrischen Kamelen und Dromedaren identifiziert. Obwohl sich die Dromedare phylogenetisch stärker von den beiden baktrischen Kamelarten unterschieden, teilten die baktrischen Hauskamele mehr Varianten mit den Dromedaren (66,73 %) als mit den wilden baktrischen Kamelen (39,31 %) (ergänzende Abb. 4), was auf die enorme Verringerung der genetischen Varianten bei den wilden baktrischen Kamelen und auf den Genfluss zwischen Dromedaren und baktrischen Hauskamelen zurückzuführen ist. Bei den baktrischen Hauskamelen wurden in den ostasiatischen und zentralasiatischen Populationen 12,68 Millionen bzw. 11,61 Millionen Varianten identifiziert (ergänzende Abb. 4). Obwohl die Anzahl der Hauskamele aus Ostasien höher war als die aus Zentralasien, zeigte die Anzahl der Varianten in den einzelnen Populationen keine signifikante Abweichung zwischen den beiden Gebieten (P-Wert = 0,77, zweiseitiger t-Test; ergänzende Tabelle 4).

Genetische Diversität und Differenzierung

Für einen detaillierteren Vergleich der genetischen Diversität zwischen den verschiedenen Populationen entfernten wir zunächst 14 Individuen, die eine enge genetische Verwandtschaft mit den übrigen zeigten (ergänzende Tabelle 5). Die paarweise Nukleotiddiversität π (Abb. 1a) der Dromedare (1,54 × 10-3) war signifikant höher als die der Trampeltiere aus allen geografischen Regionen (0,88 × 10-3-1,11 × 10-3; ergänzende Tabelle 6), was im Gegensatz zu früheren Schätzungen der Heterozygosität auf der Grundlage einzelner Genome stand10. Ein wichtiger Grund dafür könnte die Hybridisierungspraxis zwischen Dromedaren und Trampeltieren in Zentralasien sein30. Unter den Trampeltieren wies die Wildpopulation im Vergleich zu allen Hauspopulationen das niedrigste π (0,88 × 10-3) auf (Abb. 1a und ergänzende Tabelle 6). Dieses Ergebnis widerspricht zwar in vielen Fällen der Tatsache, dass Wildtiere in der Regel eine höhere genetische Vielfalt aufweisen als ihre domestizierten Artgenossen, wie z. B. Hunde24, Schweine27 und Kaninchen32, doch würde dies bei gefährdeten Wildtieren mit einer extrem kleinen Populationsgröße wie z. B. Pferden33 der Fall sein. Darüber hinaus wiesen die in Zentralasien lebenden Hauskatzenpopulationen im Allgemeinen eine höhere Diversität (1,03 × 10-3-1,11 × 10-3) auf als die in Ostasien lebenden (0,95 × 10-3-1,02 × 10-3; Abb. 1a). Diese Tendenz wurde auch durch den Watterson’s θ bestätigt (ergänzende Abb. 5). Auch hier könnte die Hybridisierung mit Dromedaren in Zentralasien für die höhere Diversität in dieser Region verantwortlich sein.

Wir haben dann die paarweise genetische Distanz zwischen den Kamelpopulationen mit Hilfe von Weir’s Fst (Abb. 1b) gemessen. Das Ergebnis stimmte gut mit der bekannten Phylogenie überein, die zeigte, dass die Dromedare die höchste Fst mit den Trampeltieren (0,54-0,64) und die wilden Trampeltiere die zweithöchste Fst mit den Hauskamelen (0,27-0,31) aufwiesen. Die Differenzierung zwischen den baktrischen Hauskamelen war wesentlich geringer, was auf ihren erst kürzlich erfolgten einheitlichen Ursprung schließen lässt. Interessanterweise wiesen unter den einheimischen Trampeltieren die Kamele aus dem IRAN die größte Divergenz zu allen anderen auf (0,05-0,06). Zur Validierung der Populationsdifferenzierung erstellten wir einen NJ-Baum für alle Individuen auf der Grundlage ihrer paarweisen Identität-nach-Zustand-Matrix (IBS) (siehe ergänzende Abbildung 6). Der NJ-Baum unterstützte ebenfalls eine monophyletische Klade aller baktrischen Hauskamele, innerhalb derer IRAN den tiefsten Zweig bildete.

Ein potenzielles Problem bei den Populationsgenomschätzungen war der Referenzbias, bei dem die Verwendung eines einzigen Referenzgenoms zu einer geringen Effizienz bei der Variantenbestimmung für Individuen führen würde, die sich stark davon unterscheiden34. Um die Verzerrung zu untersuchen, verglichen wir die fehlende Anzahl von Varianten zwischen den drei Arten, wobei wir die Sequenzierungstiefe als Kovariate verwendeten (ergänzende Tabelle 7). Die Varianzanalyse (ANOVA) zeigte, dass die domestizierten Trampeltiere zwar keinen signifikanten Unterschied zu den Wildkamelen aufwiesen (P-Wert = 0,50), aber eine geringere Anzahl fehlender Varianten als die Dromedare (P-Wert = 4,38 × 10-3). Um die Auswirkung der Verzerrung auf unsere Schätzungen zu bewerten, haben wir die genetische Vielfalt und Fst nur mit synonymen SNPs neu berechnet (ergänzende Abb. 7), da kodierende Sequenzen mit größerer Wahrscheinlichkeit artenübergreifend invariant sind. Im Ergebnis stimmten die auf den synonymen SNPs basierenden Schätzungen für alle Arten gut mit dem Gesamtgenom überein, was darauf hindeutet, dass die Referenzverzerrung nur geringe Auswirkungen auf unsere genomischen Populationsschätzungen hatte.

Populationsstruktur mit Vermischung

Um die Gesamtpopulationsstruktur mit potenzieller Vermischung aufzuzeigen, haben wir die SNPs beschnitten, indem wir diejenigen mit hohem Kopplungsungleichgewicht und potenziellen funktionellen Auswirkungen entfernt haben. Die multidimensionale Skalierungsanalyse (MDS) auf der Grundlage der beschnittenen Untergruppe ergab ähnliche Ergebnisse wie der vollständige Satz (Abb. 2a und ergänzende Abb. 8). Wie erwartet konnten die Dromedare und wilden Trampeltiere durch die erste bzw. zweite Koordinate getrennt werden. Als die dritte Koordinate in die MDS aufgenommen wurde, konnte IRAN von allen anderen Trampeltieren getrennt werden (Abb. 2a).

Um unterschiedliche Abstammungsanteile zu schätzen, führten wir eine Populationsstrukturanalyse mit Admixture35 durch, indem wir K Abstammungspopulationen annahmen (Abb. 2b). Das Kreuzvalidierungsverfahren bestätigte, dass K = 3 optimal ist (ergänzende Abb. 9), was eine klare Aufteilung zwischen Dromedaren, wilden Trampeltieren und baktrischen Hauskamelen zeigt. Zumindest bei einem Dromedar wurde eine offensichtliche Introgression von baktrischen Hauskamelen in die iranischen Dromedare beobachtet. Dementsprechend war die Abstammung von Dromedaren in den zentralasiatischen Populationen baktrischer Kamele, einschließlich IRAN, KAZA und RUS, mit einem geschätzten Anteil von 1-10 % weit verbreitet. Darüber hinaus beobachteten wir bei einigen wildlebenden Tieren eine Abstammung von baktrischen Hauskamelen mit einem Anteil von 7-15 %. Dies könnte auf Polymorphismen der Vorfahren zurückzuführen sein, aber auch auf introgressive Hybridisierung, die bei mtDNAs36 und Y-Chromosomen15 beobachtet wurde und die genetische Besonderheit der wilden Arten bedrohen könnte. Überraschenderweise trugen die wilden Kamele so gut wie nichts zur Abstammung der einheimischen Populationen bei, selbst nicht bei den mongolischen Populationen, die enge Lebensräume mit den wilden Kamelen teilen. Obwohl die meisten baktrischen Hauskamele keine Differenzierung aufwiesen, als K wuchs, war IRAN die erste Population, die sich mit einer eindeutigen Abstammung abspaltete (K = 5; Abb. 2b).

Als weitere Methode zur Untersuchung der Populationsstruktur mit Vermischung haben wir den Populationsbaum für die Kamele mit TreeMix37 abgeleitet (Abb. 2c). Wenn keine Migrationsspur hinzugefügt wurde, zeigte die Baumtopologie erneut, dass IRAN die erste Population war, die von allen einheimischen Trampeltieren getrennt wurde (ergänzende Abb. 10). Die Erhöhung der Migrationsspuren (m = 1-3) konnte die Anpassung des Modells erheblich verbessern (ergänzende Abb. 11), wodurch Genflüsse von Dromedaren zu baktrischen Hauskamelen in Zentralasien, einschließlich KAZA, RUS und IRAN, mit Migrationsgewichten zwischen 4 % und 9 % identifiziert wurden (ergänzende Tabelle 8). Es ist erwähnenswert, dass die Migrationsroute zwar am Ende des Dromedarzweigs in den IRAN führte, aber in der Mitte des Dromedarzweigs nach KAZA und RUS (Abb. 2c). Dies könnte auf eine mit dem iranischen Dromedar verwandte Geisterpopulation hindeuten, die zur Abstammung von KAZA und RUS beigetragen hat. Eine zusätzliche Migrationsspur (m = 4) konnte die Anpassung des Modells weiter verbessern, was auf eine Migration des Dromedars zu einer XJ-Rasse hindeutet (Abb. 2c). Obwohl TreeMix kein starkes Signal für eine Migration zwischen den wilden und den domestizierten Trampeltieren feststellte, zeigten die Rückstände eine mäßige Vermischung zwischen den wilden und den ostasiatischen Rassen (ergänzende Abb. 11). Anschließend verwendeten wir den weniger parametrisierten Drei- und Vier-Populations-Test (F3/F4)38 , um die statistische Signifikanz dieser Vermischungsereignisse zu bewerten. Auch hier bestätigte der F3-Test die Vermischung von Dromedaren und Trampeltieren in KAZA, RUS und IRAN (ergänzende Tabelle 9). Der empfindlichere F4-Test bestätigte ein signifikant höheres Ausmaß der Vermischung zwischen Dromedaren und Trampeltieren in Zentralasien im Vergleich zu denen in Ostasien (ergänzende Tabelle 10). Unter den letzteren wurde ein höheres Ausmaß an Vermischung mit Dromedaren in XJ als in MG/IMG festgestellt.

Beweise für zentralasiatischen Ursprung durch Entfernen der Introgression

Ost- und Zentralasien waren die beiden alternativen Domestikationsregionen für baktrische Kamele, basierend auf archäologischen Beweisen1,12,17, aber die wahrscheinlichste Region blieb ungeklärt. Obwohl wir die größte genetische Differenzierung zwischen der iranischen Population und allen anderen einheimischen Populationen beobachteten, würde die Existenz von Vermischungen zwischen Dromedaren und baktrischen Kamelen in Zentralasien die Rückschlüsse auf die Herkunft schwächen. Um diesen Effekt zu verringern, haben wir versucht, die introgressierten Segmente der Dromedare aus den Genomen der baktrischen Kamele mit Hilfe des „BABA/ABBA“-Tests39 zu entfernen. Wir gruppierten die ost- und zentralasiatischen Populationen und verglichen die gemeinsame Nutzung von Allelen zwischen den beiden Gruppen und Dromedaren (Abb. 3a). Da die Abstammung von baktrischen Kamelen bei einem Dromedar und die Abstammung von Hauskamelen bei drei Wildkamelen (Abb. 2b) Störfaktoren darstellen würden, haben wir diese vier Individuen aus der Analyse herausgenommen. Wir verwendeten die Statistik fd, eine robuste Version des Patterson’s D, um introgressierte Segmente zu lokalisieren40 und wendeten ein strenges Signifikanzniveau von Z-score = 2 an, indem wir die Jackknife-Prozedur verwendeten (ergänzende Abb. 12). In den insgesamt 21.153 nicht überlappenden 100 kb-Segmenten des Genoms gab es in den zentralasiatischen Populationen (11.711, Z-score > 2) erwartungsgemäß weit mehr Segmente mit mutmaßlichen Signalen der Introgression als in den ostasiatischen Populationen (3891, Z-score < -2). Anschließend führten wir die Admixture-Analyse auf der Grundlage der verbleibenden Segmente durch und bestätigten, dass die Introgression der Dromedare effektiv reduziert wurde (ergänzende Abb. 13). Die Neuberechnung der paarweisen Fst nach Entfernung der Introgression zeigte immer noch, dass IRAN die am stärksten differenzierte (0,04-0,06) aller einheimischen Populationen war (Abb. 3b). Um weitere Einblicke in die Phylogenie der Population zu erhalten, rekonstruierten wir den NJ-Baum auf der Grundlage der paarweisen Fst und führten den Bootstrap-Test durch (Abb. 3b). Es bestätigte sich, dass IRAN die erste Population war, die sich von allen inländischen Bactrian-Populationen trennte, gefolgt von KAZA und RUS. Die Bayes’sche binäre Markov-Chain-Monte-Carlo (MCMC)-Analyse, die auf der Phylogenie basiert, stützt sich stark auf Zentralasien als Herkunftsgebiet der baktrischen Hauskamele (Wahrscheinlichkeit = 99,78 %) und eine anschließende Ausbreitungsroute von Zentral- nach Ostasien (ergänzende Abb. 14).

Als unabhängigen Nachweis rekonstruierten wir auch den Maximum-Likelihood-Baum der mtDNAs in voller Länge auf der Grundlage der 128 Proben, die wir in dieser Studie sequenziert haben, sowie von 39 zusätzlichen Proben, die in Genbank verfügbar sind (Abb. 3c und ergänzende Tabelle 11). Die Introgression von mtDNAs konnte leicht identifiziert und aus dem Baum ausgeschlossen werden. Zum Beispiel wurden zwei neu sequenzierte mtDNAs aus KAZA und RUS mit Dromedaren geclustert. Innerhalb der Klade der baktrischen Hauskamele waren zwar die meisten Kamele aus verschiedenen geografischen Regionen gemischt, aber zwei mtDNAs aus dem IRAN bildeten die basalsten Äste der Hauspopulationen (Abb. 3c). Die Bayes’sche binäre MCMC-Analyse bestätigte erneut die zentralasiatische Herkunft der baktrischen Hauskamele (Wahrscheinlichkeit = 76,43 %).

Demografische Geschichte der baktrischen Kamele

Wir führten mehrere parametrische Modellierungsanalysen durch, um die demografische Dynamik der Kamele in der Geschichte zu ermitteln. In Übereinstimmung mit früheren Studien10 zeigten die Langzeitverläufe der baktrischen Kamele auf der Grundlage des paarweise sequentiell markovianischen Koaleszenzmodells (PSMC)41 einen enormen Rückgang der effektiven Populationsgröße der Urkamele vor mehr als einer Million Jahren (ergänzende Abb. 15). Obwohl die langfristigen Entwicklungen der wilden und der domestizierten baktrischen Kamele im Allgemeinen ähnlich waren, wichen sie offensichtlich bereits vor 0,4 Millionen Jahren voneinander ab, was die ersteren als direkte Vorfahren der letzteren ausschließt, wie frühere mtDNA-Analysen zeigen9,14.

Um die Divergenzzeit zwischen den Kamelpopulationen genauer zu untersuchen, verwendeten wir den generalisierten phylogenetischen Koaleszenz-Sampler (G-PhoCS)42. Angesichts der Phylogenie der Kamelpopulationen konnte G-PhoCS die mutationsskalierte Populationsgröße und die Divergenzzeit der Populationen auf der Grundlage der nicht verknüpften neutralen Loci in den einzelnen Genomen der einzelnen Populationen schätzen. Um die Komplexität des Modells zu reduzieren, haben wir nur Dromedare, wilde Trampeltiere und drei repräsentative Populationen (IRAN, KAZA und MG) von Trampeltieren einbezogen (ergänzende Abb. 16 und ergänzende Tabelle 12). Nach Abb. 3b waren IRAN und KAZA die ersten beiden zentralasiatischen Populationen, die sich trennten, und die Aufspaltung von MG könnte auf die Ausbreitung von Zentral- nach Ostasien hinweisen. Das Alter wurde kalibriert, indem die Divergenz zwischen baktrischem Kamel und Dromedar mit 5,73 Millionen Jahren gemäß der Timetree-Datenbank43 angenommen wurde. Wenn kein Migrationsband einbezogen wurde, konnte problemlos eine Konvergenz aller Parameterschätzungen erreicht werden (ergänzende Abb. 17 und ergänzende Tabelle 13). Ähnlich wie bei den PSMC-Ergebnissen verringerte sich die effektive Populationsgröße generell von den angestammten zu den modernen Populationen (Abb. 4). Der Zeitpunkt der Divergenz zwischen Wild- und Hauskamelen wurde auf 0,43 Mio. Jahre geschätzt (95 % Konfidenzintervall: 0,13-0,73 Mya; Abb. 4), also etwas später als der auf mtDNA basierende Wert (0,714 oder 1,1 Mya9). Unter den einheimischen Populationen wurde der IRAN vor etwa 4,45 Tausend Jahren (95% CI: 0,07-17,6 Kya) von den anderen getrennt, und dann wurden die zentral- und ostasiatischen Populationen vor etwa 2,40 Tausend Jahren (95% CI: 0,01-7,84 Kya; Abb. 4) getrennt.

Um den Genfluss zu berücksichtigen, haben wir auch versucht, Migrationsbänder von Dromedaren zu baktrischen Kamelpopulationen einzuführen (ergänzende Abb. 16 und ergänzende Tabelle 12). Die Schätzungen konnten nur konvergieren, wenn ein Migrationsband von iranischen Dromedaren nach IRAN und ein Migrationsband von einer Geisterpopulation nach KAZA eingeführt wurden (ergänzende Abb. 18 und ergänzende Tabelle 13). Obwohl die Divergenzzeit zwischen wilden und domestizierten Trampeltieren durch das Migrationsmodell nicht verändert wurde (0,46 Mya, 95 % CI: 0,24-0,71 Mya), wurde die erste Divergenzzeit der domestizierten Populationen (0,19 Mya, 95 % CI: 0,08-0,31 Mya) unrealistisch, da sie weit jenseits der bekannten Geschichte der Domestikation von Nutztieren lag (11,5 Kya44). Außerdem wurde die Gesamtmigrationsrate für das Migrationsband nach IRAN und KAZA nur auf 0,27 % bzw. 0,16 % geschätzt, was weit unter der mit Admixture geschätzten Rate (1-10 %) liegt. Ein möglicher Grund für die schlechte Schätzung könnte eine komplexere Vermischungsgeschichte sein als das von G-PhoCS angenommene kontinuierliche Migrationsmodell mit konstanten Raten.