We’re all Siberians

Counterintuitive as it seems to many, computer models have shown that, with a very high probability, the last common ancestor of all living humans has lived surprisingly recently, maybe as little as 3500 years ago. More recently, genomic testing has helped to validate some of the assumptions underlying these models, showing that a pair of randomly picked individuals from e.g. Hungary and UK share about five *genetic* ancestors as little as ~1000 years back, with translates readily to thousands, if not millions of shared genealogical ancestors.

What I’m doing here is to make an even more minimalistic model: One that does without explicitly modeling large-scale migrations. Instead, all I model is admixture between neighbouring villages. And even with such an impoverished model, pretty much everyone shares ancestors by around generation 160.

East_Siberia67

Red is villages where everyone has ancestors who lived in East Siberia at the beginning of the simulation. Yellow, villages where some people do.

 

A stochastic process determines the distance of the partner, in such a way that most of the time, they’ll be from the same village or its immediate neighbours. In the current run I used the formula `int(2**abs(random.gauss(0.0,1.5))-1)` which (with a distance limit of 20) translates to roughly the following distribution:

percentage  distance cumulative percentage (distance >= n)
49.7226908118 0 100.0
21.4748636521 1 50.2773091882
10.8451461364 2 28.8024455361
6.10030067738 3 17.9572993996
3.66890586732 4 11.8569987223
2.33591620348 5 8.18809285495
1.58101339709 6 5.85217665147
1.09519578129 7 4.27116325438
0.781781237651 8 3.17596747309
0.571969005558 9 2.39418623543
0.423234999729 10 1.82221722988
0.326552881308 11 1.39898223015
0.244112319749 12 1.07242934884
0.198378869541 13 0.828317029091
0.157760344686 14 0.62993815955
0.122557623144 15 0.472177814864
0.0946762653423 16 0.349620191719
0.0866528530251 17 0.254943926377
0.0648893471149 18 0.168291073352
0.0547597890645 19 0.103401726237
0.0486419371727 20 0.0486419371727

In words, 80% of matings happen within a distance of less than 3,less than 2% pick a mate 10 or more steps away, and the theoretical maximum of 20 steps only becomes relevant about once in 2000 matings.

For computational reason, I’m working with a limited population size. The effect of this is to further slow the spread of shared ancestry: With the same rates, a village of 1000 will receive gene flow from 20 steps away once every other generation on average, while a village of 67 will do so once every 30 generations. Since, once a certain heritage has established itself, its only a question of a few generations before everyone in the village has it, this effectively slows the spread.

But wait, you’ll say (people always do in comment sections whenever this comes up in the news), some communities are pretty isolated, what about them?

There’s a method for that too: The Village class has a parameter exogamy_taboo that allows to set an arbitrarily low likelihood of outside mating actually happening: If the above distance formula throws up anything other than 0, another random number is generated and, unless that second number is below a certain threshold, the distance calculation is repeated to hopefully yield 0 this time.

In the iteration presented, such a taboo against exogamy was present in all villages of Shangri-La, the entire Southern half of Papua, and two thirds of all Alaskan villages, with default strength (~98% of matings inside the village, the rest distributed as per above). A stronger taboo is present in every other village of the Southern India and a ridiculously strong one in one in four West African villages (the weak taboo in the Carribean, East Siberia, Southern Africa and the remainder of Alaska only slightly slows things down and doesn’t really deserve to be called a taboo):

