Canada population projections

population.cc

#include <list>
#include <iostream>
#include <boost/random/linear_congruential.hpp>
#include <boost/random/variate_generator.hpp>
#include <boost/generator_iterator.hpp>
#include <boost/random/uniform_real.hpp>

using namespace std;

#define MALE 0
#define FEMALE 1


#define SCALE 1

struct person
{
	int sex;
	int birth_year;
};


const int male_population_by_age[] = {
	194329,
	192745,
	188788,
	186052,
	183504,
	184901,
	182492,
	183194,
	185098,
	189710,
	190557,
	194152,
	201893,
	208964,
	215639,
	218032,
	223270,
	229031,
	237774,
	244779,
	240613,
	233174,
	234773,
	240165,
	243590,
	241212,
	238223,
	235928,
	236425,
	233312,
	227319,
	224145,
	225343,
	226494,
	227051,
	223023,
	225595,
	230065,
	239105,
	240998,
	238554,
	237441,
	242316,
	256627,
	275416,
	285257,
	287131,
	278327,
	278264,
	272983,
	266144,
	263181,
	257843,
	249121,
	246035,
	235919,
	224684,
	215556,
	210930,
	205526,
	200107,
	198858,
	198253,
	169793,
	158979,
	153943,
	146575,
	135384,
	128065,
	117608,
	112022,
	105820,
	99717,
	97126,
	92663,
	88343,
	86347,
	83202,
	78812,
	72666,
	65981,
	60946,
	54551,
	49816,
	44393,
	38762,
	33149,
	28964,
	24177,
	18877,
	13477,
	10490,
	8207,
	6436,
	5116,
	3653,
	2556,
	1679,
	1037,
	702
};

const int female_population_by_age[] = {
	185630,
	184279,
	179636,
	175835,
	172900,
	174365,
	172481,
	172674,
	174064,
	179672,
	180937,
	185080,
	190781,
	200241,
	205124,
	207897,
	213672,
	220285,
	225876,
	230536,
	226465,
	221089,
	223357,
	227364,
	230356,
	231544,
	231507,
	232597,
	234353,
	231525,
	227496,
	224908,
	225176,
	225004,
	225827,
	221342,
	222643,
	226002,
	235395,
	235425,
	234573,
	234955,
	238517,
	251969,
	271554,
	280100,
	282155,
	275441,
	276421,
	272075,
	265818,
	263221,
	259357,
	252511,
	250806,
	242297,
	231257,
	221001,
	216827,
	211785,
	207090,
	205935,
	205667,
	177245,
	165840,
	160947,
	154352,
	143852,
	137920,
	128231,
	123170,
	117750,
	112874,
	111312,
	108249,
	103354,
	103012,
	100695,
	98401,
	94118,
	88040,
	84621,
	79830,
	76396,
	72128,
	66237,
	59785,
	54192,
	47776,
	40378,
	31535,
	25557,
	21544,
	17651,
	14945,
	11343,
	8404,
	6113,
	4219,
	2813
};

const int births_by_mother_age[] = {
	/* Age of mother, 0  to 9 years  */ 0,
	/* Age of mother, 10 to 15 years */	104,
	/* Age of mother, 15 to 19 years */	15534,
	/* Age of mother, 20 to 24 years */	57778,
	/* Age of mother, 25 to 29 years */	116878,
	/* Age of mother, 30 to 34 years */	120734,
	/* Age of mother, 35 to 39 years */	57733,
	/* Age of mother, 40 to 44 years */	11364,
	/* Age of mother, 45 to 49 years */	605,
	/* Age of mother, not stated */		133
};

const int births_by_mother_age_range_size[] = {
	10, 5, 5, 5, 5, 5, 5, 5, 5, 50
};


const int male_deaths_by_age[] = {
	145,
	55,
	37,
	27,
	28,
	16,
	20,
	17,
	24,
	19,
	24,
	26,
	32,
	46,
	56,
	67,
	95,
	148,
	175,
	205,
	194,
	188,
	185,
	207,
	171,
	193,
	180,
	187,
	191,
	193,
	171,
	179,
	188,
	195,
	224,
	230,
	231,
	253,
	299,
	319,
	349,
	346,
	395,
	428,
	538,
	575,
	689,
	726,
	770,
	857,
	927,
	993,
	1086,
	1165,
	1238,
	1273,
	1388,
	1397,
	1448,
	1605,
	1764,
	1822,
	2000,
	1937,
	1982,
	2012,
	2198,
	2266,
	2394,
	2319,
	2432,
	2522,
	2605,
	2838,
	2955,
	3120,
	3408,
	3559,
	3839,
	3911,
	3850,
	3954,
	3916,
	4037,
	3926,
	3809,
	3585,
	3586,
	3267,
	2811,
	2253,
	1854,
	1543,
	1351,
	1144,
	853,
	652,
	459,
	308,
	195
};

