Tomorrow, we will discuss in our cours forecasting tools for demographics. But first, we will see basic static tools. Before playing with longitudinal dataset, let us use “standard” life tables. Some (French) datasets are available on the INED website, but let us use the most popular ones, the TV8890 and TD8890. Those two tables can be downloaded from wikipedia,
url="https://fr.wikipedia.org/wiki/Table_de_mortalit%C3%A9" download.file(url,"mortalite.html") library(XML) tables=readHTMLTable("mortalite.html") |
The code to get the dataset is the following (here, we should be careful since there is a space in the number, as a separator for thousands)
TV8890_0=tables[[2]] a1=as.numeric(as.character(TV8890_0[,1])) a2=as.numeric(as.character(TV8890_0[,3])) espace=function(x) gsub("[[:space:]]", "", x) b1=espace(as.character(TV8890_0[,2])) b2=espace(as.character(TV8890_0[,4])) TV8890=data.frame(x=c(a1,a2),lx=as.numeric(c(b1,b2)) |
One can also read the second life table (note that there is a typo in wikipedia since here, the two tables are exactly the same)
TD8890_0=tables[[1]] a1=as.numeric(as.character(TD8890_0[,1])) a2=as.numeric(as.character(TD8890_0[,3])) b1=espace(as.character(TD8890_0[,2])) b2=espace(as.character(TD8890_0[,4])) TD8890=data.frame(x=c(a1,a2),lx=as.numeric(c(b1,b2))) |
It is possible to use that survival function to compute some sort to life expectancy at birth
sum(TV8890$lx)/100000-1 [1] 72.01518 |
One can visualize the survival probability (up to a 100,000 scaling constant)
plot(TV8890,type="l") |
or the death probability, i.e. the probability to die at some specific age x, given that you did reach age x, also called the force of mortality
n=nrow(TV8890) px=(TV8890$lx[1:(n-1)]-TV8890$lx[2:n])/ TV8890$lx[1:(n-1)] x=TV8890$x[1:(n-1)] plot(x,px,type="l",xlab="age") |
A more popular visualization is obtained with a log scale for the probability
plot(x,px,type="l",log="y") |
Finally, we can compute the density of the age at death
pbx=TV8890$lx[1:(n-1)]*px/100000 plot(x,pbx,type="l") |
that can also be used to compute life expectancy
sum(x*pbx) [1] 72.01518 |
That is for the static case. For longitudinal tables, we can use those from the Human Mortality Database. For instance, for France, we can get information of the number of death at age x during year t, as well as the exposure (number of people alive). For France, there are datasets available here
url="http://freakonometrics.free.fr/FranceDeaths_1x1.txt" download.file(url,"FRD.txt") url="http://freakonometrics.free.fr/FranceExposures_1x1.txt" download.file(url,"FRE.txt") |
and for Canada
url="http://freakonometrics.free.fr/CanadaDeaths_1x1.txt" download.file(url,"CAD.txt") url="http://freakonometrics.free.fr/CanadaExposures_1x1.txt" download.file(url,"CAE.txt") |
The following code can be use to read those files.
FRD=read.table("FRD.txt",skip = 3,header=TRUE) tail(FRD) Year Age Female Male Total 22195 2015 105 297.96 35.86 333.82 22196 2015 106 182.95 20.39 203.34 22197 2015 107 104.87 10.50 115.37 22198 2015 108 57.27 5.07 62.34 22199 2015 109 31.93 2.59 34.52 22200 2015 110+ 33.03 1.61 34.64 |