for village in shangri_la + papua[:2 * len(papua) // 3]+random.sample(alaska,2 * len(alaska) // 3):
 village.set_exogamy_taboo(True)

for village in alaska + carribean + east_siberia + southern_africa:
 village.set_exogamy_taboo(True, 0.3)

for village in random.sample(western_africa, len(western_africa)//4):
 village.set_exogamy_taboo(True, 0.005)

for village in random.sample(india[:len(india) // 2], len(india) // 4):
 village.set_exogamy_taboo(True, 0.01)

If you look closely, you can in fact tell that places in e.g. Western Africa are lagging behind: In this still, there are as sprinkling of villages in Western Africa where none or only few of the inhabitants have any French ancestry – at a time when French genes have long spread to Eastern and Southern Africa! Similarly, if you watch the animations closely, you’ll see that it takes everyone longer than it should to spread through Papua New Guinea!

Ohne Titel

Some more gifs for purely illustrative purposes:

 

The following is from a beta before I had gotten around labelling the regions and marking the landbridges, but it’s with a larger village size. I think it shows that this increases the speed of spread!

Krasnoyarsk123

 

Here’s the code to initialise the world map (preferrably, this should be done by a parser with a csv as input, but oh well…)


import random, string
import optparse # command line options
from Region import *
from admixturePlotter import *

distance_limit = 21
villagesize = vs = 67

worldmap = [Region(37, 37, distance_limit, name="Europe"),
 Region(40, 40, distance_limit, name="Western_Africa"),
 Region(40, 40, distance_limit, name="Eastern_Africa"),
 Region(47, 47, distance_limit, name="Central_Asia"),
 Region(40, 40, distance_limit, name="Southern_Africa"),
 Region(44, 35, distance_limit, name="Middle_East"),
 Region(20, 20, distance_limit, name="France"),
 Region(15, 15, distance_limit, name="Spain"),
 Region(32, 32, distance_limit, name="India"),
 Region( 10, 20, distance_limit, name="Madagascar"),
 Region( 23, 23, distance_limit, name="South_Central_Siberia"),
 Region( 23, 23, distance_limit, name="North_Central_Siberia"),
 Region( 35, 35, distance_limit, name="Krasnoyarsk"),
 Region(50, 50, distance_limit, name="North_America"),
 Region(42, 42, distance_limit, name="East_Siberia"),
 Region(54, 54, distance_limit, name="Far_East"),
 Region(40, 40, distance_limit, name="South_America"),
 Region(20, 20, distance_limit, name="Carribean"),
 Region(18, 18, distance_limit, name="Alaska"),
 Region(25, 25, distance_limit, name="Greenland"),
 Region(27, 27, distance_limit, name="Southeast_Asia"),
 Region( 15, 15, distance_limit, name="Papua"),
 Region(40, 40, distance_limit, name="Australia"),
 Region(20, 20, distance_limit, name="Chukchi"),
 Region(14, 11, distance_limit, name="Shangri_La"),
 Region( 20, 20, distance_limit, name="Indonesia"),
 Region(17, 17, distance_limit, name="Britain"),
 Region(8, 8, distance_limit, name="Ireland"),
 Region(20, 35, distance_limit, name="Patagonia"), 
 ]

#The default initialisation of a Region() is as a ghost, with 
#size etc. parameters but no actual Villages or Persons created.
#Let's change that! Also, let's make it so we can refer to the 
# continents more easily!

for continent in worldmap:
 continent.populate()
 continent.calculate_distances()
 exec(string.lower(continent.name).replace(" ","_") + " = continent")

  Now lets’ join the regions – the Region() class has a dedicated method for this. The junctiondepth parameter determines how far into the bordering region the villages in the border region look for updating their distance matrixes: e.g a value of 8 (default) means that an Eastern African village 4 steps from Eastern Africa’s south point will only update with villages below 3 distance from the landbridge (i.e. total cross-region distances < 8) even when the distance limit inside a region is, say, 26. gap is only needed for display purposes – it determines the regions’ x and y offset and thus how far apart they will be drawn.

eastern_africa.make_anchor()
eastern_africa.joinregions(southern_africa,"se", "ne", 
  junctiondepth= distance_limit//2, gap=0)
eastern_africa.joinregions(southern_africa, "s")
#southern_africa.joinregions(eastern_africa, "n", 
  junctiondepth= distance_limit // 4)
southern_africa.joinregions(western_africa, "nw", 
  junctiondepth=distance_limit//2, gap=0)
southern_africa.joinregions(eastern_africa, "nw", "sw", 
  junctiondepth=distance_limit//2)
eastern_africa.joinregions(middle_east, "ne", 
  junctiondepth=distance_limit//2)
eastern_africa.joinregions(western_africa, "nw", "ne", 
  junctiondepth=distance_limit//2)
eastern_africa.joinregions(western_africa, "sw", "se", 
  junctiondepth=distance_limit//2)
middle_east.joinregions(europe, "nw", 
  junctiondepth=distance_limit//2)
middle_east.joinregions(shangri_la, "e","nw")
shangri_la.joinregions(india, "se", "n", gap=1)
middle_east.joinregions(india, "se", 
  junctiondepth=2 * distance_limit //3)
europe.joinregions(france, "sw", junctiondepth=distance_limit//2, 
  gap=0)
france.joinregions(spain, "sw", junctiondepth=distance_limit//3, 
  gap=0)
europe.joinregions(britain, "w", "ne", 
  junctiondepth=distance_limit//2, gap=3)
france.joinregions(britain, "n", junctiondepth=distance_limit//3)
britain.joinregions(ireland, "w", junctiondepth=distance_limit//3, 
  gap=2)
spain.joinregions(western_africa, "s", 
  junctiondepth=distance_limit//2)
eastern_africa.joinregions(madagascar, "se", "nw", 
  junctiondepth= distance_limit // 5, gap=3)
india.joinregions(far_east, "ne", junctiondepth=distance_limit//2)
india.joinregions(southeast_asia, "e", 
  junctiondepth=distance_limit//3)
southeast_asia.joinregions(indonesia, "se", "w", 
  junctiondepth = distance_limit//4, gap=5)
indonesia.joinregions(papua, "e", junctiondepth=4, gap=3)
#bridges = bridges[:-1]
#papua.y_offset *= 3
#indonesia.joinregions(papua, "e", junctiondepth=4)
papua.joinregions(australia,"s", "ne",junctiondepth=3, gap=4)
far_east.joinregions(north_central_siberia, "nw", 
  junctiondepth = distance_limit//2, gap=1)
north_central_siberia.joinregions(south_central_siberia, "s")
south_central_siberia.joinregions(east_siberia, "ne", gap=2)
south_central_siberia.joinregions(central_asia, "nw", "e")
south_central_siberia.joinregions(far_east, "ne", "nw")
south_central_siberia.joinregions(far_east, "se", "w")
central_asia.joinregions(europe, "w", "ne")
central_asia.joinregions(europe,"w", "ne")
central_asia.joinregions(middle_east, "s", "ne")
east_siberia.joinregions(krasnoyarsk, "nw", "e")
east_siberia.joinregions(chukchi, "ne", "nw")
chukchi.joinregions(alaska, "e", "nw")
alaska.joinregions(north_america, "ne", "nw")
north_america.joinregions(greenland, "ne", "nw", junctiondepth=4, 
  gap=3)
north_america.joinregions(south_america, "se")
north_america.joinregions(carribean, "e", "nw")
carribean.joinregions(south_america, "se", "n")
krasnoyarsk.joinregions(central_asia, "sw", "n")
south_america.joinregions(patagonia, "sw", "nw", gap=0)
south_america.joinregions(patagonia, "s", "ne")
india.joinregions(southeast_asia, "ne", "nw")
north_central_siberia.joinregions(east_siberia, "ne", "w")
north_central_siberia.joinregions(krasnoyarsk, "ne", "se")
north_central_siberia.joinregions(krasnoyarsk, "nw", "s")
north_central_siberia.joinregions(central_asia, "nw", "ne", 
junctiondepth=distance_limit//3 )
southeast_asia.joinregions(far_east, "ne", "s")
india.joinregions(southeast_asia, "ne", "nw")

 

In order to run the simulation, I use the following loop: The first ten lines are for displaying status to the console. Lines 14-16 do the actual work of producing offspring, lines 17-19 make the newly generated offspring grow up (and the previous adults irrelevant), and the last 3 lines make sure no village becomes totally abandoned by triggering a migration event from the most populous nearby village (distance range 4) if the population has dropped below threshold due to the stochastic nature of the model (below 1/4 of the initial population size).


allcontinents = flattenonce([continent for continent in worldmap])

for generation in range(181):
 i = generation
 print "***"
 print "generation", i
 world_pop_size = 0
 for village in allcontinents:
   world_pop_size += len(village.adults)
 print "current total population", world_pop_size
 print "***"

 random.shuffle(allcontinents)
 for village in allcontinents: 
   for person in village.adults:
     person.spawn(1) 
 for village in allcontinents:
   village.adults = village.children
   village.children = []
 for village in allcontinents:
   if len(village.adults) <= vs // 4:
   village.immigrate()

 

And finally, there’s an nifty module for displaying the results graphically, called like this:


if generation % 5 == 0:
for origin_continent in worldmap:
print "mapping ", origin_continent.name
showstatus_multiregion(
[(continent_drawn, continent_drawn.x_offset,
continent_drawn.y_offset) for continent_drawn
in worldmap
],
origin_continent , i, vs, preferred_marker=","
)

Modules: Region.py, Village.py, Person.py admixtureMapper.py, all at GitHub!

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s