yule.py

yule.py

#!/usr/bin/env python3

import matplotlib.pyplot as plt
import numpy as np

NUM_CITIES = 1000
city_pops = [1] * NUM_CITIES  # initialize NUM_CITIES cities each with population 1

for person in range(100000):
    prob_vector = np.divide(city_pops, sum(city_pops))
    chosen_city = np.random.choice(list(range(NUM_CITIES)), size=1, p=prob_vector)[0]
    city_pops[chosen_city] += 1

print(city_pops)
plt.hist(city_pops, bins=100)
plt.show()