const int female_deaths_by_age[] = {
	101,
	46,
	25,
	21,
	20,
	23,
	23,
	20,
	22,
	15,
	14,
	14,
	19,
	27,
	42,
	46,
	52,
	64,
	65,
	79,
	67,
	65,
	57,
	66,
	67,
	81,
	71,
	75,
	92,
	85,
	85,
	83,
	104,
	94,
	112,
	104,
	126,
	152,
	168,
	197,
	175,
	231,
	246,
	276,
	331,
	418,
	459,
	473,
	512,
	598,
	592,
	622,
	703,
	714,
	819,
	804,
	841,
	914,
	957,
	1033,
	1081,
	1175,
	1292,
	1227,
	1271,
	1329,
	1416,
	1447,
	1555,
	1575,
	1672,
	1805,
	1824,
	1974,
	2192,
	2267,
	2506,
	2861,
	2979,
	3283,
	3318,
	3660,
	3849,
	4144,
	4555,
	4607,
	4775,
	4910,
	4707,
	4517,
	3987,
	3556,
	3274,
	2979,
	2777,
	2236,
	1835,
	1497,
	1099,
	783
};

int person_age(const struct person &p, int year)
{
	int age = year - p.birth_year;

	if (age < 0) age = 0;
	if (age > 99) age = 99;

	return age;
}

double birth_probability(const struct person &p, int year)
{
	double probability = 0.0;

	if (p.sex == FEMALE) {
		int birth_rate;
		int age = person_age(p, year);

		if      (age <= 9)  birth_rate = births_by_mother_age[0] / births_by_mother_age_range_size[0];
		else if (age <= 14) birth_rate = births_by_mother_age[1] / births_by_mother_age_range_size[1];
		else if (age <= 19) birth_rate = births_by_mother_age[2] / births_by_mother_age_range_size[2];
		else if (age <= 24) birth_rate = births_by_mother_age[3] / births_by_mother_age_range_size[3];
		else if (age <= 29) birth_rate = births_by_mother_age[4] / births_by_mother_age_range_size[4];
		else if (age <= 34) birth_rate = births_by_mother_age[5] / births_by_mother_age_range_size[5];
		else if (age <= 39) birth_rate = births_by_mother_age[6] / births_by_mother_age_range_size[6];
		else if (age <= 44) birth_rate = births_by_mother_age[7] / births_by_mother_age_range_size[7];
		else if (age <= 49) birth_rate = births_by_mother_age[8] / births_by_mother_age_range_size[8];
		else                birth_rate = births_by_mother_age[9] / births_by_mother_age_range_size[9];

		probability = (double)birth_rate / (double)female_population_by_age[age];

	}

	return probability;
}


double death_probability(const struct person &p, int year)
{
	double probability = 0.0;
	int age = person_age(p, year);

	if (p.sex == MALE) {
		probability = (double)male_deaths_by_age[age] / (double)male_population_by_age[age];
	} else {
		probability = (double)female_deaths_by_age[age] / (double)female_population_by_age[age];
	}

	return probability;
}



int main(int argc, char **argv)
{
	int seed = atoi(argv[1]);
	const int start_year = 2009;
	const int end_year   = 2075;

	list<struct person> population;

	boost::minstd_rand generator(seed);
	boost::uniform_real<> uniform_distribution(0, 1);
	boost::variate_generator<boost::minstd_rand&, boost::uniform_real<> > rand(generator, uniform_distribution);

	// Set up an initial population
	for (int age = 0; age < 100; ++age) {

		// Make a bunch of men
		for (int i = 0; i < male_population_by_age[age] / SCALE; ++i) {
			struct person p;
			p.sex = MALE;
			p.birth_year = start_year - age;

			population.push_back(p);
		}

		// Make a bunch of women
		for (int i = 0; i < female_population_by_age[age] / SCALE; ++i) {
			struct person p;
			p.sex = FEMALE;
			p.birth_year = start_year - age;
			
			population.push_back(p);
		}
	}

	// Run the simulation
	for (int year = start_year; year <= end_year; ++year) {
		int total_population = 0;
		int total_children = 0;
		int total_working_age = 0;
		int total_retired = 0;
		int births = 0;
		int deaths = 0;

		// Update the population
		list<struct person>::iterator person_iter = population.begin();
		bool advance;


		while(person_iter != population.end())
		{
			advance = true;
			total_population += SCALE;

			int age = person_age(*person_iter, year);
			if      (age < 18) total_children    += SCALE;
			else if (age < 67) total_working_age += SCALE;
			else               total_retired     += SCALE;

			double birth_variate = rand();
			double death_variate = rand();
			double sex_variate   = rand();

/*
			cout << "person:"
				 << " sex=" << person_iter->sex
				 << " birth_year=" << person_iter->birth_year
				 << " age=" << person_age(*person_iter, year)
				 << endl;

			cout << "variates:"
				 << " birth=" << birth_variate
 				 << " death=" << death_variate
 				 << " sex="   << sex_variate
 				 << endl;

 			cout << "probabalities:"
				 << " birth=" << birth_probability(*person_iter, year)
 				 << " death=" << death_probability(*person_iter, year)
 				 << endl;
*/

			if (birth_variate < birth_probability(*person_iter, year)) {
				struct person child;
				child.sex = sex_variate < 0.5 ? MALE : FEMALE; // TODO: Make this use real data
				child.birth_year = year;

				population.push_back(child);
				births += SCALE;
			}

			if (death_variate < death_probability(*person_iter, year)) {
				person_iter = population.erase(person_iter);
				advance = false;
				deaths += SCALE;
			}

			if (advance) person_iter++;
		}

		cout << year << "\t"
		     << total_population << "\t"
   		     << total_children << "\t"
   		     << total_working_age << "\t"
   		     << total_retired << "\t"
		     << births << "\t"
		     << deaths 
		     << endl;
	}

}