# Chapter 12 # coded by Minato Nakazawa # 26 June 2014. T12.1s <- data.frame( EAG = c(0, 1, 1:17*5), L0lx = c(100000, 66769, 48876, 44523, 42044, 39182, 35280, 32226, 28373, 24198, 19733, 15307, 11332, 7900, 5132, 2975, 1399, 483, 90), L115lx = c(100000, 98206, 97911, 97774, 97646, 97430, 97131, 96767, 96282, 95587, 94507, 92748, 89872, 85330, 78320, 68070, 53953, 36999, 20119) ) T12.1s T12.2 <- data.frame( EAG = c(0, 1, 1:20*5), n = c(1, 4, rep(5, 20)), qx.F.W12 = c(132.25, 83.06, 24.28, 18.94, 25.51, 32.41, 36.61, 41.48, 46.27, 51.34, 58.67, 77.72, 102.44, 149.84, 209.17, 305.14, 430.57, 583.37, 749.66, 888.16, 968.73, 1000.00)/1000, mx.F.W12 = c(144.69, 21.97, 4.92, 3.82, 5.16, 6.58, 7.45, 8.46, 9.46, 10.53, 12.07, 16.15, 21.55, 32.29, 46.50, 71.50, 109.74, 169.45, 263.21, 405.52, 619.78, 945.62)/1000, mx.M.W12 = c(174.20, 22.22, 4.78, 3.46, 4.94, 7.03, 7.75, 8.92, 10.60, 13.15, 16.20, 21.58, 28.51, 40.65, 57.45, 83.70, 125.81, 189.47, 287.93, 434.39, 651.01, 974.73)/1000 ) library(fmsb) # qxtomx(T12.2$qx.F.W12, ax=c(0.35, 0.34, rep(0.5, 20)), n=c(1, 4, rep(5,20))) # matches with T12.2$mx.F.W12 for ages from 0 to 75, but are quite different # after age 80. Coale and Demeny must have used adjustment for elderlies. lifetable2(mx=T12.2$mx.F.W12, n=T12.2$n, ax=c(0.35, 0.34, rep(0.5, 20))) lifetable2(mx=T12.2$mx.M.W12, n=T12.2$n, ax=c(0.35, 0.34, rep(0.5, 20))) # Calculated tables are slightly different from the T12.2, probably due to # adjustment in elderlies noted above. # # Exercises # 1. Gambia data (For T12E.2, only WEST model was entered) T12E.1 <- data.frame( Age = c(1/52, 1/12, 1:5), lx = c(94.1, 90.8, 76.1, 64.3, 56.5, 52.3, 50.1)/100 ) T12E.2 <- data.frame( lx.F.W01 = c(63445,54958,51154,48696,46836), lx.F.W02 = c(66601,58514,54891,52549,50776), lx.F.W03 = c(69444,61785,58353,56135,54456), lx.F.W04 = c(72027,64811,61578,59488,57907), lx.F.W05 = c(74389,67625,64593,62634,61152), lx.F.W06 = c(76562,70251,67423,65596,64213), lx.F.W07 = c(78571,72713,70088,68391,67107), lx.F.W08 = c(80438,75028,72604,71037,69852), lx.F.W09 = c(82178,77211,74986,73547,72459), lx.F.W10 = c(83807,79276,77246,75933,74940), lx.F.W11 = c(85336,81233,79394,78206,77307), lx.F.W12 = c(86775,83092,81441,80374,79567), lx.F.W13 = c(88121,84865,83405,82462,81749), lx.F.W14 = c(89396,86646,85413,84616,84013), lx.F.W15 = c(90606,88290,87242,86559,86037), lx.F.W16 = c(91769,89864,88987,88407,87954), lx.F.W17 = c(92884,91352,90635,90153,89772), lx.F.W18 = c(93949,92759,92192,91806,91496), lx.F.W19 = c(94965,94089,93664,93372,93134), lx.F.W20 = c(95931,95347,95059,94859,94693), lx.F.W21 = c(96884,96531,96355,96231,96127), lx.F.W22 = c(97718,97507,97400,97324,97260), lx.F.W23 = c(98470,98361,98305,98264,98230), lx.F.W24 = c(99095,99048,99024,99007,98992), lx.F.W25 = c(99555,99540,99533,99527,99522), lx.M.W01 = c(58050,50262,46851,44617,42957), lx.M.W02 = c(61614,54105,50817,48663,47062), lx.M.W03 = c(64826,57643,54497,52437,50906), lx.M.W04 = c(67743,60918,57929,55972,54517), lx.M.W05 = c(70411,63965,61142,59293,57919), lx.M.W06 = c(72865,66812,64160,62424,61133), lx.M.W07 = c(75135,69481,67004,65382,64177), lx.M.W08 = c(77243,71992,69691,68185,67066), lx.M.W09 = c(79209,74360,72237,70846,69813), lx.M.W10 = c(81049,76601,74654,73378,72430), lx.M.W11 = c(82775,78726,76953,75791,74928), lx.M.W12 = c(84401,80745,79144,78096,77316), lx.M.W13 = c(85983,82816,81428,80520,79844), lx.M.W14 = c(87487,84756,83560,82777,82194), lx.M.W15 = c(88804,86446,85414,84738,84235), lx.M.W16 = c(90084,88086,87208,86632,86203), lx.M.W17 = c(91322,89716,88976,88477,88098), lx.M.W18 = c(92517,91266,90662,90244,89921), lx.M.W19 = c(93666,92736,92266,91933,91672), lx.M.W20 = c(94767,94129,93791,93547,93353), lx.M.W21 = c(95866,95460,95236,95070,94937), lx.M.W22 = c(96901,96648,96501,96391,96302), lx.M.W23 = c(97838,97699,97616,97552,97499), lx.M.W24 = c(98652,98588,98548,98517,98492), lx.M.W25 = c(99289,99266,99252,99240,99231), lx.M.S03 = c(71056,59951,54829,52083,50567), lx.M.S04 = c(73058,62645,57842,55267,53846), lx.F.S03 = c(73567,61479,55913,52974,51304) ) # Taking average of males and females life tables with weighing by # sex ratio is unnecessary. Just find best-fit model among all models # including both males and females. Fitting is usually done by # minimizing the sum of squared residuals. (MSE <- colSums((T12E.2 - c(76100,64300,56500,52300,50100))^2)) TCOL <- which.min(MSE) (Tlx <- T12E.2[,TCOL]) print(names(T12E.2)[TCOL]) # best fit lx in the lifetables given in T12E.2 T12E.3 <- data.frame( EAG = c(0, 1, 1:19*5), n = c(1, 4, rep(5, 19)), mx.F.S03 = c(319.16, 95.63, 15.44, 8.12, 11.21, 14.05, 15.14, 15.80, 16.67, 17.02, 18.08, 23.72, 32.68, 54.37, 84.46, 135.11, 209.40, 297.11, 431.19, 622.46, 898.27) / 1000, mx.M.S03 = c(359.07, 89.99, 13.77, 6.50, 9.83, 14.94, 14.95, 14.86, 16.34, 19.45, 22.85, 29.13, 38.79, 57.69, 85.10, 131.36, 200.38, 284.46, 412.23, 594.36, 856.77) ) (LT <- lifetable2(mx=T12E.3[,paste("mx.",substr(names(T12E.2)[TCOL],4,8),sep="")], n=T12E.3$n, ax=c(0.35, 0.3, rep(0.5, 19)))) (P0 <- sum(LT$Lx[1:2])/(5*LT$lx[1])) (P1 <- LT$Lx[3]/sum(LT$Lx[1:2])) # e(x) and P(x) are slightly different from the given table, but almost fit well. # e(0) = 25.0, e(5) = 42.5 # P(1) = 0.785, which means the probability of surviving from 0-4 to 5-9. # Possible reason of poor estimates are unregistered stillbirths, fitting was # only based on childhood mortality. Childhood mortality is very high in this data. # # 2. Iceland life table for males in 1856-65, which model may fit? T12E.4 <- data.frame( EAG = c(1, 5, 10, 15), n = c(4, 5, 5, 5), Observed = c(0.832, 0.776, 0.761, 0.745), West = c(0.845, 0.775, 0.756, 0.743), North = c(0.868, 0.783, 0.751, 0.734), East = c(0.829, 0.774, 0.759, 0.751), South = c(0.856, 0.771, 0.756, 0.745) ) (MSE <- colSums((T12E.4[,4:7]-T12E.4[,3])^2)) names(T12E.4)[which.min(MSE)+3] T12E.4[,which.min(MSE)+3] # East model fitted best, probably due to the level # of infant mortality. Some kind of infectious disease # outbreak may have occurred in 19C Iceland.