Skip to content

Instantly share code, notes, and snippets.

@tbates
Created November 20, 2012 13:11
Show Gist options
  • Save tbates/4117872 to your computer and use it in GitHub Desktop.
Save tbates/4117872 to your computer and use it in GitHub Desktop.
persons = rnorm(100, 20,10)
nYearsToRun = 500
births = deaths = pop = rep(NA,nYearsToRun)
for (y in 1:nYearsToRun) {
populationSize = length(persons)
persons = persons+1 # age everybody
# death
persons <- persons[persons < rnorm(populationSize, 77, 15)]
nDeaths = populationSize - length(persons)
populationSize = length(persons)
# births
giveBirth = rbinom(populationSize,1, .02) # .1 chance of giving birth per year
newPeople <- persons[giveBirth & persons < 60 ]
nBirths = length(newPeople)
newPeople <- rep(0, nBirths)
persons = c(persons,newPeople)
populationSize = length(persons)
births[y] = nBirths
deaths[y] = nDeaths
pop[y] = populationSize
}
x = c(1:nYearsToRun)
plot(pop ~ x, ylim=c(0,max(pop)))
lines(x,deaths,col="red")
lines(x,births,col="green")
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment