# Answering code for the exercises of Chapter 6. options(stringsAsFactors=FALSE) # It's necessary to replace Age's label. T6E1 <- data.frame( Age = c("0-4", "0", "1", "2", "3", "4", "5-9", "10-14", "15-19", "20-24", "25-29", "30-34", "35-39", "40-44", "45-49", "50-54", "55-59", "60-64", "65-69", "70-74", "75-79", "80-84", "85-89", "90+"), Px.Males = c(1571400, 317000, 323000, 325900, 314900, 290500, 1557600, 1947500, 2121200, 1942800, 1708200, 1764700, 1734500, 1417200, 1368500, 1381000, 1382000, 1277400, 1088600, 900100, 578400, 274500, 97200, 32700), Px.Females = c(1491400, 300800, 308400, 309500, 298100, 274600, 1473800, 1847000, 2014500, 1893200, 1683800, 1749100, 1709800, 1396800, 1351900, 1398300, 1445500, 1429100, 1323100, 1233300, 972500, 620300, 297900, 130100), Dx.Males = c(4566, 3914, 287, 149, 121, 95, 391, 546, 1669, 1668, 1409, 1735, 2246, 3280, 5647, 10497, 18820, 27701, 39171, 51908, 52096, 37844, 19875, 9119), Dx.Females = c(3346, 2861, 208, 140, 75, 60, 253, 353, 588, 672, 702, 1079, 1576, 2132, 3639, 6351, 10854, 16897, 24598, 37623, 49866, 54879, 44632, 31655) ) # Exercise 1 # Calculating qx for Females # At first, omitting unnecessary data and unite age-class 1-4 T6E1x <- T6E1[c(2,3,7:24),] T6E1x$Age[2] <- "1-4" # for this, stringsAsFactors=FALSE is needed. T6E1x[2,2:5]<-T6E1x[2,2:5]+colSums(T6E1[4:6,2:5]) # Above process is easier in Excel-like spreadsheet. # Checking the shape of population pyramid library(pyramid) par(cex=0.7) pyramids(Left=T6E1x$Px.Males, Right=T6E1x$Px.Females, Center=T6E1x$Age, Laxis=0:5*500000, AxisFM="d", Cadj=-0.02) # Set n and nax T6E1x$n <- c(1, 4, rep(5,18)) T6E1x$nax <- c(0.1, 0.4, rep(0.5,17), NA) # Calculate ASDR for Femaeles as Mx.Females T6E1x$Mx.Females <- T6E1x$Dx.Females/T6E1x$Px.Females T6E1x$qx.Females <- ifelse(is.na(T6E1x$nax), 1, T6E1x$n*T6E1x$Mx.Females/(1+T6E1x$n*(1-T6E1x$nax)*T6E1x$Mx.Females)) # You can see the resulted table by below # The answers given in the book is wrong (rounding error?) print(T6E1x) # # Exercise 2 data includes lx as Proportion # T6E2 <- data.frame( Month = 0:24, Proportion = c(1.000, 0.988, 0.972, 0.955, 0.937, 0.908, 0.880, 0.844, 0.795, 0.735, 0.680, 0.578, 0.485, 0.392, 0.366, 0.334, 0.306, 0.232, 0.169, 0.108, 0.105, 0.101, 0.099, 0.072, 0.047)) # 6d12 = l6-l12 T6E2$Proportion[T6E2$Month==12]-T6E2$Proportion[T6E2$Month==18] # p20 = l21/l20 T6E2$Proportion[T6E2$Month==21]/T6E2$Proportion[T6E2$Month==20] # l12 T6E2$Proportion[T6E2$Month==12] # 3L20 = 3*(l20+l23)/2 in rough calculation 3*(T6E2$Proportion[T6E2$Month==20]+T6E2$Proportion[T6E2$Month==23])/2 # If you like to calculate accurately, 3L20=L20+L21+L22 T6E2$dx <- T6E2$Proportion - c(T6E2$Proportion[-1],0) T6E2$Lx <- T6E2$Proportion - T6E2$dx/2 T6E2$Lx[T6E2$Month==20]+T6E2$Lx[T6E2$Month==21]+T6E2$Lx[T6E2$Month==22] # 6q6 = 6d6/l6 = (l6-l12)/l6 (T6E2$Proportion[T6E2$Month==6]-T6E2$Proportion[T6E2$Month==12])/T6E2$Proportion[T6E2$Month==6] # Answer = 3L6/3L3*100, a kind of survivorship ratio (Difficult English) sum(T6E2$Lx[T6E2$Month %in% 6:8])/sum(T6E2$Lx[T6E2$Month %in% 3:5])*100 # # Exercise 3 Graunt's life table, data includes x and lx # T6E3 <- data.frame( x = c(0, 6, 16, 26, 36, 46, 56, 66, 76), # Exact Age lx = c(100, 64, 40, 25, 16, 10, 6, 3, 1)) # Survivors # (a) 26q0 = 1 - (l26/l0) 1-T6E3$lx[T6E3$x==26]/T6E3$lx[T6E3$x==0] # (b) 20p26 = l46/l26 T6E3$lx[T6E3$x==46]/T6E3$lx[T6E3$x==26] # (c) e0x requires n, dx, ax, Lx, Tx T6E3$n <- c(6, rep(10, 8)) # Assuming 76+ for 10 yrs # Assuming ax as 0.2 for age 0-6, # 0.5 for age 76+ (it means 10*0.5 years survival for 76+) T6E3$ax <- c(0.2, rep(0.5, 8)) T6E3$dx <- T6E3$lx - c(T6E3$lx[-1], 0) T6E3$Lx <- T6E3$n*(T6E3$lx - (1-T6E3$ax)*T6E3$dx) T6E3$Tx <- rev(cumsum(rev(T6E3$Lx))) T6E3$ex <- T6E3$Tx/T6E3$lx print(T6E3) print(T6E3$ex[1]) # nLx is the number of person-lived between exact ages x and x+n.