-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathPop_sim.py
More file actions
388 lines (325 loc) · 13.8 KB
/
Pop_sim.py
File metadata and controls
388 lines (325 loc) · 13.8 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
import numpy as np
from faker import Faker
import random
import matplotlib.pyplot as plt
import networkx as nx
from collections import Counter
fake = Faker()
debug = False
"""
Implement genetics system:
inherit good or bad traits
eugenics pushes for good traits
but small gene pool makes emergence of bad traits more likely
"""
class Person:
ID = 0
def __init__(self, age):
Person.ID += 1
self.ID = Person.ID
self.age = age
self._surname = None
# Define is_child based on age
self.is_child = True if self.age < 18 else False
# 50/50 gender probability
self.gender = random.choice(["male", "female"])
if self.gender == "male":
self.first_name = fake.first_name_male()
elif self.gender == "female":
self.first_name = fake.first_name_female()
self.name = self.first_name
@property
def surname(self):
return self._surname
@surname.setter
def surname(self, value):
self._surname = value
if self._surname is not None:
self.name = self.first_name + " " + self._surname
# Recent stats: 3% of population is homosexual and 4% is bisexual
orientation_prob = random.random()
if orientation_prob < 0.03:
orientation = "homosexual"
elif orientation_prob < 0.07:
orientation = "bisexual"
else:
orientation = "heterosexual"
self.orientation = orientation
if self.orientation == "heterosexual":
if self.gender == "male":
self.spouse_gender = "female"
else:
self.spouse_gender = "male"
elif self.orientation == "homosexual":
if self.gender == "male":
self.spouse_gender = "male"
else:
self.spouse_gender = "female"
else:
self.spouse_gender = ["male", "female"]
if self.age < 18:
self.is_child = True
else:
self.is_child = False
self.spouse = None
self.parents = []
self.children = []
class Family:
def __init__(self, parents, surname):
self.parents = parents
self.surname = surname
self.children = []
def add_child(self, child):
self.children.append(child)
for parent in self.parents:
parent.children.append(child)
def simulate_population(population_size, initial_family_count, debug, time_steps):
fake = Faker()
Faker.seed(0) # Set seed for consistent data generation
population = []
families = []
step = 0
AVG_LIFESPAN = 80 + (step // 10)
# Initialize the population
print("Initializing the population...")
print("Initializing age distribution...")
age_distribution = [random.randint(0, 15) for _ in range(int(population_size * 0.55))] + [random.randint(35, 60) for
_ in range(int(population_size * 0.45))]
surnames = set()
# Populate array of unique surnames
while len(surnames) < population_size / 2:
surname = fake.last_name()
if surname not in surnames:
surnames.add(surname)
if debug:
print(f"{len(surnames)} total surnames generated. ")
for surname in surnames:
print(surname)
print("Initializing raw population...")
debug_counter = 0
for i in range(population_size):
# Initialize the population with age according to the age distribution
age = age_distribution[i]
# Initialize person
person = Person(age)
# Add person to population
population.append(person)
# Give adults a surname
if person.is_child is False:
try:
print(list(surnames))
person.surname = random.choice(list(surnames))
# Remove surname from pool once assigned to assure uniqueness
surnames.remove(person.surname)
except:
print("Error: Not enough surnames to assign to population. Len: " + str(len(surnames)))
exit(1)
else:
person.surname = None
if debug:
debug_counter+=1
print(f"{debug_counter} Surname assigned: " + person.name)
counter = Counter(person.surname for person in population)
duplicates = [item for item, count in counter.items() if count > 1 and item is not None and item is not True]
if debug:
print("duplicates = " + str(duplicates))
if len(duplicates) > 0:
print("Error: Duplicate surnames detected: " + str(len(duplicates)) + ' ' + str(duplicates))
exit(1)
for person in population:
if person.is_child is False:
print(person.name)
# DEBUG: Print population statistics
# if debug is True:
# # for person in population:
# # print(person.name + ' ' + str(person.age))
# max_age = max(person.age for person in population)
# print(f"Minimum age: {min(person.age for person in population)}")
# print(f"Maximum age: {max_age}")
# age_histogram = [0] * (max_age + 1)
# orientation_histogram = [0] * 3
#
# for person in population:
# if person.age <= max_age:
# age_histogram[person.age] += 1
# else:
# age_histogram.append(1)
# max_age += 1
#
# plt.bar(range(len(age_histogram)), age_histogram)
# plt.xlabel('Age')
# plt.ylabel('Population Count')
# plt.title('Population Demographics by Age, Initial Population')
# plt.show()
#
# total_males = 0
# for person in population:
# if person.gender == "male":
# total_males += 1
# print(f"Total males: {total_males}\nTotal females: {500 - total_males}\nGender ratio: {total_males / 500}")
#
# for person in population:
# if person.orientation == "heterosexual":
# orientation_histogram[0] += 1
# elif person.orientation == "homosexual":
# orientation_histogram[1] += 1
# else:
# orientation_histogram[2] += 1
#
# labels = ["heterosexual", "homosexual", "bisexual"]
# x = range(len(orientation_histogram))
# plt.bar(x, orientation_histogram, tick_label=labels)
# plt.show()
# Simulate population dynamics for the given number of time steps (years)
print("Simulating population dynamics...")
max_age = max(person.age for person in population)
age_histogram = [0] * (max_age + 1)
for step in range(time_steps):
print(f"Simulating year {step}...")
# Update ages and track demographics
print("Aging population...")
for person in population:
person.age += 1
if person.age <= max_age:
age_histogram[person.age] += 1
else:
age_histogram.append(1)
max_age += 1
# Determine who dies (based on age and life expectancy)
# Define Gompertz-Makeham parameters
alpha = 0.05 # Baseline mortality rate
beta = 0.001 # Rate of age-related increase in mortality
print("Determining mortality due to old age...")
# Determine who dies based on age and Gompertz-Makeham distribution
for family in families:
for parent in family.parents:
age_difference = parent.age - AVG_LIFESPAN
if age_difference >= 0:
survival_probability = np.exp(alpha + beta * age_difference)
if random.random() > survival_probability:
family.parents.remove(parent)
population.remove(parent)
print(f"{parent.first_name} {parent.surname} died at age {parent.age}")
print("Generating families...")
# Allow couples to form families and have children
eligible_parents = [person for person in population if 20 <= person.age <= 50]
print("Generating parenting couples...")
parents = []
print(", " .join(person.name for person in eligible_parents))
while len(eligible_parents) >= 2:
while True:
print("Choosing first parent...")
parents.append(random.choice(eligible_parents))
print(f"Parent: {parents[0].name}, {parents[0].age}, {parents[0].orientation}\n\n")
print("Eligible parents: " + str(len(eligible_parents)))
while True:
parents.append(random.choice(eligible_parents))
for person in parents:
print(f"Person: {person.name}, {person.age}, {person.orientation}, {person.spouse_gender}")
if debug is True:
if len(parents) == 2:
print("First parent gender:", parents[0].gender)
print("First parent spouse gender:", parents[0].spouse_gender)
print("Second parent gender:", parents[1].gender)
print("Second parent spouse gender:", parents[1].spouse_gender)
if parents[1].gender in parents[0].spouse_gender and parents[0].gender in parents[1].spouse_gender:
print("Eligible\n")
eligible_parents.remove(parents[0])
eligible_parents.remove(parents[1])
break
else:
print("Not eligible\n")
parents.remove(parents[1])
if parents[0].gender in parents[1].spouse_gender and parents[1].gender in parents[0].spouse_gender:
print("Couple formed!")
# If hetero, assign male surname, otherwise doesn't matter
if parents[0].gender == "male":
parents[1].surname = parents[0].surname
else:
parents[0].surname = parents[1].surname
for parent in parents:
print(f"{parent.name}, {parent.age}, {parent.orientation}")
family_surname = parents[0].surname
print("Remaining eligible parents: " + str(len(eligible_parents)))
if len(parents) != 2:
print("Parents need to be two. Resetting...")
parents = []
continue
break
# Assign spouses
if len(parents) != 2:
exit(1)
parents[0].spouse = parents[1]
parents[1].spouse = parents[0]
family = Family(parents, family_surname)
families.append(family)
# Assigned any children without parents to this family, up to 3
print("Assigning children to family...")
for person in population:
# If child without parents, assign to family
if person.is_child is True and person.parents == []:
person.parents = parents
# Patriarchal surname assignment
if parents[0].gender == "male":
person.surname = parents[0].surname
else:
person.surname = parents[1].surname
# Add child to family
family.add_child(person)
print(f"Assigned {person.name} to family {family_surname}")
# If family has 3 children, stop assigning children to this family
if len(family.children) >= 3:
break
# if population capacity not yet reached, have children
if len(population) < population_size:
print("Adding children to the population...")
child_surname = family_surname
child = Person(0, child_surname)
family.add_child(child)
population.append(child)
eligible_parents.append(child)
print(f"{parents[0].name} and {parents[1].name} had a baby named {child.first_name}!")
# Print family lineages
if step % 10 == 0:
print(f"Family Lineages at Year {step}:")
for i, family in enumerate(families):
if i >= 10:
break
print(f"Couples: {', '.join([parent.name for parent in family.parents])}")
for child in family.children:
print(f" Child: {child.name} ({child.gender})")
# Plot age histogram
plt.bar(range(len(age_histogram)), age_histogram)
plt.xlabel('Age')
plt.ylabel('Population Count')
plt.title('Population Demographics by Age')
plt.show()
# View lineage of a person
def view_lineage(person):
G = nx.DiGraph()
# Create nodes for each person
for p in population:
G.add_node(p.name)
# Create edges for family relationships
for family in families:
for child in family.children:
parent1_name = family.parents[0].name
parent2_name = family.parents[1].name
child_name = child.name
G.add_edge(parent1_name, child_name)
G.add_edge(parent2_name, child_name)
# Draw family tree
pos = nx.spring_layout(G)
nx.draw_networkx(G, pos=pos, with_labels=True)
plt.axis('off')
plt.show()
# Select a person and view their lineage
if len(population) > 0:
selected_person = random.choice(population)
view_lineage(selected_person)
# Run the simulation
population_size = 500
initial_family_count = 112
family_size = 4.5
time_steps = 100
simulate_population(population_size, initial_family_count, family_size, time_steps)