learn.R

5. Plotting a Lexis Surface

# Here, we will use the R function filled.contour() to plot a surface.
# Specifically, we'll use data from the human mortality database to
# construct central death rates for the United States, 1959-1999, and
# plot a Lexis surface of the logarithm of these rates.

# I assume that you've downloaded the necessary data from HMD. I have
# called them "usa_deaths_1x1.txt" and "usa_exposures_1x1.txt" and
# placed them in my working directory -- you will need to insert
# proper path as necessary

# read the data into the R workspace

usa.death <- read.table("usa_deaths_1x1.txt", header=TRUE, skip=1)
usa.expose <- read.table("usa_exposures_1x1.txt", header=TRUE, skip=1)

# construct the central death rates and build a data frame to hold
# them along with age

fmx <- data.frame(Year=usa.death$Year, Mx=usa.death$Female/usa.expose$Female)

# In order to make a surface plot, you need three things: (1) a 1xk
# vector of x-axis values (e.g., age), (2) a 1xp vector of y-axis
# variables (e.g., years), and (3) a kxp matrix of the thing you're
# interested in (log(Mx)).

# The data are listed rectangularly, with observations for ages 0:110
# for each year stacked on top of each other. We thus need to pull
# out each year's mortality rates and stack them next to each other to
# make the z matrix of plotting heights

# initialize a matrix of the correct size

zz <- matrix(0,nrow=111,ncol=41)

# insert the columns

for(i in 0:40){
zz[,i+1] <- fmx$Mx[fmx$Year==i+1959]
}

age <- 0:110
year <- 1959:1999

# need to take transpose of matrix zz to orient it in a standard Lexis
# format (of course, you could have simply constructed your matrix zz
# row-wise rather than column-wise)

filled.contour(x=year,y=age,z=log(t(zz)),col=topo.colors(20),xlab="Year", ylab="Age")

# I used topo.colors. You can try other color schemes. Some build-in
# ones include: heat.colors(n), terrain.colors(n), topo.colors(n),
# cm.colors(n).

Posted by jhj1 at August 19, 2005 09:56 PM

Comments

Post a comment